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.

test1
function test1
% -x^2u''-x/6u'+u=19/6x^(5/2),u(0)=u(1)=0
% exact solution: x^(3/2)-x^(5/2)
%

display('Chebpack');
n=80;dom=[0 1];kind=1;tic;x = pd(n,dom,kind);X=mult(n,dom);D=deriv(n,dom);
A=-X^2*D^2-X*D/6+speye(n);b = x2t((x.^(5/2))*19/6,kind);T=cpv(n,dom,dom);
A(n-1,:)=T(1,:);b(n-1)=0;A(n,:)=T(2,:);b(n)=0;
c=A\b;sol = t2x(c,kind);solex=x.^(3/2)-x.^(5/2);toc;err=norm(sol-solex)
% figure(1);semilogy(x,abs(sol-solex));
% figure(2);plot(x,sol);

display('bvp4c');
solinit=bvpinit(linspace(0,1,5),[0;0]);S=[0,1;1,5/6];
options=bvpset('SingularTerm',S,'RelTol',1.e-5);
tic;sol=bvp4c(@ex1_bvp4code,@ex1_bvp4cbc,solinit,options);toc;
x=sol.x;solex=x.^(3/2)-x.^(5/2);err=norm(solex-sol.y(1,:))

display('bvp5c');
solinit=bvpinit(linspace(0,1,5),[0;0]);S=[0,1;1,5/6];
options=bvpset('SingularTerm',S,'RelTol',1.e-4);
tic;sol=bvp5c(@ex1_bvp4code,@ex1_bvp4cbc,solinit,options);toc;
x=sol.x;solex=x.^(3/2)-x.^(5/2);err=norm(solex-sol.y(1,:))

display('dms');
n=64;tic;g=[1 0 0;1 0 0];[t,D2t,D1t,phip,phim]=cheb2bc(n,g);
f=19/6*((t+1)/2).^(5/2);
p=-(t+1).^2.*phip(:,2)-(t+1)/6.*phip(:,1);
m=-(t+1).^2.*phim(:,2)-(t+1)/6.*phim(:,1);
D=-diag((t+1).^2)*D2t-diag(t+1)/6*D1t+eye(size(D1t));
sol=D\(f-p-m);solex=((t+1)/2).^(3/2)-((t+1)/2).^(5/2);toc;
err=norm(sol-solex)

display('Chebfun');
tic;[d,x]=domain(0,1);I=eye(d);D=diff(d);
A=-diag(x.^2)*D^2-diag(x)/6*D+I;A.bc='dirichlet';
f=19/6*x.^(5/2);
sol=A\f;solex=x.^(3/2)-x.^(5/2);toc;
err=max(abs(sol-solex))

function dydx=ex1_bvp4code(x,y)
dydx=[0
        -19/6*x^(3/2)];
end
function res=ex1_bvp4cbc(ya,yb)
    res=[ya(1)
            yb(1)];
end
end

Contact us