We solve the following PDE in polar coordinates \(u=u(r,q)\):
\begin{align} (r\,u_{r})_{r} +\frac{1}{4r}\,(u_{q})_{q}=0,\quad 1<r<3,\quad 0<q<\pi \end{align}
with the boundary conditions \(u(1,q)=\sin(2\,q)\), \(u(3,q)=0\), and \(u(r,0)=u(r,\pi)=0\).
Thanks to Rene at reneif07@gmail.com
![]() | ![]() |
download example: 2DLPSemiAnnulus.edp or return to 2D examples
//We solve the following PDE in polar coordinates // (r*u_r)_r +(0.25/r)*(u_q)_q=0, for 1<r<3, 0<q<pi //Boundary conditions u(1,q)=sin(2*q), u(3,q)=0, u(r,0)=u(r,pi)=0 load "Element_P3"; real beta=pi/2, r1=1, r2=3,a=1./4.; border L1(t=-beta,beta){ x=r1*sin(t); y=r1*cos(t);label=4;} border L2(t=-beta,beta){ x=-r2*sin(t); y=r2*cos(t);label=3;} border L3(t=0,1){x=-r2*sin(beta)+(-r1*sin(beta)+r2*sin(beta))*t; y=r2*cos(beta)+(r1*cos(beta)-r2*cos(beta))*t;label=2;} border L4(t=0,1){x=r1*sin(beta)+(r2*sin(beta)-r1*sin(beta))*t; y=r1*cos(beta)+(r2*cos(beta)-r1*cos(beta))*t;label=1;} for (int i=1; i<=40; i++) { int n=5*i; mesh T0h = buildmesh(L1(n)+L2(n)+L3(n)+L4(n)); //plot(T0h,wait=1); fespace V0h(T0h,P3); V0h u,v,r,Ur,Uq,q; u=0;v=0;r=sqrt(x^2+y^2);q=atan2(y,x); Uq=sin(2.*q); Ur=((9-(r)^(2.))/(8.*r)); // macro Grad(u) [(dx(u)*cos(q)+dy(u)*sin(q))/sqrt(r), (dy(u)*cos(q)-dx(u)*sin(q))*sqrt(a)/sqrt(r)] // EOM macro Lap(u,v) (Grad(u)'*Grad(v)) // macro Jac()(r) //Jacobian // problem static(u,v)=int2d(T0h)(Lap(u,v)*Jac) + on(1,2,3,u=0) + on(4,u=Uq); //exact solution V0h f=Uq*Ur; // static; V0h err=u-f; V0h Int=(err)^2.*Jac; real l2e=sqrt(int2d(T0h)(Int)); cout << " L2-error " << " = "<<l2e<< endl; plot(u,wait=1,value=1,fill=1,cmm="Max(u_numeric)="+u[].max+", Min(u_numeric)="+u[].min+", L^2-error="+l2e); plot(f,value=1,wait=1,fill=1,cmm="Max(u_exact)="+f[].max+", Min(u_exact)="+f[].min+", L^2-error="+l2e); }