zgq1004118 发表于 2007-2-23 08:00

在Matlab中解非线性方程组遇到的问题!

刚刚接触MATLAB不久,属于 菜菜鸟 级的。在求解一组非线性方程时无法收敛!是否有其他解法?请高人指教!!
function non_linearfun
x0=;
answer_x_y_z=fsolve(@fun,x0)
function y=fun(x)
%----------------------Variable-------------------
%x=sigma_r./A;
%y=sigma_a./A
%z=Theta/A
%

Kn      =   434912.7857;
A       =   0.9922;
B       =   0.05;
D       =   19.844;
Ri      =   40.4027718;
dm      =   80;
n       =   1.5;
alpha0=   40/180*pi;
%
%---------------------------------------------------
sinao   =   sin(alpha0);
cosao   =   cos(alpha0);

phi0    =   0;
phi1    =   2/11*pi;
phi2    =   2*2/11*pi;
cosphi0 =   cos(phi0);
cosphi1 =   cos(phi1);
cosphi2 =   cos(phi2);

%------------------------------------------------
B0=sinao+x(1)+Ri*x(3)*cosphi0;
B1=sinao+x(1)+Ri*x(3)*cosphi1;
B2=sinao+x(1)+Ri*x(3)*cosphi2;
%------------------------------------
%------------------------------------
C0=cosao+x(2)*cosphi0;
C1=cosao+x(2)*cosphi1;
C2=cosao+x(2)*cosphi2;
%----------------------------------
%----------------------------------
D0=(cosao+x(2)*cosphi0)*cosphi0;
D1=(cosao+x(2)*cosphi1)*cosphi1;
D2=(cosao+x(2)*cosphi2)*cosphi2;

%-----------------------------------
Fa      =   10000;
Fr      =   0;
My      =   0;                                    
%---------------------------------------

y(1)=   Fa-Kn*(A^n)*(((B0^2+C0^2)^.5-1)^1.5*B0/(B0^2+C0^2)^0.5+2*((B1^2+C1^2)^0.5-1)^n*B1/(B1^2+C1^2)^0.5+2*((B2^2+C2^2)^0.5-1)^n*B2/(B2^2+C2^2)^.5);
y(2)=   Fr-Kn*(A^n)*(((B0^2+C0^2)^.5-1)^1.5*D0/(B0^2+C0^2)^0.5+2*((B1^2+C1^2)^0.5-1)^n*D1/(B1^2+C1^2)^0.5+2*((B2^2+C2^2)^0.5-1)^n*D2/(B2^2+C2^2)^.5);
y(3)=   My-0.5*dm*Kn*(A^n)*(((B0^2+C0^2)^.5-1)^1.5*B0*cosphi0/(B0^2+C0^2)^0.5+2*((B1^2+C1^2)^0.5-1)^n*B1*cosphi1/(B1^2+C1^2)^0.5+2*((B2^2+C2^2)^0.5-1)^n*B2*cosphi2/(B2^2+C2^2)^.5);

zgq1004118 发表于 2007-2-23 23:41

说明

我再解释一下:@) ,三个未知数x(1),x(2),x(3).
B0-B2,C0-C2,D0-D2是含有未知数的变量。
其它都是常量,其数值已经给出。这是一个很复杂的方程组,已被我简化写出。
高人啊,快来帮帮兄弟啊:@L

xjzuo 发表于 2007-2-25 18:38

请把你的"无法收敛"贴出来看看.
另:我运行了一下,似乎可以得到结果.

zgq1004118 发表于 2007-2-28 01:13

> Warning: Cannot determine from calling sequence whether to use new
(2.0 or later) FSOLVE function or grandfathered FSOLVE function.
Assuming new syntax; if call was grandfathered FSOLVE syntax,
this may give unexpected results.
To force new syntax and avoid warning, add options structure argument:
      x = fsolve(@sin,3,optimset('fsolve'));
To force grandfathered syntax and avoid warning, add options array argument:
      x = fsolve(@sin,3,foptions);

> In E:\Matlab\toolbox\optim\fsolve.m (parse_call) at line 389
In E:\Matlab\toolbox\optim\fsolve.m at line 101
In E:\Matlab\work\non_linearfunz6.m at line 3
Optimization terminated successfully:
Norm of relative change in X is less than max(options.TolX^2,eps) and
sum-of-squares of function values is less than sqrt(options.TolFun)

answer_x_y_z =

Column 1

          0.261630241614871 +    0.0432136503412044i

Column 2

         -0.147501542528093 -    0.0461512387637423i

Column 3

       -0.00424579357031622 -   0.00155455859937376i
这是我运行的结果,是完全数解,得不到实数解,请问高人的解怎样?

xjzuo 发表于 2007-2-28 08:32

回复

由于你的问题背景没有讲清楚,所以不知道你要求什么样的解.
我试了试以下指令,就你给的这段简化程序而言,似乎并非你说的不收敛:
%%%-------------------------------------------------------%%%
x0=;
x = fsolve(@myfunt1,x0,optimset('fsolve'))
%%%-------------------------------------------------------%%%
x0=;
x = fsolve(@myfunt1,x0,optimset('fsolve'))
%%%-------------------------------------------------------%%%

zgq1004118 发表于 2007-3-1 02:27

