-- FreeFem++ v4.6 (Thu Apr  2 15:47:38 CEST 2020 - git v4.6)
 Load: lg_fem lg_mesh lg_mesh3 eigenvalue 
    1 : // membrane_error.edp
    2 : //
    3 : //  Discussion:
    4 : //
    5 : //    Solve the membrane problem using a manufactured solution.
    6 : //
    7 : //    We use a circular domain C with boundary Gamma1 + Gamma2.
    8 : //
    9 : //    We have 
   10 : //      u = sin ( x^2 + y^2 - 1 )
   11 : //    for which we need to specify a nonhomogeneous Neumann condition on Gamma2.
   12 : //
   13 : //  Location:
   14 : //
   15 : //    http://people.sc.fsu.edu/~jburkardt/freefem++/membrane_error/membrane_error.edp
   16 : //
   17 : //  Licensing:
   18 : //
   19 : //    This code is distributed under the GNU LGPL license.
   20 : //
   21 : //  Modified:
   22 : //
   23 : //    28 July 2015
   24 : //
   25 : //  Reference:
   26 : //
   27 : //    Frederic Hecht,
   28 : //    Freefem++,
   29 : //    Third Edition, version 3.22
   30 : //
   31 : cout << "\n";
   32 : cout << "MEMBRANE_ERROR:\n";
   33 : cout << "  FreeFem++ version.\n";
   34 : //
   35 : //  Define Gamma1 and Gamma2, the boundaries.
   36 : //  A and B are the lengths of the semimajor and semiminor axes:
   37 : //  THETA specifies the angle at which we switch from GAMMA1 to GAMMA2.
   38 : //
   39 : real a = 1.0;
   40 : real b = 1.0;
   41 : real theta = 4.0 * pi / 3.0;
   42 : border Gamma1 ( t = 0, theta )    { x = a * cos(t); y = b*sin(t); }
   43 : border Gamma2 ( t = theta, 2 * pi ) { x = a * cos(t); y = b*sin(t); }
   44 : //
   45 : // Define the exact solution.
   46 : //
   47 : func uexact = sin ( x^2 + y^2 - 1.0 );
   48 : //
   49 : //  Make space for the errors.
   50 : //
   51 : real[int] L2error(2);
   52 : 
   53 : int n;
   54 : 
   55 : for ( n = 0; n < 2; n++ )
   56 : {
   57 : //
   58 : //  Define Th, the triangulation of the "left" side of the boundaries.
   59 : //
   60 :   mesh Th = buildmesh ( Gamma1(20*(n+1)) + Gamma2(10*(n+1)) );
   61 : //
   62 : //  Define Vh, the finite element space defined over Th, using P2 basis functions.
   63 : //
   64 :   fespace Vh ( Th, P2 );
   65 : //
   66 : //  Define U and V, piecewise P2 continuous functions over Th.
   67 : //
   68 :   Vh u;
   69 :   Vh v;
   70 : //
   71 : //  Define F, the right hand side.
   72 : //
   73 :   func f = - 4.0 * ( cos ( x^2 + y^2 - 1.0 ) 
   74 :     - ( x^2 + y^2 ) * sin ( x^2 + y^2 - 1.0 ) );
   75 : //
   76 : //  Solve the problem.
   77 : //
   78 :   solve Laplace ( u, v ) 
   79 :     = int2d ( Th ) ( dx(u)*dx(v) + dy(u)*dy(v) )
   80 :     - int2d ( Th ) ( f * v ) 
   81 :     - int1d ( Th, Gamma2 ) ( 2.0 * v )
   82 :     + on ( Gamma1, u = 0.0 );
   83 : //
   84 : //  Plot the solution.
   85 : //
   86 :   plot ( u, fill = true, wait = true, ps = "membrane_error_u_"+n+".ps" );
   87 : //
   88 : //  Plot the mesh.
   89 : //
   90 :   plot ( Th, wait = true, ps = "membrane_error_mesh_"+n+".ps" );
   91 : //
   92 : //  Compute the error.
   93 : //
   94 :   L2error[n] = sqrt ( int2d ( Th ) ( ( u - uexact )^2 ) );
   95 :   cout << "  L2error " << n << " = " << L2error[n] << endl;
   96 : }
   97 : 
   98 : cout << "  Convergence rate = " 
   99 :      << log ( L2error[0] / L2error[1] ) / log ( 2.0 ) << endl;
  100 : //
  101 : //  Terminate.
  102 : //
  103 : cout << "\n";
  104 : cout << "MEMBRANE_ERROR:\n";
  105 : cout << "  FreeFem++ version.\n";
  106 : 
  107 :  sizestack + 1024 =3272  ( 2248 )


MEMBRANE_ERROR:
  FreeFem++ version.
  --  mesh:  Nb of Triangles =    162, Nb of Vertices 97
  -- Solve : 
          min -0.831386  max 0.0182362
  L2error 0 = 0.0183715
  --  mesh:  Nb of Triangles =    640, Nb of Vertices 351
  -- Solve : 
          min -0.83894  max 0.00464906
  L2error 1 = 0.0046527
  Convergence rate = 1.98133

MEMBRANE_ERROR:
  FreeFem++ version.
times: compile 0.004613s, execution 1.78816s,  mpirank:0
 CodeAlloc : nb ptr  3706,  size :479840 mpirank: 0
Ok: Normal End