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,DF,T,DFphys]=deriv_frac(n,dom,q)
function [x,DF,T,DFphys]=deriv_frac(n,dom,q)
% DF: the matrix of the fractional Caputo derivative 
% of order q>0, dom=[0,b], in spectral form
% DFphys: T*DF/T, the same but in physical form
% call: [x,DF,T,DFphys]=deriv_frac(64,[0,1],0.5);
if (~isa(n,'double')||n<=2||floor(n)-n~=0), error('n should be a natural number n>2');end
if (~isa(dom,'double')||length(dom)~=2||dom(1)~=0||dom(2)<=0), error('dom should be a vector [0,b] with b>0');end
if (~isa(q,'double')||q<=0), error('q should be a number q>0');end
x=pd(n,dom,2);T=cpv(n,x,dom);D=deriv(n,dom);I=zeros(n);
ex=floor(q);qq=q-ex;
if qq==0, DF=D^ex;return;end
alpha=2/dom(2);
I(:,1)=x.^(1-qq)/(1-qq)/2;
I(:,2)=alpha*x.^(2-qq)/(2-qq)/(1-qq)-x.^(1-qq)/(1-qq);
I(:,3)=4*alpha^2/(3-qq)/(2-qq)/(1-qq)*x.^(3-qq)-4*alpha/(2-qq)/(1-qq)*x.^(2-qq)+x.^(1-qq)/(1-qq);
for k=4:n
    I(:,k)=2*(alpha*x-1).*I(:,k-1)+((1-qq)/(k-3)-1).*I(:,k-2)-2*(-1)^(k-1)/(k-1)/(k-3).*x.^(1-qq);
    I(:,k)=I(:,k)/(1+(1-qq)/(k-1));
end
DF=T\I*D^(ex+1)/gamma(1-qq);
DFphys=I*D^(ex+1)/T/gamma(1-qq);
end

Contact us