Poisson's Equation with Point Source and Adaptive Mesh Refinement

This example solves a Poisson's equation with a delta-function point source on the unit disk using the adaptmesh function in the Partial Differential Equation Toolbox™.

Specifically, we solve the Poisson's equation

$$-\Delta u = \delta(x,y)$$

on the unit disk with zero Dirichlet boundary conditions. The exact solution expressed in polar coordinates is

$$u(r,\theta) = \frac{\log(r)}{2\pi},$$

which is singular at the origin.

By using adaptive mesh refinement, the Partial Equation Toolbox™ can accurately find the solution everywhere away from the origin.

Problem Definition

The following variables will define our problem:

  • g: A specification function that is used by initmesh. For more information, please see the documentation page for circleg and pdegeom.

  • c, a: The coefficients of the PDE.

  • f: A function that captures a point source at the origin. It returns 1/area for the triangle containing the origin and 0 for other triangles.

g = @circleg;
c = 1;
a = 0;
f = @circlef;

Boundary Conditions

Plot the geometry and display the edge labels for use in the boundary condition definition.

figure;
pdegplot(g, 'edgeLabels', 'on');
axis equal
title 'Geometry With Edge Labels Displayed';
% Create a pde entity for a PDE with a single dependent variable
numberOfPDE = 1;
pb = pde(numberOfPDE);
% Create a geometry entity
pg = pdeGeometryFromEdges(g);
% Solution is zero at all four outer edges of the circle
pb.BoundaryConditions = pdeBoundaryConditions(pg.Edges(1:4), 'u', 0);

Generate Mesh

adaptmesh solves elliptic PDEs using adaptive mesh generation. The 'tripick' parameter lets the user specify a function that returns which triangles will be refined in the next iteration. circlepick returns triangles with computed error estimates greater a given tolerance. The tolerance is provided to circlepick using the 'par' parameter.

[u,p,e,t] = adaptmesh(g,pb,c,a,f,'tripick','circlepick','maxt',2000,'par',1e-3);
Number of triangles: 258
Number of triangles: 515
Number of triangles: 747
Number of triangles: 1003
Number of triangles: 1243
Number of triangles: 1481
Number of triangles: 1705
Number of triangles: 1943
Number of triangles: 2155

Maximum number of triangles obtained.

Plot Finest Mesh

figure;
pdemesh(p,e,t);
axis equal

Plot Error

x = p(1,:)';
y = p(2,:)';
r = sqrt(x.^2+y.^2);
uu = -log(r)/2/pi;
figure;
pdeplot(p,e,t,'xydata',u-uu,'zdata',u-uu,'mesh','off');

Plot FEM Solution on Finest Mesh

figure;
pdeplot(p,e,t,'xydata',u,'zdata',u,'mesh','off');

Was this topic helpful?