Code covered by the BSD License

Chebpack

Damian Trif (view profile)

15 Jul 2011 (Updated )

The MATLAB package Chebpack solves specific problems for differential or integral equations.

[x,solnum]=ibvp_sys_test(n,dom,kind,y0)
```function [x,solnum]=ibvp_sys_test(n,dom,kind,y0)
% Example: y1'=y2;y2'+4x^2 y1=2cos(x^2), x in [0,40]
%          y1(0)=y2(0)=0
% From: M. K. El Daou, Proceedings of the World Congress on Engineering
%       2008, Vol II, London, U.K.
% Exact solutions: y1ex=sin(x^2), y2ex=2xcos(x^2)
% call: y0=[0,0];for k=1:20,[x,solnum]=ibvp_sys_test(128,[2*k-2 2*k],2,y0);
%       y0=solnum(end,:);end;
tic;m=length(y0);A=spalloc(m*n,m*n,5*m*n);b=zeros(m*n,1);
x=pd(n,dom,kind);X=mult(n,dom);D=deriv(n,dom);
myDE;myBC;sol=reshape(A\b,n,m);% the solution in spectral space
solnum = t2x(sol,kind);% the solution in physical space
toc;myOUT;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function myDE
% describes the differential system
A=[D -speye(n);4*X^2 D];
r=[zeros(n,1),2*cos(x.^2)];b=reshape(x2t(r,kind),m*n,1);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function myBC
% implements the initial conditions
T=cpv(n,dom(1),dom);
for k=1:m, A(k*n,:)=0;A(k*n,(k-1)*n+1:k*n)=T;b(k*n)=y0(k);end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function myOUT
xx=linspace(dom(1),dom(2),10000);
f1 = barycheb(xx,solnum(:,1),x,kind);f2 = barycheb(xx,solnum(:,2),x,kind);
y1ex=sin(x.^2);y2ex=2*x.*cos(x.^2);err=norm(solnum-[y1ex,y2ex]);
subplot(2,1,1);semilogy(x,abs(sol),'.');
title('Absolute value of the Chebyshev coefficients of y_1 and y_2');
grid;hold on;
subplot(2,1,2);plot(xx,f1,'b',xx,f2,'r');legend('y_1(x)','y_2(x)');
title(['Numerical solution y_{1,2}(x)  (error = ',num2str(err),')']);
xlabel('x');ylabel('y_{1,2}(x)');grid;hold on;pause(0.25)
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
end```