convect-3d.edp

FreeFem++可以生成高质量的有限元网格。可以用于流体力学

load "msh3"

int nn=8;

mesh Th2=square(nn,nn,[x*2.-1.,y*2.-1.]);
fespace Vh2(Th2,P1);
int[int] rup=[0,2],  rdown=[0,1], rmid=[1,1,2,1,3,1,4,1];
real zmin=-1,zmax=1.;

mesh3 Th=buildlayers(Th2,nn,
  zbound=[zmin,zmax],
  // reftet=r1, 
  reffacemid=rmid, 
  reffaceup = rup,
  reffacelow = rdown);
func  real hill(real r2){return exp(-10.*(r2));};

fespace Vh(Th,P13d);

macro Grad(u) [dx(u),dy(u),dz(u)]// EOM
macro div(u1,u2,u3) (dx(u1)+dy(u2)+dz(u3)) //EOM

Vh v,vo;
Vh2 v2;
real x0=0.3,y0=0.3,z0=0;
vo=hill(square(x-x0)+square(y-y0)+square(z-z0));
real t=0;
v2=vo(x,y,0);
plot(v2,cmm=" cut y = 0.5, time ="+t,wait=1);
real dt=0.1;
func u1=1.;
func u2=2.;
func u3=3.;
v=convect([u1,u2,u3],-dt,vo);
v2=v(x,y,0);
t += dt;
plot(v2,cmm=" cut y = 0.5, time ="+t,wait=1);