牛顿法解非线性方程组结果为NaN（已解决）

原野05071515

function f=f1(x0)
x=x0(1);
y=x0(2);
z=x0(3);
f1=x^2+y^2+z^2-1;
f2=2*x^2+y^2-4*z;
f3=3*x^2-4*y+z^2;
f=[f1 f2 f3];

function f=df1(x0)
x=x0(1);
y=x0(2);
z=x0(3);
f=[2*x 2*y 2*z
   4*x 2*y -4
   6*x -4 2*z];

function [x,n,data]=mulNewton(x0,tol)
if nargin==1
    tol=1.0e-4;
end
x1=x0-f1(x0)/df1(x0);
n=1;
while (norm(x1-x0)>tol)
    x0=x1;
    x1=x0-f1(x0)/df1(x0);
    n=n+1;
    data(:,n)=x1;
end
x=x1;

x0=[1e10 1e10 1e10];
[x,n,data]=mulNewton(x0);
disp(x)
disp(n)


为什么求出的解为以下形式：
>> clear
new_main
   NaN   NaN   NaN

   244
希望大侠们能帮下忙~

原野05071515

我把程序改了一下就出结果：
function f=f1(x0)
x=x0(1);
y=x0(2);
z=x0(3);
f1=x^2+y^2+z^2-1;
f2=2*x^2+y^2-4*z;
f3=3*x^2-4*y+z^2;
f=[f1;f2;f3];

function F=df1(x0)
x=x0(1);
y=x0(2);
z=x0(3);
F=[2*x 2*y 2*z
4*x 2*y -4
6*x -4 2*z];

function [x,n,data]=new_ton(x0,tol)
if nargin==1
    eps=1.0e-4;
end
x1=x0-inv(df1(x0))*f1(x0);
n=1;
tol=1;
while (tol>eps)
    x0=x1;
    x1=x0-inv(df1(x0))*f1(x0);
    tol=norm(x1-x0);
    n=n+1;
    data(:,n)=x1;
end
x=x1;


主函数：
x0=[7.8e-1;4.9e-1;3.7e-1];
[x,n,data]=new_ton(x0);
disp('迭代结果')
x
disp('迭代次数')
n

结果如下：

迭代结果

x =

    0.7852
    0.4966
    0.3699

迭代次数

n =
    2