% expm4.m - for Table 8 % This computes \norm T_t\norm_\infty % using a trigonometric basis close('all') clear all % Set up basic data %a=length of interval a=10 %b=diffusion paramter b=20 % number of trigonometric basis elements N=200 % set up tridiagonal matrix M=zeros(N,N); for r=1:N for s=1:N if abs(r-s)>0 M(r,s)=-(2*r*s/a)*(1-(-1)^(r+s))/(r^2-s^2); end; end; end; for r=1:N M(r,r)=-pi^2*r^2*a^(-2)/b; end; %set up alpha alp=zeros(N,1); for r=1:N alp(r,1)=(sqrt(2*a)/(pi*r))*(1-(-1)^r); end; % compute exp(Mt)*alp in the trigonometric basis format long bet=zeros(N,10); for t=1:10 bet(:,t)=expm(t*M)*alp; end; % convert to the standard basis xnum=200; vectors=zeros(xnum, N); for rx=1:xnum for n=1:N vectors(rx,n)=sqrt(2/a)*sin(pi*rx*n/(xnum)); end; end; % compute the answer values=vectors*bet; Tnorm=max(values); Tnorm'