-- FreeFem++ v  3.320001 (date Ven  7 nov 2014 14:28:55 CET)
 Load: lg_fem lg_mesh lg_mesh3 eigenvalue 
    1 : //  Discussion:
    2 : //
    3 : //  Boundary layer
    4 : //  -eps * Laplace(u) + 2 du/dx + du/dy = f on [-1,+1]x[-1,+1]
    5 : //  u = u0 on the boundary.
    6 : //  u0 = (1-exp(-(1-x)/eps)*(1-exp(-(1-y)/eps)*cos(pi(x+y))
    7 : //  f = -eps*(u0'')+2 du0/dx+du0/dy
    8 : //
    9 : //  Discussion:
   10 : //
   11 : //    Suggested parameter values:
   12 : //    * eps = 0.1
   13 : //    * eps = 0.001
   14 : //
   15 : //  Location:
   16 : //
   17 : //    http://people.sc.fsu.edu/~jburkardt/examples/mitchell_freefem++/test06a.edp
   18 : //
   19 : //  Licensing:
   20 : //
   21 : //    This code is distributed under the GNU LGPL license.
   22 : //
   23 : //  Modified:
   24 : //
   25 : //    18 December 2014
   26 : //
   27 : //  Author:
   28 : //
   29 : //    John Burkardt
   30 : //
   31 : //  Reference:
   32 : //
   33 : //    Frederic Hecht,
   34 : //    Freefem++,
   35 : //    Third Edition, version 3.22
   36 : //
   37 : //    William Mitchell,
   38 : //    A collection of 2D elliptic problems for testing adaptive
   39 : //    grid refinement algorithms,
   40 : //    Applied Mathematics and Computation,
   41 : //    Volume 220, 1 September 2013, pages 350-364.
   42 : //
   43 : border bottom ( t = -1.0, +1.0 )   { x = t;    y = -1.0; label = 1; }
   44 : border right  ( t = -1.0, +1.0 )   { x = +1.0; y = t;    label = 1; }
   45 : border top    ( t = +1.0, -1.0 )   { x = t;    y = +1.0; label = 1; }
   46 : border left   ( t = +1.0, -1.0 )   { x = -1.0; y = t;    label = 1; }
   47 : //
   48 : //  Define Th, the triangulation of the "left" side of the boundaries.
   49 : //
   50 : int n = 40;
   51 : mesh Th = buildmesh ( bottom ( n ) + right ( n ) + top ( n ) + left ( n ) );
   52 : //
   53 : //  Define Vh, the finite element space defined over Th, using P1 basis functions.
   54 : //
   55 : fespace Vh ( Th, P1 );
   56 : //
   57 : //  Define U, V, and F, piecewise continuous functions over Th.
   58 : //
   59 : Vh u;
   60 : Vh v;
   61 : //
   62 : //  Define the problem parameter.
   63 : //
   64 : real eps = 0.1;
   65 : //
   66 : //  Define the right hand side function F.
   67 : //
   68 : func f = exp ( - 2.0 / eps ) / pow ( eps, 2 ) * 
   69 :   (
   70 :     - exp ( ( 1.0 - y ) / eps ) * 
   71 :     ( 
   72 :       -2.0 * exp ( 1.0 / eps ) * pow ( pi * eps, 2 )
   73 :      + exp ( x / eps ) * ( 1.0 + 2.0 * pow ( pi * eps, 2 ) ) 
   74 :     ) 
   75 :     * cos ( pi * ( x + y ) ) - 
   76 :     ( 
   77 :       3.0 * exp ( 2.0 / eps ) - exp ( ( 1.0 - x ) / eps )
   78 :       - exp ( ( 1.0 - y ) / eps ) - exp ( ( x + y ) / eps ) 
   79 :     )
   80 :     * eps * pi * sin ( pi * ( x + y ) ) 
   81 :   );
   82 : //
   83 : //  Define the Dirichlet boundary condition function G.
   84 : //
   85 : func g = ( 1.0 - exp ( - ( 1.0 - x ) / eps ) ) * 
   86 :          ( 1.0 - exp ( - ( 1.0 - y ) / eps ) ) * 
   87 :          cos ( pi * ( x + y ) );
   88 : //
   89 : //  Solve the variational problem.
   90 : //
   91 : solve Laplace ( u, v ) 
   92 :   = int2d ( Th ) ( eps*dx(u)*dx(v) + eps*dy(u)*dy(v) + 2*dx(u)*v + dy(u)*v )
   93 :   - int2d ( Th ) ( f * v ) 
   94 :   + on ( 1, u = g );
   95 : //
   96 : //  Plot the solution.
   97 : //
   98 : plot ( u, wait = 1, fill = true, ps = "test06a_u.eps" );
   99 : //
  100 : //  Plot the mesh.
  101 : //
  102 : plot ( Th, wait = 1, ps = "test06a_mesh.eps" );
  103 : //
  104 : //  Save the mesh file.
  105 : //
  106 : savemesh ( Th, "test06a.msh" );
  107 : 
  108 :  sizestack + 1024 =6936  ( 5912 )

  --  mesh:  Nb of Triangles =   3782, Nb of Vertices 1972
  -- Solve : 
          min -12052.4  max 5729.13
  number of required edges : 0
times: compile 0s, execution 0.14s,  mpirank:0
 CodeAlloc : nb ptr  2804,  size :347176 mpirank: 0
Ok: Normal End