Crank-Nicholson MATLAB

function [t,x,U] = heat1 (T,a,m,n,beta,c,f,g)
%----------------------------------------------------------------------
% Usage: [t,x,U] = heat1 (T,a,m,n,beta,c,f,g)
%
% Description: Use the Crank-Nicolson method to solve the following
% partial differential equation called the one-dimensional
% heat equation (also the diffusion equation):
%
% (d/dt)u(t,x) = beta*(d/dx)*(d/dx)u(t,x)
%
% The heat equation is to be solved over the region:
%
% 0 <= t <= T,
% 0 < x < a
%
% subject to the following initial and boundary conditions:
%
% u(0,x) = f(x)
% (1 - c(1))*u(t,0) - c(1)*(d/dt)u(t,0) = g(t,1)
% (1 - c(2))*u(t,a) + c(2)*(d/dt)u(t,a) = g(t,2)
%
% where c(1) and c(2) are in the ran