mongolism 发表于 2008-3-28 23:00

请教方程组求解

最近在求解一个非常复杂的非线性四个未知数的方程组,但是我在这个方程组的外面是采用循环的的方式进行多次的求解,因为这个方程组的部分系数是变化的,由循环的此时确定,但是在某一次既定的循环中,所有的系数都是已知的,因此方程是可解的,问题就出来了,我的目的是进行10000次循环,在每一次循环中我都会观察计算结果,结果应该是四个数值,在前面的52次(这个次数不是很重要)循环中都得到了正确的结果,但是在第53次循环中出现了这样一个问题,没有给出数值解,是有一堆凌乱的公式,和前面的结果明显不一样,所以我知道了这个计算程序已经不能给出正确的结果了,但是如果我单独把这次的循环中的各个系数提取出来,在一个新的窗口中计算这个方程组,则是可以得到一个准确的结果的,为什么?下面是我的程序和运行错误信息,请指教:
restart;
>
> restart;

> for l from 1 to 1 by 1 do

> for m from 1 to 100 by 1 do

>   k:=l*0.004;

>   zeta:=m*0.001;

>   center:=1+k^2-2*zeta^2;
>   
eqn1:=(sqrt(center+x)^4-sqrt(center+x)^2*c^2-2*sqrt(center+x)^2*zeta*c^2*d-sqrt(center+x)^2-sqrt(center+x)^2*k^2+c^2)^2+(-2*sqrt(center+x)^3*zeta-sqrt(center+x)^3*c^2*d+2*sqrt(center+x)*zeta*c^2+c^2*d*sqrt(center+x)+sqrt(center+x)*c^2*d*k^2)^2-((-sqrt(center+x)^2+c^2)^2+(c^2*d*sqrt(center+x))^2)*((sqrt(center+x)^2-1-k^2)^2+(-2*sqrt(center+x)*zeta)^2);

>eqn2:=(sqrt(center-x)^4-sqrt(center-x)^2*c^2-2*sqrt(center-x)^2*zeta*c^2*d-sqrt(center-x)^2-sqrt(center-x)^2*k^2+c^2)^2+(-2*sqrt(center-x)^3*zeta-sqrt(center-x)^3*c^2*d+2*sqrt(center-x)*zeta*c^2+c^2*d*sqrt(center-x)+sqrt(center-x)*c^2*d*k^2)^2-((-sqrt(center-x)^2+c^2)^2+(c^2*d*sqrt(center-x))^2)*((sqrt(center-x)^2-1-k^2)^2+(-2*sqrt(center-x)*zeta)^2);

>eqn3:=(sqrt(1+k^2-2*zeta^2)^4-sqrt(1+k^2-2*zeta^2)^2*c^2-2*sqrt(1+k^2-2*zeta^2)^2*zeta*c^2*d-sqrt(1+k^2-2*zeta^2)^2-sqrt(1+k^2-2*zeta^2)^2*k^2+c^2)^2+(-2*sqrt(1+k^2-2*zeta^2)^3*zeta-sqrt(1+k^2-2*zeta^2)^3*c^2*d+2*sqrt(1+k^2-2*zeta^2)*zeta*c^2+c^2*d*sqrt(1+k^2-2*zeta^2)+sqrt(1+k^2-2*zeta^2)*c^2*d*k^2)^2-((-sqrt(1+k^2-2*zeta^2)^2+c^2)^2+(c^2*d*sqrt(1+k^2-2*zeta^2))^2)*((sqrt(center-x)^2-1-k^2)^2+(-2*sqrt(center-x)*zeta)^2):

>z:=fsolve({eqn1,eqn2,eqn3},{x,c,d},{x=0..2,c=0..2,d=0..1});

> print(z);
restart;
> od;

