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.

ibvp_int_newton_ex1.m
function [x,solnum]=ibvp_int_newton_ex1(n,dom,kind)
% Example: y(x)=1+sin^2(x)-int_0^x{3sin(x-t)y^2(t)}dt, x in [0,1]
% From Razzaghi, M., Ordokhani, Y., Math. Probl. Eng. 7, 205219 (2001)
% Exact solution: solex=cos(x)
% call: [x,solnum]=ibvp_int_newton_ex1(32,[0 1],2);
%
tic;x=pd(n,dom,kind);T=cpv(n,x,dom);myINIT;
dy=1;cont=1;% y initial values of the solution
[~,Vphys]=volt(n,dom,kind);
while dy(cont) > 1.e-12
    if cont>100, break, end; 
    cont=cont+1;KW=zeros(n,1);
    DKW=Vphys*diag(2*yv);KW=Vphys*(yv.^2);
    K=eye(n)+DKW;F=-yv-KW+myfun(x);cor=K\F;
    dy(cont)=norm(cor);y=cor+yv;yv=y;
end;toc;solnum=y;myOUT;
function [Vspec,Vphys]=volt(n,dom,~)
% Volterra operator with kernel K(x,t)
[~,J0]=prim(n,dom);V=T*J0/T;
[YY,XX]=meshgrid(x,x);KW=mykernel(XX,YY);
Vphys=V.*KW;Vspec=T\Vphys*T;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function myINIT
    yv=zeros(n,1);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function col=mykernel(x,t)
    % x,t matrices
    col=3*sin(x-t);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function f=myfun(x)
    f=1+sin(x).^2;% values of f(x) given by the user
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function myOUT
    solex=cos(x);figure(1);
    subplot(2,1,1);semilogy(abs(y),'.');grid;
    title('Absolute value of the coefficients of the solution');
    subplot(2,1,2);plot(x,solnum);grid;xlabel('x');ylabel('y(x)');
    title('The numerical solution');
    figure(2);
    semilogy(1:cont,dy,'.');grid;title('History of iterations - corrections');
    figure(3);plot(x,abs(solnum-solex));grid; 
    title('The error');xlabel('x');ylabel('err');
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
end


Contact us