hom_nonliear.m

来自「同伦方程组求非线性方程组所有解」· M 代码 · 共 28 行

%% N 为积分部数，h 为积分步长，其值视计算精度而定
%% x0 为初始值，可任意给定
%% myfun(x)为非线性方程组表达式（列向量），
%% detamyfun(x)为雅可比矩阵
N=100;
h=1/N;
x0=[rand(1) rand(1) rand(1) rand(1)]';
x=x0;
%format long
f=myfun(x);
b=-h*f;
for i=1:N
  A=detamyfun(x);
  k1=inv(A)*b;
  A=detamyfun(x+0.5*k1);
  k2=inv(A)*b;
  A=detamyfun(x+0.5*k2);
  k3=inv(A)*b;
  A=detamyfun(x+0.5*k3);
  k4=inv(A)*b;
  x=x+(k1+2*k2+2*k3+k4)/6;
end
disp('The Solution is:')
format short;
disp('x=');disp(x');
f=myfun(x)'