-- FreeFem++ v4.6 (Thu Apr 2 15:47:38 CEST 2020 - git v4.6) Load: lg_fem lg_mesh lg_mesh3 eigenvalue 1 : // Discussion: 2 : // 3 : // This example solves the Poisson equation repeatedly, using each 4 : // solution to adapt the mesh before the next solution is computed. 5 : // 6 : // Location: 7 : // 8 : // http://people.sc.fsu.edu/~jburkardt/freefem_src/poisson_adaptive/poisson_adaptive.edp 9 : // 10 : // Licensing: 11 : // 12 : // This code is distributed under the GNU LGPL license. 13 : // 14 : // Modified: 15 : // 16 : // 23 June 2015 17 : // 18 : // Reference: 19 : // 20 : // Frederic Hecht, 21 : // Freefem++, 22 : // Third Edition, version 3.22 23 : // 24 : cout << "\n"; 25 : cout << "poisson_adaptive\n"; 26 : cout << " FreeFem++ version\n"; 27 : cout << " Solve the Poisson equation in the L-shaped regi ... : on\n"; 28 : cout << " using adaptive refinement to reduce the error.\n"; 29 : // 30 : // Define the border. 31 : // 32 : border aaa (t=0,1) {x=t;y=0;}; 33 : border bbb (t=0,0.5) {x=1;y=t;}; 34 : border ccc (t=0,0.5) {x=1-t;y=0.5;}; 35 : border ddd (t=0.5,1) {x=0.5;y=t;}; 36 : border eee (t=0.5,1) {x=1-t;y=1;}; 37 : border fff (t=0,1) {x=0;y=1-t;}; 38 : // 39 : // Define Th, the triangulation of the domain to the left of the border. 40 : // 41 : mesh Th = buildmesh ( aaa(6) + bbb(4) + ccc(4) + ddd(4) + eee(4) + fff(6) ); 42 : // 43 : // Plot the initial mesh. 44 : // 45 : plot ( Th, wait = 1, ps = "poisson_adaptive_mesh_0.ps" ); 46 : // 47 : // Define Vh, the finite element space defined over Th, using P1 basis functions. 48 : // 49 : fespace Vh ( Th, P1 ); 50 : // 51 : // Define U and V as elements of the space VH. 52 : // Moreover, U is the zero function. 53 : // 54 : Vh u = 0.0; 55 : Vh v; 56 : // 57 : // Define F, the PDE right hand side function, and 58 : // G, used to define the Dirichlet boundary condition. 59 : // 60 : func f = 1.0; 61 : func g = 0.0; 62 : // 63 : // Define the weak form of the problem. 64 : // 65 : problem Problem1 ( u, v, solver = CG, eps = -1.0e-6 ) 66 : = int2d (Th) ( dx(u)*dx(v) + dy(u)*dy(v) ) 67 : + int2d (Th) ( v * f ) 68 : + on ( aaa, bbb, ccc, ddd, eee, fff, u = g ); 69 : // 70 : // Carry out an iteration in which the mesh is adaptively refined. 71 : // 72 : real error = 0.1; 73 : real coef = 0.1^(1./5.); 74 : 75 : for ( int it = 0; it < 10; it++ ) 76 : { 77 : Problem1; 78 : plot ( u, wait = 1, ps = "poisson_adaptive_uh_"+it+".ps" ); 79 : Th = adaptmesh ( Th, u, inquire = 1, err = error ); 80 : plot ( Th, wait = 1, ps = "poisson_adaptive_mesh_"+(it+1)+".ps" ); 81 : error = error * coef; 82 : } 83 : // 84 : // Terminate. 85 : // 86 : cout << "\n"; 87 : cout << "poisson_adaptive:\n"; 88 : cout << " Normal end of execution.\n"; 89 : sizestack + 1024 =1776 ( 752 ) poisson_adaptive FreeFem++ version Solve the Poisson equation in the L-shaped region using adaptive refinement to reduce the error. -- mesh: Nb of Triangles = 86, Nb of Vertices 58 GC: converge in 13 g=2.46767e-13 rho= 1.05789 gamma= 0.0950509 -- Solve : min -0.034581 max -1.41269e-62 number of required edges : 0 -- adaptmesh Regulary: Nb triangles 232 , h min 0.0395624 , h max 0.170766 area = 0.75 , M area = 99.8584 , M area/( |Khat| nt) 0.994022 infiny-regularity: min 0.489413 max 1.81052 anisomax 4.47211, beta max = 1.38412 min 0.794756 -- mesh: Nb of Triangles = 232, Nb of Vertices 140 GC: converge in 25 g=4.37019e-13 rho= 1.25989 gamma= 0.304349 -- Solve : min -0.0365388 max -1.93471e-33 number of required edges : 0 -- adaptmesh Regulary: Nb triangles 294 , h min 0.0277197 , h max 0.195416 area = 0.75 , M area = 145.823 , M area/( |Khat| nt) 1.14546 infiny-regularity: min 0.48192 max 1.9535 anisomax 4.78154, beta max = 1.21437 min 0.725526 -- mesh: Nb of Triangles = 294, Nb of Vertices 175 GC: converge in 28 g=7.50518e-13 rho= 1.24029 gamma= 0.249857 -- Solve : min -0.0368433 max 9.21779e-35 number of required edges : 0 -- adaptmesh Regulary: Nb triangles 439 , h min 0.01575 , h max 0.191932 area = 0.75 , M area = 215.717 , M area/( |Khat| nt) 1.1348 infiny-regularity: min 0.416777 max 2.23598 anisomax 7.0479, beta max = 1.26955 min 0.684531 -- mesh: Nb of Triangles = 439, Nb of Vertices 251 GC: converge in 37 g=9.24017e-13 rho= 1.18883 gamma= 0.32351 -- Solve : min -0.0370816 max 1.34934e-33 number of required edges : 0 -- adaptmesh Regulary: Nb triangles 680 , h min 0.00836428 , h max 0.214321 area = 0.75 , M area = 330.466 , M area/( |Khat| nt) 1.12232 infiny-regularity: min 0.412082 max 2.11116 anisomax 10.0641, beta max = 1.26174 min 0.67878 -- mesh: Nb of Triangles = 680, Nb of Vertices 378 GC: converge in 47 g=9.08359e-13 rho= 1.62367 gamma= 0.63713 -- Solve : min -0.0371607 max 1.90832e-33 number of required edges : 0 -- adaptmesh Regulary: Nb triangles 1140 , h min 0.00401219 , h max 0.160219 area = 0.75 , M area = 529.765 , M area/( |Khat| nt) 1.07319 infiny-regularity: min 0.396418 max 2.01354 anisomax 9.96459, beta max = 1.32451 min 0.681254 -- mesh: Nb of Triangles = 1140, Nb of Vertices 617 GC: converge in 58 g=6.30588e-13 rho= 1.57927 gamma= 0.599056 -- Solve : min -0.0372259 max 2.72104e-34 number of required edges : 0 -- adaptmesh Regulary: Nb triangles 1813 , h min 0.00220022 , h max 0.149013 area = 0.75 , M area = 851.221 , M area/( |Khat| nt) 1.08429 infiny-regularity: min 0.383251 max 2.29233 anisomax 10.7669, beta max = 1.32473 min 0.677053 -- mesh: Nb of Triangles = 1813, Nb of Vertices 968 GC: converge in 73 g=9.17893e-13 rho= 1.66321 gamma= 0.645302 -- Solve : min -0.037282 max 1.17912e-34 number of required edges : 0 -- adaptmesh Regulary: Nb triangles 2836 , h min 0.000982267 , h max 0.146597 area = 0.75 , M area = 1344.3 , M area/( |Khat| nt) 1.09469 infiny-regularity: min 0.354304 max 2.02686 anisomax 11.1307, beta max = 1.297 min 0.68333 -- mesh: Nb of Triangles = 2836, Nb of Vertices 1496 GC: converge in 94 g=9.07393e-13 rho= 1.7615 gamma= 0.85011 -- Solve : min -0.0373127 max 5.94154e-35 number of required edges : 0 -- adaptmesh Regulary: Nb triangles 4509 , h min 0.000565851 , h max 0.130189 area = 0.75 , M area = 2135.14 , M area/( |Khat| nt) 1.09357 infiny-regularity: min 0.383775 max 2.05661 anisomax 16.0154, beta max = 1.32215 min 0.679605 -- mesh: Nb of Triangles = 4509, Nb of Vertices 2355 GC: converge in 108 g=9.22451e-13 rho= 1.64455 gamma= 0.789363 -- Solve : min -0.0373178 max 4.91373e-35 number of required edges : 0 -- adaptmesh Regulary: Nb triangles 7029 , h min 0.000305076 , h max 0.118647 area = 0.75 , M area = 3383.79 , M area/( |Khat| nt) 1.11176 infiny-regularity: min 0.422411 max 2.3355 anisomax 15.8433, beta max = 1.30809 min 0.689281 -- mesh: Nb of Triangles = 7029, Nb of Vertices 3643 GC: converge in 133 g=9.72059e-13 rho= 1.70232 gamma= 0.93504 -- Solve : min -0.0373255 max 2.0637e-35 number of required edges : 0 -- adaptmesh Regulary: Nb triangles 11325 , h min 0.000142597 , h max 0.0972546 area = 0.75 , M area = 5359.23 , M area/( |Khat| nt) 1.09286 infiny-regularity: min 0.37457 max 2.22639 anisomax 22.1658, beta max = 1.4393 min 0.699301 -- mesh: Nb of Triangles = 11325, Nb of Vertices 5826 poisson_adaptive: Normal end of execution. times: compile 0.004742s, execution 0.617774s, mpirank:0 CodeAlloc : nb ptr 3703, size :481696 mpirank: 0 Ok: Normal End