Question

matlab解方程问题。

function F=jiefangcheng(x,e,m,b,c,vv12,vv13,vv21,vv23,vv31,vv32,vv44,vv55,vv66,dd12,dd13,dd21,dd23,dd31,dd32,dd44,dd55,dd66)
syms v12 v13 v21 v23 v31 v32 v44 v55 v66 d12 d13 d21 d23 d31 d32 d44 d55 d66
A=[2,-1,-1,0,0,0;-1,2,-1,0,0,0;-1,-1,2,0,0,0;0,0,0,3,0,0;0,0,0,0,3,0;0,0,0,0,0,3];
a=[x(1);x(2);0;0;0;0];
u=1/3*A*a;
V=[0,-v12,-v13,0,0,0;-v21,0,-v23,0,0,0;-v31,-v32,0,0,0,0;0,0,0,v44,0,0;0,0,0,0,v55,0;0,0,0,0,0,v66];
W=[0,-d12,-d13,0,0,0;-d21,0,-d23,0,0,0;-d31,-d32,0,0,0,0;0,0,0,d44,0,0;0,0,0,0,d55,0;0,0,0,0,0,d66];
V=subs(V,{v12,v13,v21,v23,v31,v32,v44,v55,v66},{vv12,vv13,vv21,vv23,vv31,vv32,vv44,vv55,vv66});
W=subs(W,{d12,d13,d21,d23,d31,d32,d44,d55,d66},{dd12,dd13,dd21,dd23,dd31,dd32,dd44,dd55,dd66});
m=V*u;
d=W*u;
sxx1=m(1);
syy1=m(2);
szz1=m(3);
syz1=m(4);
szx1=m(5);
sxy1=m(6);
sxx2=d(1);
syy2=d(2);
szz2=d(3);
syz2=d(4);
szx2=d(5);
sxy2=d(6);
H1=(sxx1+syy1+szz1)/3;
H2=(syz1^2+szx1^2+sxy1^2-syy1*szz1-szz1*sxx1-sxx1*syy1)/3;
H3=(2*syz1*szx1*sxy1+sxx1*syy1*szz1-sxx1*syz1^2-syy1*szx1^2-szz1*sxy1^2)/2;
H4=(sxx2+syy2+szz2)/3;
H5=(syz2^2+szx2^2+sxy2^2-syy2*szz2-szz2*sxx2-sxx2*syy2)/3;
H6=(2*syz2*szx2*sxy2+sxx2*syy2*szz2-sxx2*syz2^2-syy2*szx2^2-szz2*sxy2^2)/2;
p1=H1+H2;
q1=(2*H1^3+3*H1*H2+2*H3)/2;
o1=acos(q1/p1^(3/2));
p2=H4+H5;
q2=(2*H4^3+3*H4*H5+2*H6)/2;
o2=acos(q2/p2^(3/2));
s1=2*(H1^2+H2)^(1/2)*cos(o1/3)+H1;
s2=2*(H1^2+H2)^(1/2)*cos((o1+4*pi)/3)+H1;
s3=2*(H1^2+H2)^(1/2)*cos((o1+2*pi)/3)+H1;
s4=2*(H4^2+H5)^(1/2)*cos(o2/3)+H4;
s5=2*(H4^2+H5)^(1/2)*cos((o2+4*pi)/3)+H4;
s6=2*(H4^2+H5)^(1/2)*cos((o2+2*pi)/3)+H4;
fangcheng1=(abs(s1-s4))^8+(abs(s1-s5))^8+(abs(s1-s6))^8+(abs(s2-s4))^8+(abs(s2-s5))^8+(abs(s2-s6))^8+(abs(s3-s4))^8+(abs(s3-s5))^8+(abs(s3-s6))^8;
fangcheng2=diff(fangcheng1,x(2));
fangcheng3=diff(fangcheng1,x(1));
F(1)=fangcheng1-507.7*(0.005+x(3))^(0.280);
F(2)=sqrt(x(1)^2*(cos(e/2)^2-cos(e/2)^4))-200;
F(3)=fangcheng2/fangcheng3-m;
这是m文件

然后输入程序 哪位大神能帮我修改一下
options=optimset('display','off');
x0=[100,100,0.1];
n=1:134;
m=1-1.5*(n-1)/200;
a=0.000001;
b=m*a;
c=-a-b;
e=acos((a+c)./(a-c));
vv12=1;vv13=1;vv21=1;vv23=1;vv31=1;vv32=1;vv44=1;vv55=1;vv66=1;dd12=1;dd13=1;dd21=1;dd23=1;dd31=1;dd32=1;dd44=1;dd55=1;dd66=1;
for i=1:length(n)
mm=m(i);
bb=b(i);
cc=c(i);
ee=e(i);
x=fsolve(@(x) jiefangcheng(x,ee,mm,cc,bb,vv12,vv13,vv21,vv23,vv31,vv32,vv44,vv55,vv66,dd12,dd13,dd21,dd23,dd31,dd32,dd44,dd55,dd66),x0,options);
xx(i)=x(1);
yy(i)=x(2);
zz(i)=x(3);
end
ww=zz.*(2*yy-xx)
vv=2*sqrt(xx.^2-xx.*yy+yy.^2)
uu=zz.*(2*xx-yy)
xxx=uu./vv
yyy=ww./vv
plot(yyy,xxx,'.')

??? Error using ==> unknown
Matrix dimensions must agree.

Error in ==> jiefangcheng at 47
F(3)=fangcheng2/fangcheng3-m;

Error in ==> fsolve at 180
fuser = feval(funfcn{3},x,varargin{:});

Answer 1

发现解在+/-1.6附近 因此通过fsolve 来解非线性方程，不能得到其解析解，只能计算其数值解，通过限制其精度，达到自己的要求。 fsolve('funx',-1.7)即可求得 -1.6376 - 0.0000i fsolve('funx',1.5) 1.6376 + 0.0000i function [fx ] = funx( n ) %FUNX Summary of this function goes here % Detailed explanation goes here lanmuda=1.55; n1=1.6375696780514445567155543374853; n2=1.6443618413936927969908352393948; n3=1.6375696780514445567155543374853; c1=n2^2/n1^2; c2=n2^2/n3^2; k0=2*pi/lanmuda; r1=sqrt(n.^2-n1^2)*k0; r2=sqrt(n2^2-n.^2)*k0; r3=sqrt(n.^2-n3^2)*k0; fx=r2.*5-atan(c1.*r1./r2)-atan(c2.*r3./r2)-pi; end