with initial and boundary conditions given by
To do this we used a finite difference approximation of the spatial derivatives on the right hand side of the equation to obtain a system of ordinary differential equations
where L(t) is a tridiagonal matrix and Y(t) and F(t) are vectors.
The code pdeuler.c solves this system using Euler's method
and the code pdetrap.c solves this system using the trapezoid method
Note that these numeric schemes do not yield convergent approximations for all possible choices of v(t,x), p(t,x), q(t,x) and f(t,x). In particular, not all partial differential equations of this form even have unique differentiable solutions. Moreover, for numerical stability the ratio r=k/h2 must also be small.
To test our programs we choose
to obtain the heat equation.
We used a Makefile to compile the programs and create the output files pdeuler.out and pdetrap.out. The gnuplot commands
set style data lines plot "pdeuler.out" index 0display the initial condition
The gnuplot commands
plot "pdeuler.out" index 100, "pdetrap.out" index 100create a plot of y(T,x) versus x which compares the approximation obtained by Euler's method to the approximation obtained by the trapezoid method.
The graphs of the two approximations overlay each other and represent a smoothed version of the initial condition. Since it is known that diffusion of heat physically smoothes the temperature distribution, this is a partial check of the two methods.