Source code for pyro.multigrid.examples.mg_test_simple
#!/usr/bin/env python3"""an example of using the multigrid class to solve Laplace's equation. Here, wesolve:: u_xx + u_yy = -2[(1-6x**2)y**2(1-y**2) + (1-6y**2)x**2(1-x**2)] u = 0 on the boundarythis is the example from page 64 of the book `A Multigrid Tutorial, 2nd Ed.`The analytic solution is u(x,y) = (x**2 - x**4)(y**4 - y**2)"""importmatplotlib.pyplotaspltimportnumpyasnpimportpyro.util.io_pyroasiofrompyro.multigridimportMGfrompyro.utilimportcompare,msg# the analytic solution
[docs]deftest_poisson_dirichlet(N,store_bench=False,comp_bench=False,bench_dir="tests/",make_plot=False,verbose=1,rtol=1e-12):# test the multigrid solvernx=Nny=nx# create the multigrid objecta=MG.CellCenterMG2d(nx,ny,xl_BC_type="dirichlet",yl_BC_type="dirichlet",xr_BC_type="dirichlet",yr_BC_type="dirichlet",verbose=verbose)# initialize the solution to 0a.init_zeros()# initialize the RHS using the function frhs=f(a.x2d,a.y2d)a.init_RHS(rhs)# solve to a relative tolerance of 1.e-11a.solve(rtol=1.e-11)# alternately, we can just use smoothing by uncommenting the following# a.smooth(a.nlevels-1,50000)# get the solutionv=a.get_solution()# compute the error from the analytic solutionb=true(a.x2d,a.y2d)e=v-bprint(" L2 error from true solution = %g\n rel. err from previous cycle = %g\n num. cycles = %d"%(e.norm(),a.relative_error,a.num_cycles))# plot itifmake_plot:plt.figure(num=1,figsize=(5.0,5.0),dpi=100,facecolor='w')plt.imshow(np.transpose(v[a.ilo:a.ihi+1,a.jlo:a.jhi+1]),interpolation="nearest",origin="lower",extent=[a.xmin,a.xmax,a.ymin,a.ymax])plt.xlabel("x")plt.ylabel("y")print("Saving figure to mg_test.png")plt.savefig("mg_test.png")# store the output for later comparisonbench="mg_poisson_dirichlet"my_data=a.get_solution_object()ifstore_bench:my_data.write(f"{bench_dir}/{bench}")# do we do a comparison?ifcomp_bench:compare_file=f"{bench_dir}/{bench}"msg.warning("comparing to: %s "%(compare_file))bench_data=io.read(compare_file)result=compare.compare(my_data,bench_data,rtol)ifresult==0:msg.success(f"results match benchmark to within relative tolerance of {rtol}\n")else:msg.warning("ERROR: "+compare.errors[result]+"\n")returnresultreturnNone