bilap.edp

来自「FreeFem++可以生成高质量的有限元网格。可以用于流体力学」· EDP 代码 · 共 25 行

int n=100,nn=n+10;
real[int] xx(nn),yy(nn);

mesh Th=square(40,40);  // mesh definition of $\Omega$
fespace Vh(Th,P1);      // finite element space
macro laplacien(u,v) (dx(u)*dx(v)+dy(u)*dy(v)) // fin macro 
real f=1;
Vh u,uu,v,vv;

solve bilap([u,uu],[v,vv],solver=LU,eps=1.0e-6) =
    int2d(Th)(  laplacien(u,vv)+uu*vv //  - Delta u + uu =0 (vv)
    	         + laplacien(uu,v)  )     //  - Delta uu = 1 (v)
    	        
   - int2d(Th)(f*v)
   + on(1,2,3,4,u=0); // => v=0 also on 1,2,3,4
   
plot(u,wait=1);

for (int i=0;i<=n;i++)
 {
   xx[i]=real(i)/n;
   yy[i]=u(0.5,real(i)/n); // value of uh at point (0.5, i/10.) 
 }
 plot([xx(0:n),yy(0:n)],wait=1);