Finite Difference Method for PDE using MATLAB (m-file)

In mathematics, finite-difference methods (FDM) are numerical methods for solving differential equations by approximating them with difference equations, in which finite differences approximate the derivatives. FDMs are thus discretization methods. FDMs convert a linear (non-linear) ODE/PDE into a system of linear (non-linear) equations, which can then be solved by matrix algebra techniques. The reduction of the differential equation to a system of algebraic equations makes the problem of finding the solution to a given ODE ideally suited to modern computers, hence the widespread use of FDMs in modern numerical analysis.



MATLAB Program:

% Finite difference method
 % Approximate the solution of y"=(-2/x)y'+(2/x^2)y+ sin(lnx)/x^2
 % for 1<=x<=2 with y(1)=1 and y(2)=2.

 p = @(x) (-2/x); 
 q = @(x) (2/x^2);
 r = @(x) (sin(log(x))/x^2);

 aa = 1; bb = 2; alpha = 1; beta = 2; n=29;      % h = (bb-aa)/(n+1);   h=0.1
 
 fprintf('   x           w   \n');
 h = (bb-aa)/(n+1);
 a = zeros(1,n+1);
 b = zeros(1,n+1);
 c = zeros(1,n+1);
 d = zeros(1,n+1);
 l = zeros(1,n+1);
 u = zeros(1,n+1);
 z = zeros(1,n+1);
 w = zeros(1,n+1);
 x = aa+h;
 a(1) = 2+h^2*q(x);
 b(1) = -1+0.5*h*p(x);
 d(1) = -h^2*r(x)+(1+0.5*h*p(x))*alpha;
 m = n-1;

 for i = 2 : m
   x = aa+i*h;
   a(i) = 2+h^2*q(x);
   b(i) = -1+0.5*h*p(x);
   c(i) = -1-0.5*h*p(x);
   d(i) = -h^2*r(x);
 end

 x = bb-h;
 a(n) = 2+h^2*q(x);
 c(n) = -1-0.5*h*p(x);
 d(n) = -h^2*r(x)+(1-0.5*h*p(x))*beta;
 l(1) = a(1);
 u(1) = b(1)/a(1);
 z(1) = d(1)/l(1);

 for i = 2 : m
   l(i) = a(i)-c(i)*u(i-1);
   u(i) = b(i)/l(i);
   z(i) = (d(i)-c(i)*z(i-1))/l(i);
 end

 l(n) = a(n)-c(n)*u(n-1);
 z(n) = (d(n)-c(n)*z(n-1))/l(n);
 w(n) = z(n);

 for j = 1 : m
   i = n-j;
   w(i) = z(i)-u(i)*w(i+1);
 end
 i = 0;
 fprintf('%5.4f    %11.8f\n', aa, alpha);
 for i = 1 : n
   x = aa+i*h;
   fprintf('%5.4f    %11.8f\n', x, w(i));
 end
 i = n+1;
 fprintf('%5.4f    %11.8f\n', bb, beta);



OUTPUTs:

>> Untitled_c
   x           w   
1.0000     1.00000000
1.0333     1.03067949
1.0667     1.06155254
1.1000     1.09262608
1.1333     1.12390423
1.1667     1.15538887
1.2000     1.18708018
1.2333     1.21897697
1.2667     1.25107701
1.3000     1.28337728
1.3333     1.31587416
1.3667     1.34856355
1.4000     1.38144105
1.4333     1.41450201
1.4667     1.44774163
1.5000     1.48115505
1.5333     1.51473732
1.5667     1.54848354
1.6000     1.58238883
1.6333     1.61644834
1.6667     1.65065735
1.7000     1.68501118
1.7333     1.71950528
1.7667     1.75413522
1.8000     1.78889666
1.8333     1.82378540
1.8667     1.85879736
1.9000     1.89392857
1.9333     1.92917520
1.9667     1.96453354
2.0000     2.00000000

>> 

3 comments:

'; (function() { var dsq = document.createElement('script'); dsq.type = 'text/javascript'; dsq.async = true; dsq.src = '//' + disqus_shortname + '.disqus.com/embed.js'; (document.getElementsByTagName('head')[0] || document.getElementsByTagName('body')[0]).appendChild(dsq); })();