Code covered by the BSD License  

Highlights from
Chebpack

image thumbnail

Chebpack

by

 

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

Contact us