声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 1290|回复: 0

[编程技巧] 用matlab解超越方程组出现问题,请大家帮忙看看

[复制链接]
发表于 2011-4-1 17:29 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
一共八个方程,构成方程组。未知数是 x,c2,c3,c4,c5,c6,c7,c8。因为c1可以赋值为1.所以不用解。
我也是初学matlab,程序也是参考matlab的帮助写的,但是一直有问题,请大家帮忙看看。不胜感激。

% solve x c2 c3 c4 c5 c6 c7 c8
clear;
clc;

% =============== Define some parameters ================== %
Lambda = 1.55*1E-6;
Epsilon = 8.85*1E-12;
Mu = 4*pi*1E-7;
k = (2*pi)/Lambda;
r1 = 2*1E-5;
r2 = 3*1E-5;
n1 = 1.6;
n2 = 1.49;
syms x c2 c3 c4 c5 c6 c7 c8;
m=0;
Epsilon1 = Epsilon*n1^2;
Epsilon2 = Epsilon*n2^2;
u = (k^2*n1^2 - x^2)^0.5;
w = (x^2 - k^2*n2^2)^0.5;
f = 6*pi/Lambda*1E+8;
c1=1;
% ================= End of Definition ====================== %
% ================definition transform ===================== %
Dbslj = inline('1/2*(besselj(m-1,x)-besselj(m+1,x))');
Dbsly = inline('1/2*(bessely(m-1,x)-bessely(m+1,x))');
Dbsli = inline('1/2*(besseli(m-1,x)+besseli(m+1,x))');
Dbslk = inline('-1/2*(besselk(m-1,x)-besselk(m+1,x))');
% =============end of definition transform====================%
%================definition function==========================%
F1 = c1*besseli(m,w*r1)-c3*besselj(m,u*r1)-c4*bessely(m,u*r1);
F2 = c2*besseli(m,w*r1)-c5*besselj(m,u*r1)-c6*bessely(m,u*r1);
F3 = c3*besselj(m,u*r2)+c4*bessely(m,u*r2)-c7*besselk(m,w*r2);
F4 = c5*besselj(m,u*r2)+c6*bessely(m,u*r2)-c8*besselk(m,w*r2);
F5 = c1*(m*x)/(w^2*r1)*besseli(m,w*r1)+c2*(i/w)*Mu*f*Dbsli(m,w*r1)...
     -c3*(m*x)/(u^2*r1)*besselj(m,u*r1)-c4*(m*x)/(u^2*r1)*bessely(m,u*r1)...
     -c5*(i/u)*Mu*f*Dbslj(m,u*r1)-c6*(i/u)*Mu*f*Dbsly(m,u*r1);
F6 = -c1*(i*f*Epsilon2/w)*Dbsli(m,w*r1)+c2* (m*x)/(w^2*r1)*besseli(m,w*r1)...
     +c3*(i*f*Epsilon1/u)*Dbslj(m,u*r1)+c4*(i*f*Epsilon1/u)*Dbsly(m,u*r1)...
     -c5*(m*x)/(u^2*r1)*besselj(m,u*r1)-c6*(m*x)/(u^2*r1)*bessely(m,u*r1);
F7 = c3*(m*x)/(u^2*r2)*besselj(m,u*r2)+c4*(m*x)/(u^2*r2)*bessely(m,u*r2)...
     +c5*(i*f*Mu/u)*Dbslj(m,u*r2)+c6*(i*f*Mu/u)*Dbsly(m,u*r2)...
     -c7*(m*x)/(w^2*r2)*besselk(m,w*r2)-c8*(i*f*Mu/w)*Dbslk(m,w*r2);
F8 = -c3*(i*f*Epsilon1/u)*Dbslj(m,u*r2)-c4*(i*f*Epsilon1/u)*Dbsly(m,u*r2)...
     +c5*(m*x)/(u^2*r2)*besselj(m,u*r2)+c6*(m*x)/(u^2*r2)*bessely(m,u*r2)...
     +c7*(i*f*Epsilon2/w)*Dbslk(m,w*r2)-c8*(m*x)/(w^2*r2)*besselk(m,w*r2);
%==============================end of function======================%
%============================== program==========================%
function y = myFF(xc);
         y=[F1;F2;F3;F4;F5;F6;F7;F8];
         xc0 = input('[]');
         xc = fsolve(@myFF,xc0,options)


如果把program写到一个m文件里,一直显示function定义有问题,我也尝试吧 program 新建一个m文件,但是解出来的结果根本不对。正确的结果是:x的值应该是介于6000000~6500000之间,有很多个值,也就是说,有很多组解。

求解方程.doc (31 KB, 下载次数: 2)
回复
分享到:

使用道具 举报

您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-11-16 09:41 , Processed in 0.075561 second(s), 21 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表