function U = finedif(f,g,a,b,c,n,m)
% U = finedif(f,g,a,b,c,n,m)
% Finite difference solution to the wave equation.
% f is a function, input.
% g is a function, input.
% a is the width of [0 a], input.
% b is the width of [0 b], input.
% c is the constant in the wave equation, input.
% n is the number of grid points over [0 a], input.
% m is the number of grid points over [0 b], input.
% U is the solution matrix, output.
h = a/(n-1);
k = b/(m-1);
r = c*k/h;
r2 = r^2;
r22 = r^2/2;
s1 = 1 - r^2;
s2 = 2 - 2*r^2;
U = zeros(n,m);
for i=2:(n-1),
U(i,1) = feval(f,h*(i-1));
U(i,2) = s1*feval(f,h*(i-1)) + k*feval(g,h*(i-1)) ...
+ r22*(feval(f,h*(i)) + feval(f,h*(i-2)));
end
for j=3:m,
for i=2:(n-1),
U(i,j) = s2*U(i,j-1) + r2*(U(i-1,j-1) + U(i+1,j-1)) - U(i,j-2);
end
end