-- FreeFem++ v4.6 (Thu Apr 2 15:47:38 CEST 2020 - git v4.6) Load: lg_fem lg_mesh lg_mesh3 eigenvalue 1 : // Location: 2 : // 3 : // http://people.sc.fsu.edu/~jburkardt/freefem_src/poisson_square/poisson_square.edp 4 : // 5 : // Discussion: 6 : // 7 : // Solve the steady Poisson equation in the unit square. 8 : // 9 : // - uxx - uyy = f(x,y) 10 : // 11 : // with boundary condition 12 : // 13 : // u = 0 14 : // 15 : // Licensing: 16 : // 17 : // This code is distributed under the GNU LGPL license. 18 : // 19 : // Modified: 20 : // 21 : // 12 February 2021 22 : // 23 : // Reference: 24 : // 25 : // John Chrispell, Jason Howell, 26 : // Finite Element Approximation of Partial Differential Equations 27 : // Using FreeFem++, or, How I Learned to Stop Worrying and Love 28 : // Numerical Analysis. 29 : // 30 : // Frederic Hecht, 31 : // New development in FreeFem++, 32 : // Journal of Numerical Mathematics, 33 : // Volume 20, Number 3-4, 2012, pages 251-265. 34 : // 35 : cout << "\n"; 36 : cout << "poisson_square:\n"; 37 : cout << " FreeFem++ version.\n"; 38 : cout << " Solve the Poisson equation in a square region.\n"; 39 : // 40 : // n: number of nodes in each spatial direction. 41 : // 42 : int n = 16 + 1; 43 : // 44 : // h: spacing between nodes. 45 : // 46 : real h = 1.0 / ( n - 1 ); 47 : // 48 : // Th: the triangulation of the unit square. 49 : // 50 : mesh Th = square ( n, n ); 51 : // 52 : // Define Vh, the finite element space. 53 : // Here, "P1" requests piecewise linear functions. 54 : // 55 : fespace Vh ( Th, P1 ); 56 : // 57 : // uh: the trial function that approximates the solution. 58 : // 59 : Vh uh; 60 : // 61 : // vh: the test functions that multiply the equation. 62 : // 63 : Vh vh; 64 : // 65 : // f: the right hand side function f(x,y). 66 : // 67 : func f = 68 : + sin(5*pi*x*(1-x)) * sin(4*pi*y*(1-y)) * pow(5*pi*(1-x)-5*pi*x,2) 69 : + 10 * cos(5*pi*x*(1-x)) * pi * sin(4*pi*y*(1-y)) 70 : + sin(5*pi*x*(1-x)) * sin(4*pi*y*(1-y)) * pow(4*pi*(1-y)-4*pi*y,2) 71 : + 8 * sin(5*pi*x*(1-x)) * pi * cos(4*pi*y*(1-y)); 72 : // 73 : // g: the boundary condition function g(x,y). 74 : // 75 : func g = 0.0; 76 : // 77 : // define "poisson()" to be the problem to be solved. 78 : // 79 : problem poisson ( uh, vh ) 80 : = int2d ( Th ) 81 : ( 82 : dx(uh) * dx(vh) 83 : + dy(uh) * dy(vh) 84 : ) 85 : - int2d ( Th ) 86 : ( 87 : f * vh 88 : ) 89 : + on ( 1, 2, 3, 4, uh = g ); 90 : // 91 : // Solve the problem for uh. 92 : // 93 : poisson; 94 : // 95 : // Plot the computed solution. 96 : // 97 : plot ( uh, fill = 1, value = 1, ps = "poisson_square_fem.ps", 98 : cmm = "poisson_square: FEM solution" ); 99 : // 100 : // For this problem, we know the exact solution. 101 : // u: the exact solution. 102 : // 103 : Vh u = sin(5*pi*x*(1-x)) * sin(4*pi*y*(1-y)); 104 : // 105 : // Plot the exact solution. 106 : // 107 : plot ( u, fill = 1, value = 1, ps = "poisson_square_exact.ps", 108 : cmm = "poisson_square: exact" ); 109 : // 110 : // Compute the L2 error. 111 : // 112 : real l2error = sqrt ( int2d ( Th ) ( ( u - uh )^2 ) ); 113 : cout << endl; 114 : cout << " H = " << h << endl; 115 : cout << " L2 error = " << l2error << endl; 116 : // 117 : // Terminate. 118 : // 119 : cout << "\n"; 120 : cout << "poisson_square:\n"; 121 : cout << " Normal end of execution.\n"; 122 : 123 : exit ( 0 ); 124 : sizestack + 1024 =7864 ( 6840 ) poisson_square: FreeFem++ version. Solve the Poisson equation in a square region. -- Square mesh : nb vertices =324 , nb triangles = 578 , nb boundary edges 68 -- Solve : min -0.667279 max 0.949072 H = 0.0625 L2 error = 0.0101365 poisson_square: Normal end of execution. current line = 123 exit(0) err code 0 , mpirank 0 CodeAlloc : nb ptr 3747, size :479456 mpirank: 0 Ok: Normal End