FreeFEM++ 例題9-6 2次元ストークス流
uは流速、pは圧力である。まず、キャビティ流れを扱おう。この問題はボックスの三方を固定壁で囲まれ、上側の側面がu=1で移動することで、領域内の流体が対流を起こすという問題である。境界条件として、上側で流速u=1、その他の壁面ではすべりなし条件(no slip条件:u,v=0)を与える。
コード名: cavity.edp (u1,u2,pは流速と圧力、v1,v2,qは重み関数)
mesh Th=square(8,8);    //四角に切る
fespace Xh(Th,P2);
fespace Mh(Th,P1);
Xh u2,v2;
Xh u1,v1;
Mh p,q;
solve Stokes ([u1,u2,p],[v1,v2,q],solver=Crout) =
    int2d(Th)( ( dx(u1)*dx(v1) + dy(u1)*dy(v1)   +  dx(u2)*dx(v2) + dy(u2)*dy(v2) )
            - p*q*(0.000001)
            - p*dx(v1)- p*dy(v2)  - dx(u1)*q- dy(u2)*q )
  + on(3,u1=1,u2=0)  + on(1,2,4,u1=0,u2=0);    //境界条件
 
label=1
label=2
label=3
label=4
(u=1,v=0)
ペナルティ法: 非圧縮の条件に微少量εと圧力pを掛けたものを加える。
ここでは流速の重み関数をv、圧力の重み関数をqと表し