/////////////////////////////////////////////////////////////
//                                                         //
//  Solve poisson equation with Omega = circle of radius 1 //
//                                                         //
//  | - ( u_xx + u_yy ) = f  in Omega                      //
// <                                                       //
//  |                 u = g  on dOmega                     //
//                                                         //
/////////////////////////////////////////////////////////////

verbosity=0; // disable freefem output

// Parametrization of the circle
border dOmega(t=0,2*pi){x=cos(t); y=sin(t);}

// Mesh generation
int nodes = 100;
mesh Omega = buildmesh(dOmega(nodes));

// Finite element space
fespace X(Omega,P1);

// Finite element variables
X u,v;

// Functions
func f=((x-.1)^2+(y-.2)^2<.25) ? 10 : 0;  // source   f=10*(characteristic function of a circle)
func g=sin(4*atan2(y,x));                 // dirichlet g=sin(4*theta) with theta=arctan(y/x)

// Weak formulation of the problem
problem poisson(u,v) =  int2d(Omega)( dx(u)*dx(v) + dy(u)*dy(v) ) //  bilinear form
                      - int2d(Omega)( f*v )                       //  linear form
                      + on( dOmega, u=g );                        //  dirichlet boundary condition 
      
poisson;  // solve the problem 

plot(u, fill=true, value=true, wait=true, dim=2); // plot solution 2d with filled level sets and legend, waiting for a key press
plot(u, fill=true, value=true, wait=true, dim=3); // plot solution 3d with filled level sets and legend, waiting for a key press