0

尽管我阅读了算法，但我对非线性有限差分法的编码有问题。这是《数值分析》 （Burden RL, Faires JD, 2010）一书中的一个例子，第 11.4 节，第 694 页。

下面是我从书中得到的算法。我需要为此编写代码：

非线性有限差分

我成功地为线性问题编写了相同的代码（上面提到的书中的第 11.3 节）：

function LinearFiniteDifference() 
  clear
  clc
 
 p = @(x) (-2/x);  
 q = @(x) (2/x^2);
 r = @(x) (sin(log(x))/x^2);
 
 ga = 1; gb = 2; alpha = 1; beta = 2; N = 9;      
 % h = (gb - ga)/(N + 1);   h = 0.1
  
 fprintf('   x           w   \n');
 h = (gb - ga)/(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 = ga + 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;
 
 for i = 2 : N - 1 
   x = ga + 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);
 endfor
 
 x = gb - 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 : N - 1 
   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);
 endfor
 
 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 i = N - 1 : -1 : 1 
   w(i) = z(i) - u(i)*w(i+1);
 endfor
 i = 0;
 fprintf('%5.4f    %11.8f\n', ga, alpha);
 for i = 1 : N 
   x = ga + i*h;
   fprintf('%5.4f    %11.8f\n', x, w(i));
 endfor
 i = N + 1;
 fprintf('%5.4f    %11.8f\n', gb, beta);
  
endfunction

我认为非线性也是一样的，但事实并非如此。请给我一些帮助，我真的很感激。太感谢了！

4

0 回答 0