-- FreeFem++ v4.14 (mer. 06 mars 2024 16:59:04 CET - git v4.14-1-g2b2052ae)
   file : driven_cavity_stokes.edp
 Load: lg_fem lg_mesh lg_mesh3 eigenvalue 
    1 : //  driven_cavity_stokes.edp
    2 : //
    3 : //  Discussion:
    4 : //
    5 : //    Steady Stokes problem in a cavity.
    6 : //
    7 : //      - uxx - uyy + px = uf
    8 : //      - vxx - vyy + py = vf
    9 : //        ux  + vy       = pf
   10 : //
   11 : //      u = 1 along top wall.
   12 : //
   13 : //  Licensing:
   14 : //
   15 : //    This code is distributed under the MIT license.
   16 : //
   17 : //  Modified:
   18 : //
   19 : //    22 May 2020
   20 : //
   21 : //  Reference:
   22 : //
   23 : //    John Chrispell, Jason Howell,
   24 : //    Finite Element Approximation of Partial Differential Equations
   25 : //    Using FreeFem++, or, How I Learned to Stop Worrying and Love
   26 : //    Numerical Analysis.
   27 : //
   28 : //    Frederic Hecht,
   29 : //    New development in FreeFem++,
   30 : //    Journal of Numerical Mathematics,
   31 : //    Volume 20, Number 3-4, 2012, pages 251-265.
   32 : //
   33 : cout << "\n";
   34 : cout << "driven_cavity_stokes():\n";
   35 : cout << "  FreeFem++ version.\n";
   36 : //
   37 : //  eps: the penalty coefficient for the continuity equation.
   38 : //
   39 : real eps = 1.0E-10;
   40 : //
   41 : //  n: the number of nodes in each spatial direction.
   42 : //
   43 : int n = 30;
   44 : //
   45 : //  Th: the triangulation of the square.
   46 : //
   47 : mesh Th = square ( n, n );
   48 : //
   49 : //  Vh: the finite element space for 2D velocity vector (piecewise linear "bubble") 
   50 : //
   51 : fespace Vh ( Th, P1b );
   52 : // 
   53 : //  Qh: the finite element space for the scalar pressure (piecewise linear).
   54 : //
   55 : fespace Qh ( Th, P1 );
   56 : //
   57 : //  u1, u2: velocity trial functions.
   58 : //
   59 : Vh u1, u2;
   60 : //
   61 : //  p: pressure trial function.
   62 : //
   63 : Qh p;
   64 : //
   65 : //  v1, v2: velocity test functions.
   66 : //
   67 : Vh v1, v2;
   68 : //
   69 : //  q: pressure test function.
   70 : //
   71 : Qh q;
   72 : //
   73 : //  Define and solve the problem.
   74 : //
   75 : solve Stokes ( [u1,u2,p], [v1,v2,q] ) = int2d ( Th, qforder = 3 ) 
   76 :   ( 
   77 :     dx(u1) * dx(v1) 
   78 :   + dy(u1) * dy(v1)
   79 :   + dx(p)  *    v1
   80 :   + dx(u2) * dx(v2) 
   81 :   + dy(u2) * dy(v2)
   82 :   + dy(p)  *    v2
   83 :   + dx(u1) *    q
   84 :   + dy(u2) *    q
   85 :   - eps * p   * q
   86 :   ) 
   87 :   + on ( 1, 2,    4, u1 = 0 )
   88 :   + on (       3,    u1 = 1 )
   89 :   + on ( 1, 2, 3, 4, u2 = 0 );
   90 : //
   91 : //  Plot the velocity vectors, using 1 color for the vectors.
   92 : //
   93 : plot ( [ u1, u2 ], nbarrow = 1, ps = "driven_cavity_stokes_velocity.ps" );
   94 : //
   95 : //  Plot the pressure.
   96 : //
   97 : plot ( p, ps = "driven_cavity_stokes_pressure.ps" );
   98 : //
   99 : //  Terminate.
  100 : //
  101 : cout << "\n";
  102 : cout << "driven_cavity_stokes():\n";
  103 : cout << "  Normal end of execution.\n";
  104 : 
  105 : exit ( 0 );
  106 : 
  107 :  sizestack + 1024 =1904  ( 880 )


driven_cavity_stokes():
  FreeFem++ version.
  -- Square mesh : nb vertices  =961 ,  nb triangles = 1800 ,  nb boundary edges 120
  -- Solve : 
          min -0.190866  max 1
          min -0.326409  max 0.313182
          min -487.151  max 399.511

driven_cavity_stokes():
  Normal end of execution.
  current line = 105
exit(0)
 err code 0 ,  mpirank 0
 CodeAlloc : nb ptr  3932,  size :524080 mpirank: 0
Ok: Normal End