用牛顿迭代法解三元三次方程的MATLAB程序!

用牛顿迭代法解三元三次方程的MATLAB程序!
用牛顿迭代法解下面的三元三次方程,用MATLAB实现,急用,
(1-c)*(1-y)*(1-m)*10.61+y(1-c)*(1-m)*30.56+c*(1-y)*(1-m)*60.67+m*(1-c)*(1-y)*0.95+y*m*(1-c)*28.45+y*c*(1-m) *8.19+c*m*(1-y) *2.51+c* y*m* 76.69-3.21=0
(1-c)*(1-y)*(1-m)*15.21+y*(1-c)*(1-m)*14.58+c*(1-y)*(1-m)*65.2+m*(1-c)*(1-y)*0.96+y*m*(1-c)*14.23+y*c*(1-m) *17.41+c*m*(1-y) *1.62+m*y *c*80.32=0
(1-c)*(1-y)*(1-m)*44.77+y*(1-c)*(1-m)*9.9+c*(1-y)*(1-m)*3.97+m*(1-c)*(1-y)*0.78+y*m*(1-c)*0.8+y*c*(1-m) *3.3+c*m*(1-y) *7.95+ m*c*y*66.79=0
c和m和y分别都是未知数,解的范围都是从0-1.
小妹财富不多,
jimtien的解答,我的方程漏掉了一些数据,在这里我把完整的写出来,此方程中X、Y、Z是已知数据,我的意思是通过这个方程,将一组X、Y、Z数据,转换成c、m、y数据.能给出一个能加载 三维数据的 最终程序么?
(1-c)*(1-y)*(1-m)*10.61+y(1-c)*(1-m)*30.56+c*(1-y)*(1-m)*60.67+m*(1-c)*(1-y)*0.95+y*m*(1-c)*28.45+y*c*(1-m) *8.19+c*m*(1-y) *2.51+c* y*m* 76.69-X=0
(1-c)*(1-y)*(1-m)*15.21+y*(1-c)*(1-m)*14.58+c*(1-y)*(1-m)*65.2+m*(1-c)*(1-y)*0.96+y*m*(1-c)*14.23+y*c*(1-m) *17.41+c*m*(1-y) *1.62+m*y *c*80.32-Y=0
(1-c)*(1-y)*(1-m)*44.77+y*(1-c)*(1-m)*9.9+c*(1-y)*(1-m)*3.97+m*(1-c)*(1-y)*0.78+y*m*(1-c)*0.8+y*c*(1-m) *3.3+c*m*(1-y) *7.95+ m*c*y*66.79-Z=0
补充：X Y Z的参考数据
X Y Z
38.98 38.38 32.29
7.07 9.73 10.97
0.67 0.85 0.79
69.85 68.99 56.01
1.11 1.7 1.6
52.05 54.56 49.08
3.21 2.15 3.76
44.6 46.98 41.77
1.21 1.56 1.46
63.46 63.94 51
3 4.05 4.41
13.01 13.6 11.03

用fsolve可解出来：
先构造函数：
function output=solveproblem(X)
c=X(1);
m=X(2);
y=X(3);
output(1)=(1-c)*(1-y)*(1-m)*10.61+y*(1-c)*(1-m)*30.56+c*(1-y)*(1-m)*60.67+m*(1-c)*(1-y)*0.95+y*m*(1-c)*28.45+y*c*(1-m) *8.19+c*m*(1-y) *2.51+c* y*m* 76.69-3.21
output(2)=(1-c)*(1-y)*(1-m)*15.21+y*(1-c)*(1-m)*14.58+c*(1-y)*(1-m)*65.2+m*(1-c)*(1-y)*0.96+y*m*(1-c)*14.23+y*c*(1-m) *17.41+c*m*(1-y) *1.62+m*y *c*80.32
output(3)=(1-c)*(1-y)*(1-m)*44.77+y*(1-c)*(1-m)*9.9+c*(1-y)*(1-m)*3.97+m*(1-c)*(1-y)*0.78+y*m*(1-c)*0.8+y*c*(1-m) *3.3+c*m*(1-y) *7.95+ m*c*y*66.79
然后求
x0=[1 1 1]; %初值
ans1 = ones(20,3);
options=optimset('fsolve');
for k=1:20
[x,fval,exitflag,output]=fsolve(@solveproblem,x0,options);
x0=x;
ans1(k,:)=x;
end
ans1 =
-0.2892 0.9143 0.1893
带入三个方程验证,分别等于0.0012,4.8363e-004,-0.0016
再问： 能看下我的问题补充，非常感谢~
再答： function output=solveproblem(X1,Y1) c=X1(1); m=X1(2); y=X1(3); X=Y1(1); Y=Y1(2); Z=Y1(3); output(1)=(1-c)*(1-y)*(1-m)*10.61 + y*(1-c)*(1-m)*30.56 + c*(1-y)*(1-m)*60.67 + m*(1-c)*(1-y)*0.95 + y*m*(1-c)*28.45+y*c*(1-m) *8.19 + c*m*(1-y)*2.51 + c*y*m*76.69 - X output(2)=(1-c)*(1-y)*(1-m)*15.21 + y*(1-c)*(1-m)*14.58 + c*(1-y)*(1-m)*65.2 + m*(1-c)*(1-y)*0.96 + y*m*(1-c)*14.23+y*c*(1-m) *17.41+ c*m*(1-y)*1.62 + m*y*c*80.32 - Y output(3)=(1-c)*(1-y)*(1-m)*44.77 + y*(1-c)*(1-m)*9.9 + c*(1-y)*(1-m)*3.97 + m*(1-c)*(1-y)*0.78 + y*m*(1-c)*0.8+y*c*(1-m) *3.3 + c*m*(1-y)*7.95 + m*c*y*66.79 - Z ----------------------------------------- XYZ=[38.98 38.38 32.29; 7.07 9.73 10.97; 0.67 0.85 0.79; 69.85 68.99 56.01; 1.11 1.7 1.6; 52.05 54.56 49.08; 3.21 2.15 3.76; 44.6 46.98 41.77; 1.21 1.56 1.46; 63.46 63.94 51; 3 4.05 4.41; 13.01 13.6 11.03]; [m,n]=size(XYZ); x0=[1 1 1]; %初始值 ans1 = [zeros(m,3) XYZ]; %每行1~3列保存结果，4~6列保存X,Y,Z值 options = optimset('fsolve'); for k=1:m [x,fval,exitflag,output]=fsolve(@solveproblem,x0,options,XYZ(k,:)); ans1(k,1:3)=x; end