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_ode_ex4(n,dom,kind)
function [x,solnum]=ibvp_ode_ex4(n,dom,kind)
% Example: y''-3/x y'+(4x^2+3/x^2)y=8 x^3, x in [1,25]
%          y(1)=2+sin(1), y'(1)=2+2*cos(1)+sin(1) 
% From: M. K. El Daou, International Journal of Computational and 
%       Mathematical Sciences,3:4 2009, pp. 153-159
% Exact solution: solex=x sin(x^2)+2x
% call: [x,solnum]=ibvp_ode_ex4(512,[1 25],2);
% 
tic;
A=zeros(n);b=zeros(n,1);
x=pd(n,dom,kind);X=mult(n,dom);D=deriv(n,dom);
myDE;myBC;sol=A\b;% the solution in spectral form
toc;
solnum = t2x(sol,kind);% the solution in physical form
myOUT
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function myDE
    f=8*x.^5;
    A=X^2*D^2-3*X*D+(4*X^4+3*speye(n));b=x2t(f,kind);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function myBC
    T=cpv(n,dom(1),dom);% T(1,:) for xc=1
    A(n-1,:)=T(1,:);b(n-1)=2+sin(1);
    A(n,:)=T(1,:)*D;b(n)=2+2*cos(1)+sin(1);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function myOUT
    xx=linspace(1,25,10000);
    fxnum = barycheb(xx,solnum,x,kind);
    solexx=xx.*sin(xx.^2)+2*xx;solex=x.*sin(x.^2)+2*x;
    subplot(2,1,1);semilogy(abs(sol),'.');grid;
    title('The absolute value of the coefficients of the solution');
    subplot(2,1,2);plot(xx,fxnum,'r',xx,solexx,'b');grid;
    title(['The solution (err = ',num2str(norm(solnum-solex)),')']);
    legend('numeric','exact');xlabel('x');ylabel('y(x)');
    display(solnum(end)-solex(end));% the error at x=25
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
end

Contact us