多谢版主关心。我现在把方程换个写法。三个未知数:x(1),x(2),x(3);三个方程是y(1),y(2),y(3);
欲求出一组实数解!!这样写出来便于理解方程组的结构,就是方程本身太复杂,不便于在matlab中编程,因此就没有写出原始表达式。
不知版主上次选取x0=做初值,在运行时的结果是否跟我的一样也是完全数解?

function non_linearfun
x0=;
answer_x_y_z=fsolve(@fun,x0)
function y=fun(x)
%----------------------function-------------------
y(1)=10000-429834.2416*((((.642787609+x(1)+40.304027718*x(3))^2+(x(2)+.766044443)^2)^.5-1)^1.5*(.642787609+x(1)+40.304027718*x(3))/((.642787609+x(1)+40.304027718*x(3))^2+(x(2)+.766044443)^2)^.5+2*(((.642787609+x(1)+33.9059057*x(3))^2+(.841253532*x(2)+.766044443)^2)^.5-1)^1.5*(.642787609+x(1)+33.9059057*x(3))/((.642787609+x(1)+33.9059057*x(3))^2+(.841253532*x(2)+.766044443)^2)^.5+2*(((.642787609+x(1)+16.7428982*x(3))^2+(.415415013*x(2)+.766044443)^2)^.5-1)^1.5*(.642787609+x(1)+16.7428982*x(3))/((.642787609+x(1)+16.7428982*x(3))^2+(.415415013*x(2)+.766044443)^2)^.5);
y(2)=0-429834.2416*((((.642787609+x(1)+40.304027718*x(3))^2+(x(2)+.766044443)^2)^.5-1)^1.5*(.766044443+x(2))/((.642787609+x(1)+40.304027718*x(3))^2+(x(2)+.766044443)^2)^.5+2*(((.642787609+x(1)+33.9059057*x(3))^2+(.841253532*x(2)+.766044443)^2)^.5-1)^1.5*(.766044443+.841253532*x(2))*.841253532/((.642787609+x(1)+33.9059057*x(3))^2+(.841253532*x(2)+.766044443)^2)^.5+2*(((.642787609+x(1)+16.7428982*x(3))^2+(.415415013*x(2)+.766044443)^2)^.5-1)^1.5*(.766044443+.415415013*x(2))*.415415013/((.642787609+x(1)+16.7428982*x(3))^2+(.415415013*x(2)+.766044443)^2)^.5);
y(3)=0-40*429834.2416*((((.642787609+x(1)+40.304027718*x(3))^2+(x(2)+.766044443)^2)^.5-1)^1.5*(.642787609+x(1)+40.304027718*x(3))/((.642787609+x(1)+40.304027718*x(3))^2+(x(2)+.766044443)^2)^.5+2*(((.642787609+x(1)+33.9059057*x(3))^2+(.841253532*x(2)+.766044443)^2)^.5-1)^1.5*(.642787609+x(1)+33.9059057*x(3))*.841253532/((.642787609+x(1)+33.9059057*x(3))^2+(.841253532*x(2)+.766044443)^2)^.5+2*(((.642787609+x(1)+16.7428982*x(3))^2+(.415415013*x(2)+.766044443)^2)^.5-1)^1.5*(.642787609+x(1)+16.7428982*x(3))*.415415013/((.642787609+x(1)+16.7428982*x(3))^2+(.415415013*x(2)+.766044443)^2)^.5);

xjzuo 发表于 2007-3-1 08:56

回复

我认为还是一楼的表达式比较清楚.
就你给的表达式,似乎没有实数解.
另:你说的"完全数",应该不是指数论中的完全数吧?

[ 本帖最后由 xjzuo 于 2007-3-1 08:59 编辑 ]

dingd 发表于 2007-3-1 09:52

第一个公式中:

"10000-429834.2416*(...)"变为“1000000-429834.2416*()"或”100000-42983.42416*()'就会有实数解了,否则,似乎无实数解。确认代码没输错吗?

zgq1004118 发表于 2007-3-2 01:48

8楼的高手我前不久在本论坛读过你的关于非线性方程组解法的帖子,从中收获很多,在这谢了。
下面是把"10000-429834.2416*(...)"变为“1000000-429834.2416*()"后的结果,
仍然不是实数解。不知高手的结如何,能否贴出来?另外是否有其它解法?:@)
answer_x_y_z =

Column 1

         3.64616251216352 -   0.644111407183738i

Column 2

         -0.742140780488419 +   0.147938243941897i

Column 3

         -0.100008046487466 +    0.0360334493768073i

dingd 发表于 2007-3-2 10:23

"10000-429834.2416*(...)"变为“1000000-429834.2416*(...)",结果:

x1: 5.44662954384634
x3: -0.210122988803927
x2: -1.06537137387729

dingd 发表于 2007-3-2 10:28

再精确点:

x1: 5.44663122246003
x3: -0.210123048907879
x2: -1.06536950937869

zgq1004118 发表于 2007-3-5 04:39

多谢这位高人,但10楼和11楼的解也不在取值范围内,另外我用MATLAB是得不到这个解的。郁闷啊:@( :@o

grta 发表于 2007-3-7 22:44

用solve函数试试那
页: [1]
查看完整版本: 在Matlab中解非线性方程组遇到的问题!