函数文件1:
1 function b=F(f,x0,h,N) 2 % b(1,1)=x0(1)-h*x0(2)-u(1); 3 % b(2,1)=x0(2)+h*x0(1)^2-u(2)-h*f; 4 b=zeros(N,1); 5 b(1,1)=4*x0(1)-x0(2); 6 b(2,1)=h^2*x0(1)^2-2*x0(1)+x0(2)-h^2*f(1) 7 for i=3:N 8 b(i,1)=x0(i-2)+h^2*x0(i-1)^2-2*x0(i-1)+x0(i)-h^2*f(i-1); 9 end
函数文件2:
1 function g=Jacobian(x0,h,N) 2 % g(1,1)=1; 3 % g(1,2)=-h; 4 % g(2,1)=2*h*x0(1); 5 % g(2,2)=1; 6 g=zeros(N); 7 g(1,1)=4; 8 g(1,2)=-1; 9 g(2,1)=2*h^2*x0(1)-2; 10 g(2,2)=1; 11 for i=3:N 12 g(i,i-2)=1; 13 g(i,i-1)=2*h^2*x0(i-1)-2; 14 g(i,i)=1; 15 end
函数文件3:
1 function x=newton_Iterative_method(f,u,h,N) 2 % u:上一节点的数值解或者初值 3 % x0 每次迭代的上一节点的数值 4 % x1 每次的迭代数值 5 % tol 允许误差 6 % f 右端函数 7 x0=u; 8 tol=1e-11; 9 x1=x0-Jacobian(x0,h,N)F(f,x0,h,N); 10 while (norm(x1-x0,inf)>tol) 11 %数值解的2范数是否在误差范围内 12 x0=x1; 13 x1=x0-Jacobian(x0,h,N)F(f,x0,h,N); 14 end 15 x=x1;%不动点
脚本文件:
1 tic; 2 clear 3 clc 4 N=100; 5 h=1/N; 6 x=0:h:1; 7 y=inline('x.^2.*sin(x).^2+2*cos(x)-x.*sin(x)'); 8 fun1=inline('x.*sin(x)'); 9 Accurate=fun1(x); 10 f=y(x(2:N)); 11 Numerical=zeros(N+1,1); 12 Numerical(2:end,1)=newton_Iterative_method(f,Numerical(2:end,1),h,N); 13 error=Numerical-Accurate'; 14 subplot(1,3,1) 15 plot(x,Accurate); 16 xlabel('x'); 17 ylabel('Accurate'); 18 grid on 19 subplot(1,3,2) 20 plot(x,Numerical); 21 xlabel('x'); 22 ylabel('Numerical'); 23 grid on 24 subplot(1,3,3) 25 plot(x,error); 26 xlabel('x'); 27 ylabel('error'); 28 grid on 29 toc;
效果图: