返回信息流想用matlab解一组非线性方程组（超越方程，含有Bessel函数）
一直编不好，得不出解，我用的书上说的Broyden法（秩1的拟Newton法）
请各位大牛帮我看看怎么改，我比较菜，现学现用，大家别见笑
主程序部分：
Lamda=0.3e-6;%波长
c=3.0e8;%真空中的光速
n1=1.0;%空气折射率
n2=1.46;%石英折射率
r=0.3e-6;%空气孔半径
Delta=2.3e-6;%孔间距
R=0.525*Delta;%单元胞外半径
Omega=2*pi*c/Lamda;
broyden(1)
函数部分
function y=broyden(x0)
%y=broyden(x0) x0为初值
a=eye(length(x0));
x1=x0-fc(x0)/a;
n=1;
while (norm(x1-x0)>=1.0e-6)&(n<=100000000)
x0=x1;
x1=x0-fc(x0)/a;
p=x1-x0;
q=fc(x1)-fc(x0);
a=a+(q-p*a)'*p/norm(p);
n=n+1;
end
y=x1;
n
function y=fc(x)
%u=x(1);
%w=x(2);
g=1/(x(2)*r)*(BESSELJ(0,x(1)*r)*BESSELY(1,x(1)*R)-BESSELY(0,x(1)*r)*BESSELJ(1,x(1)*R))/(BESSELJ(1,x(1)*r)*BESSELY(1,x(1)*R)-BESSEL(1,x(1)*r)*BESSELJ(1,x(1)*R))-1/(x(1)^2*r^2);
f=1/r^4*(1/x(1)^2+1/x(2)^2)*(n2^2/x(1)^2+n1^2/x(2)^2);
y(1)=BESSELI(2,x(2)*r)/BESSELI(1,x(2)*r)+1/(x(2)*r)+x(2)*r/2*(1+n2^2/n1^2)*g+x(2)*r*sqrt(1/4*(1-n2^2/n1^2)*g^2+f/n1^2);
y(2)=x(2)^2+x(1)^2-Omega^2*(n2^2-n1^2)/c^2;
y=[y(1) y(2)];
这是一条镜像帖。来源：北邮人论坛 / matlab / #4507同步于 2008/12/19
该镜像源已超过 30 天没有更新，可能在源站已被删除。
Matlab机器人发帖
非线性方程组请教
Elim
2008/12/19镜像同步0 回复
订阅后，新回复会通过你的通知中心匿名送达。
0 条回复
暂无回复 · 你可以订阅本帖等待新回复。