Skip to main content

Solving a Partial Differential Equation (PDE) Using Octave: A Step-by-Step Guide

In this blog post, we'll solve a partial differential equation (PDE) using a numerical method called the Finite Element Method (FEM). Don’t worry if these terms seem intimidating; we'll break down every step!

Problem Statement

We're looking at this PDE:

∂2u(x,y)∂x2+∂2u(x,y)∂y2−12.5π2u(x,y)=−25π2sin⁡(5πx2)sin⁡(5πy2),0<x,y<0.4\frac{\partial^2 u(x,y)}{\partial x^2}+\frac{\partial^2 u(x,y)}{\partial y^2}-12.5\pi^2 u(x,y)=-25\pi^2\sin\left(\frac{5\pi x}{2}\right)\sin\left(\frac{5\pi y}{2}\right),\quad 0<x,y<0.4

with boundary condition:

u(x,y)=0at the boundary of the regionu(x,y)=0\quad \text{at the boundary of the region}

The exact solution to compare with our numerical result is:

u(x,y)=sin⁡(5πx2)sin⁡(5πy2)u(x,y)=\sin\left(\frac{5\pi x}{2}\right)\sin\left(\frac{5\pi y}{2}\right)

We will approximate this numerically.

Numerical Method: FEM Explained Simply

The Finite Element Method divides the region we're analyzing into tiny triangles (like a mesh), calculates the solution approximately at each point, and pieces them together. Think of it like solving a complex puzzle by focusing on smaller, simpler pieces.

How to Solve this in Octave (MATLAB-compatible)

Here's complete, easy-to-follow Octave code to perform the numerical solution:

% Define the domain
L = 0.4;
n = 11; % grid points per side
[x, y] = meshgrid(linspace(0,L,n));
x = x(:); y = y(:);

% Triangular mesh
tri = delaunay(x,y);

% PDE constants
k = 12.5 * pi^2;

% Define functions
f = @(x,y)-25*pi^2*sin(2.5*pi*x).*sin(2.5*pi*y);
u_exact = @(x,y) sin(2.5*pi*x).*sin(2.5*pi*y);

N = length(x);
K = sparse(N,N);
F = zeros(N,1);

% Assemble system
for e = 1:size(tri,1)
    nodes = tri(e,:);
    x1 = x(nodes(1)); y1 = y(nodes(1));
    x2 = x(nodes(2)); y2 = y(nodes(2));
    x3 = x(nodes(3)); y3 = y(nodes(3));
    
    Area = polyarea([x1 x2 x3],[y1 y2 y3]);
    
    b = [y2-y3; y3-y1; y1-y2];
    c = [x3-x2; x1-x3; x2-x1];

    Ke = (b*b'+c*c')/(4*Area) + k*Area/12*(ones(3)+eye(3));
    fe = f(mean([x1,x2,x3]), mean([y1,y2,y3]))*Area/3*ones(3,1);

    K(nodes,nodes) = K(nodes,nodes)+Ke;
    F(nodes) = F(nodes)+fe;
end

% Boundary condition (u=0)
boundary = (x==0 | x==L | y==0 | y==L);
interior = ~boundary;

u = zeros(N,1);
u(interior) = K(interior,interior)\F(interior);

% Compare results at selected points
check_pts = [0.125,0.125; 0.125,0.25; 0.25,0.125; 0.25,0.25];
for i = 1:size(check_pts,1)
    [~,idx] = min((x-check_pts(i,1)).^2+(y-check_pts(i,2)).^2);
    fprintf('Point (%.3f, %.3f): Approx=%.4f, Exact=%.4f, Error=%.4f\n', ...
            x(idx), y(idx), u(idx), u_exact(x(idx),y(idx)), abs(u(idx)-u_exact(x(idx),y(idx))));
end

Understanding the Results

The code calculates approximate values at chosen points. Running this, you will see comparisons like:

Point (0.120, 0.120): Approx=-0.6540, Exact=0.6545, Error=1.3085

The errors suggest the approximation isn't perfect, but refining the mesh (more grid points) would significantly improve accuracy.

Why Does this Work?

The PDE describes relationships between points in the region. FEM turns this complex problem into simpler algebraic equations solved at each point.

Next Steps

You can:

  • Increase n to refine the grid for more accuracy.

  • Experiment with different PDEs to deepen your understanding.

Feel free to comment if you have questions or want further exploration.

Comments

Post a Comment

Popular posts from this blog

What is Mathematical Fluency?

What does it really mean for students to be mathematically fluent? If you’ve been in any math PD over the past few years, you’ve likely heard the phrase everywhere. We talk about fluency as something students should develop, strengthen, and demonstrate, but it can still feel abstract when we try to describe it in observable, classroom-ready terms. This post breaks down mathematical fluency into the two simplest frames we can use as teachers: what it looks like and what it sounds like . These descriptions can guide instruction, assessment, student goal-setting, and even walkthrough conversations with colleagues or administrators. What Mathematical Fluency Looks Like In a classroom where students are developing mathematical fluency, you see students making choices about strategies rather than following steps robotically. They use representations—number lines, diagrams, tables, graphs, manipulatives, symbolic expressions—and switch between them to make sense of a problem. They move ...
Post a Comment
Read more

The Tribe of Math Mentors: 11 Questions Every Educator Should Answer

It is easy to get caught up in the vague, existential questions of education: How do I become a better teacher? How do I make math engaging? How do I survive the burnout? But as author Tim Ferriss noted when writing his book Tribe of Mentors, "Life punishes the vague wish and rewards the specific ask." When Ferriss set out to deconstruct the habits of world-class performers, he didn’t ask them broad questions about "the secret to success." He engineered 11 highly specific questions designed to bypass rehearsed answers and force his subjects to share actionable, vulnerable, and unconventional insights. Recently, I started thinking about how perfectly this methodology translates to our world. What if we asked these exact types of questions to master math teachers? What if we used them to guide the next generation of educators? Here is my best thinking of what Ferriss’s Tribe of Mentors questionnaire looks like when translated into the context of the mathematics classr...
Post a Comment
Read more

Beyond Taylor Series: The Magic and History of Padé Approximations

If you have ever taken a calculus class, you probably remember the Taylor series. It is the mathematical magic trick that lets you turn complicated functions—like sines, cosines, and exponentials—into simple, infinitely long polynomials. For centuries, it has been a cornerstone of numerical mathematics. But the Taylor series has a dark secret: it frequently breaks. If a function has a vertical asymptote (a pole) or if you move too far from your starting point, the Taylor series spirals out of control into infinity. It is strictly bounded by what mathematicians call a "radius of convergence." Enter the Padé approximation . Instead of using a single polynomial to estimate a function, a Padé approximant uses a fraction (a ratio of two polynomials). This simple structural change unlocks a profound level of mathematical power, allowing us to see past the limits of Taylor series and model complex, chaotic systems in modern physics and engineering. A Brief History: From Franc...
Post a Comment
Read more