声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 2312|回复: 7

[综合讨论] 用fsolve 解非线性方程组出现问题

[复制链接]
发表于 2008-1-23 16:12 | 显示全部楼层 |阅读模式

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

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

x
代码为:
format long
Ul11=9.92;
Ul21=8.72;
Ul12=13.2;
Ul22=12.6;
Ul13=19.2;
Ul23=20.0;
options=optimset('LargeScale','off','NonlEqnAlgorithm','gn','MaxIter',100000,'MaxFunEvals',10000000);
x0=[0.0001,0.0001,60];
x=fsolve(@myfun2,x0,options);
alpha1=x(1)
beta1=x(2)
Zc1=x(3)
Ul11=9.92;
Ul21=8.48;
Ul12=12.4;
Ul22=12.2;
Ul13=17.6;
Ul23=19.6;
x0=[0.0001,0.0001,60];
x=fsolve(@myfun2,x0,options);
alpha2=x(1)
beta2=x(2)
Zc2=x(3)

w1=2*pi*100*10^3;
w2=2*pi*200*10^3;
galpha1=alpha1;
gbeta1=beta1;
galpha2=alpha2;
gbeta2=beta2;
x0=[2.5*1.000000e-011; 2.5*1.000000e-005; 0.5*1.000000e-006 ;1.000000e-008];
[x,fval] = fsolve(@mymean,x0,options);
C0=x(1);
G0=x(2);
L0=x(3);
R0=x(4);
C0
G0
L0
R0
omga1=w2;
[Z0,Y0,gama,alpha,beta,Zc]=trans_line(R0,G0,C0,L0,omga1)
V=omga1/beta     
函数为:
function F=myfun2(x)
%Ul11,Ul21分别是终端接阻抗Z1的条件下中点和终端处的电压;
%Ul12,Ul22分别是终端接阻抗Z2的条件下中点和终端处的电压;
%Ul13,Ul23分别是终端开路的条件下中点和终端处的电压;
%x(1)、x(2)分别频率为衰减常数和相移常数,x(3)为特性阻抗
global Ul11 Ul21 Ul12 Ul22 Ul13 Ul23
l=44.84;
Z1=50;
Z2=100;
F=[Ul11-Ul21*abs(cosh(l*(x(1)+j*x(2)))+x(3)/Z1*sinh(x(1)+j*x(2)));
    Ul12-Ul22*abs(cosh(l*(x(1)+j*x(2)))+x(3)/Z2*sinh(x(1)+j*x(2)));
    Ul13-Ul23*abs(cosh(l*(x(1)+j*x(2))))];


function [Z0,Y0,gama,alpha,beta,Zc]=trans_line(R0,G0,C0,L0,omga)
%R0:两根导线每单位长度具有的电阻;
%L0:两根导线每单位长度具有的电感;C0: 每单位长度导线之间具有的电容;
%G0:每单位长度导线之间具有的电导;
Z0=R0+j*omga*L0 ;             %单位长度的阻抗
Y0=G0+j*omga*C0;              %单位长度的导纳
gama=sqrt(Z0*Y0);            %gama为传播常数  
alpha=real(gama);             %aerf为衰减常数   
beta=imag(gama);             %beta 为相移常数
Zc=sqrt(Z0/Y0);              % Zc为特性阻抗

function CGRL=mymean(x)
global w1 w2 gbeta1 galpha1 gbeta2 galpha2
    CGRL=[(x(4)/2*sqrt(x(1)/x(3))+x(2)/2*sqrt(x(3)/x(1)))*(1-(x(4)/x(3)-x(2)/x(1))^2/(8*w1^2))-galpha1;
      w1*sqrt(x(3)*x(1))*(1+(x(4)/x(3)-x(2)/x(1))^2/(8*w1^2))-gbeta1;
      (x(4)/2*sqrt(x(1)/x(3))+x(2)/2*sqrt(x(3)/x(1)))*(1-(x(4)/x(3)-x(2)/x(1))^2/(8*w2^2))-galpha2;
   w2*sqrt(x(3)*x(1))*(1+(x(4)/x(3)-x(2)/x(1))^2/(8*w2^2))-gbeta2];


代码为上面所示,运行时 前两个fsolve函数,初值给什么值结果就是什么值,而且提示Optimization terminated: directional derivative along search direction less than TolFun and infinity-norm of gradient less than 10*(TolFun+TolX).
第三个fsolve出现提示Optimization terminated: search direction less than TolX.


哪位高手告诉我程序错在哪?另外怎么确定初值?

[ 本帖最后由 eight 于 2008-1-23 18:40 编辑 ]
回复
分享到:

使用道具 举报

发表于 2008-1-23 18:40 | 显示全部楼层
这个只是 warning,不是 error。请 help fsolve 一下,或在 matlab 中根据提示的某些关键词搜索一下。还可以搜索版面

[ 本帖最后由 eight 于 2008-1-23 18:41 编辑 ]
发表于 2008-1-24 08:55 | 显示全部楼层
试一下1stOpt,不需要猜初值。
 楼主| 发表于 2008-1-24 10:09 | 显示全部楼层
但是结果很显然不对    Zc的大小应该在45-----100之间,
1stopt 我装了一个  但不会用   相关的教程没找到
 楼主| 发表于 2008-1-25 09:59 | 显示全部楼层
请问用1stopt解方程组的时候,结果中目标函数值不为0,得到的结果可信可靠吗?
发表于 2008-1-25 13:49 | 显示全部楼层
不为0或不接近于0,表明结果有差异。多算两次。还有2.0以前的版本,对非线性方程组,求解能力要弱不少。
发表于 2008-1-25 15:15 | 显示全部楼层
这个问题用matlab可以解决,但用1stopt 可能更合适一些。后者的语法可以参考论坛的相关帖子,以及1stopt 自带的例子。其语法基本上是傻瓜型地套用,自己稍微看看就会明白,不过有时候自己先(化简)处理一下问题,可能会解决得更快。

另:问问题时,将公式贴出,并描述一下问题,是一个好习惯,也是一个可以更快得到解答的方法。

[ 本帖最后由 xjzuo 于 2008-1-25 15:18 编辑 ]
 楼主| 发表于 2008-1-25 17:00 | 显示全部楼层
谢谢大家!
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

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

GMT+8, 2024-11-18 10:26 , Processed in 0.068704 second(s), 17 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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