Partial Differential Equation Toolbox

Poisson's Equation on Rectangular Domain Using a Fast Poisson Solver

This example solves a Poisson's equation using the poisolv function in the Partial Equation Toolbox™.

The specific Poisson's problem is

$$-\Delta u = 3x^2 $$

on a square domain with Dirichlet boundary conditions,

A Fast Poisson Solver, like the one found in poisolv, is able to numerically solve a Poisson's Equation on a regular grid is less time than the more general solver found in assempde.

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 squarereg and pdegeom.

  • b: A boundary file used by assempde. For more information on this file, please see the documentation pages for squareb4 and pdebound.

  • c, a, f: The coefficients of the PDE.

g = @squareg;
c = 1;
a = 0;
f = '3*x.^2';

Boundary Conditions

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

figure;
pdegplot(g, 'edgeLabels', 'on');
axis([-1.1 1.1 -1.1 1.1]);
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);
innerBCFunc = @(thePde, loc, state) 0.2*cos(pi/2*loc.y);
b1 = pdeBoundaryConditions(pg.Edges(2), 'u', innerBCFunc);
b2 = pdeBoundaryConditions(pg.Edges([1 3 4]), 'u', 0);
pb.BoundaryConditions = [b1 b2];

Mesh for Fast Poisson Solver

The fast Poisson solver requires a regular rectangular grid. A mesh meeting this requirement that has been generated by poimesh is shown below.

[p,e,t] = poimesh(g,16);
figure;
pdemesh(p,e,t);
axis equal

Solve Using Both Fast Poisson Solver and Standard Solver

We solve using both the fast Poisson solver that is implemented in poisolv and the standard solver that is implemented in assempde.

fprintf('%-5s|%15s|%15s\n','n','Fast Poisson','Standard');
fprintf('-----|---------------|---------------\n');
for n = [16 32 64 128]
    [p,e,t] = poimesh(g,n);
    tic;
    u = poisolv(pb,p,e,t,f);
    tfast = toc;
    tic;
    u1 = assempde(pb,p,e,t,c,a,f);
    tstandard = toc;
    fprintf('%-5d|%15.5g|%15.5g\n',n,tfast,tstandard);
end
n    |   Fast Poisson|       Standard
-----|---------------|---------------
16   |         0.1732|       0.035352
32   |       0.048946|       0.035086
64   |        0.10908|       0.090939
128  |        0.23142|        0.48139

Plot FEM Solution on Finest Mesh

figure;
pdesurf(p,t,u);