> od;
   {c = 1.000000343, d = 0.005656862445, x = 0.004394748134}
   {d = 0.005656832890, x = 0.005534166657, c = 0.9999886865}
   {c = 0.9999730302, d = 0.005656773773, x = 0.006476053075}
   {d = 0.005656648276, x = 0.007297256086, c = 0.9999533744}
   {x = 0.008034792781, c = 0.9999297197, d = 0.005656354291}
   {c = 0.9999020668, d = 0.005655844108, x = 0.008709852888}
   {c = 0.9998704165, d = 0.005655048101, x = 0.009335896289}
   {c = 0.9998347703, d = 0.005653874224, x = 0.009922060245}
   {c = 0.9997951295, d = 0.005652256234, x = 0.01047490051}
   {c = 0.9997514960, x = 0.01099927785, d = 0.005650120785}
   {c = 0.9997038725, d = 0.005647308214, x = 0.01149882736}
   {c = 0.9996522622, d = 0.005643706383, x = 0.01197640559}
   {c = 0.9995966678, d = 0.005639230319, x = 0.01243430139}
   {c = 0.9995370934, d = 0.005633748277, x = 0.01287430511}
   {c = 0.9994735443, d = 0.005627063940, x = 0.01329778253}
   {c = 0.9994060251, d = 0.005619082078, x = 0.01370593370}
   {c = 0.9993345413, d = 0.005609659459, x = 0.01409968603}
   {x = 0.01447968170, c = 0.9992591006, d = 0.005598567937}
   {c = 0.9991797097, d = 0.005585679622, x = 0.01484652481}
   {c = 0.9990963784, d = 0.005570727597, x = 0.01520051673}
/                                    -7                  \
{ d = 0.2669176555, x = 8.495593356 10, c = 0.001658312396 }
\                                                          /
   {c = 0.9989179288, d = 0.005534049763, x = 0.01587120454}
   {c = 0.9988228327, d = 0.005511856317, x = 0.01618806304}
   {c = 0.9987238396, d = 0.005486749433, x = 0.01649251654}
   {c = 0.9986209626, d = 0.005458505500, x = 0.01678445332}
   {d = 0.005426865389, x = 0.01706363663, c = 0.9985142167}
   {c = 0.9984036193, d = 0.005391513081, x = 0.01732967899}
   {c = 0.9982891886, d = 0.005352152505, x = 0.01758214880}
   {c = 0.9981709434, d = 0.005308499985, x = 0.01782056379}
   {c = 0.9980489077, d = 0.005260138263, x = 0.01804416250}
   {c = 0.9979231049, d = 0.005206713333, x = 0.01825219817}
   {c = 0.9977935617, d = 0.005147820446, x = 0.01844376691}
   {c = 0.9976603078, d = 0.005083003682, x = 0.01861778395}
   {c = 0.9975233750, d = 0.005011797600, x = 0.01877305142}
   {c = 0.9973827985, d = 0.004933694252, x = 0.01890818383}
   {c = 0.9972386202, d = 0.004848049515, x = 0.01902141845}
   {c = 0.9970908859, d = 0.004754179283, x = 0.01911075812}
   {c = 0.9969396412, d = 0.004651477802, x = 0.01917417862}
   {x = 0.01920863894, c = 0.9967849498, d = 0.004538930760}
   {c = 0.9966268715, d = 0.004415727178, x = 0.01921121496}
   {d = 0.004280661563, c = 0.9964654850, x = 0.01917783850}
   {c = 0.9963008756, d = 0.004132483620, x = 0.01910388564}
   {x = 0.01898319765, c = 0.9961331523, d = 0.003969510716}
   {c = 0.9959624362, d = 0.003789924217, x = 0.01880847708}
   {c = 0.9957888871, d = 0.003591208983, x = 0.01856967465}
   {c = 0.9956126971, d = 0.003370349304, x = 0.01825382208}
   {c = 0.9954341276, d = 0.003123109890, x = 0.01784231289}
   {c = 0.9952535277, d = 0.002843761390, x = 0.01730831803}
   {c = 0.9950714120, x = 0.01660932168, d = 0.002523642439}
   {c = 0.9948885869, d = 0.002148907648, x = 0.01567231444}
   {c = 0.9947065267, d = 0.001693638091, x = 0.01434942639}
   {d = 0.001097123183, x = 0.01224812569, c = 0.9945286044}
      /
      |
fsolve\

   /
{
   \

/            2                   2                         2
\(0.994398 + x)- (0.994398 + x) c- 0.106 (0.994398 + x) cd

                                  2\ 2
   - 0.9944139104 - 1.000016 x + c /^

MathPen001 发表于 2008-4-2 16:00

确实复杂,看得有点晕。
页: [1]
查看完整版本: 请教方程组求解