24小时热门版块排行榜

>论坛更新日志 (3541)
>休闲灌水 (197)
>虫友互识 (195)
>文献求助 (176)
>导师招生 (140)
>考博 (103)
>硕博家园 (93)
>基金申请 (72)
>论文投稿 (64)
>考研 (51)
>博后之家 (50)
>公派出国 (32)
>招聘信息布告栏 (31)
>论文道贺祈福 (23)
>教师之家 (23)
>找工作 (20)

返回列表

edisonchou12

捐助贵宾 (初入文坛)

应助: 0 (幼儿园)
金币: 215.5
帖子: 2
在线: 1.4小时
虫号: 2425655
注册: 2013-04-19
专业: 原子和分子物理

[求助] MATLAB求解多元多次方程求得的解竟然不能是原方程为0

程序如下：
NE = 9;
Delta00 = 0.5;
Deltat = 0.2;
steps = 101;

Classicalcorrel00=zeros(1,steps);
Timelist00=zeros(1,steps);

for i=1:steps
Currt=(i-1)*Deltat;
Timelist00(1,i)=Currt;
p=(2/3+(1-NE*(Delta00^2)*(Currt^2)/4)*exp(-NE*(Delta00^2)*(Currt^2)/8)/3);

xx3=-2/3+4*p/3;
T003=1/3-2*p/3;
T006=-2/3+4*p/3;
T007=-2/3;
T306=1/3;
T315=1-2*p;
syms a3 b3 c3 T11 T33 real
Eqn100 = a3 - (c3*T007 - 2*T006*T11 + b3*T306 + T315*T33 + xx3)/(1 + 2*b3^2 + 2*(T11^2) + T33^2);
Eqn101 = b3 - (T006 + a3*T306)/(1 + a3^2);
Eqn102 = c3 - (T003 + a3*T007)/(1 + a3^2);
Eqn103 = T11 - ( T007 - a3*T006)/(1 + a3^2);
Eqn104 = T33 - ( T007 + a3*T315)/(1 + a3^2);
%    [a3 b3 c3 T11 T33]=solve(Eqn100,Eqn101,Eqn102,Eqn103,Eqn104,a3,b3,c3,T11, T33);
[a3 b3 c3 T11 T33]=solve(Eqn100,Eqn101,Eqn102,Eqn103,Eqn104);
Classicalcorrel00(1,i)=((1 + a3^2)*(1 + b3^2 + c3^2 + 2*T11^2 + T33^2) + (26/9)*(1 - 2*p)^2 + ...
      8/9 - 2*(a3*xx3 + (2*T11 + T33 + a3*c3)*T007 + b3*a3*T306 + ...
      c3*T003 + (b3 - 2*a3*T11)*T006 + a3*T33*T315))/8;
end

MATLAB虽然可以很快得出结果，但是把解代入原方程，Eqn100，Eqn101，Eqn102，Eqn103，Eqn104不为0
比如，以i=1为例
MATLAB得到的解为：
K>> [a3 b3 c3 T11 T33]

ans =

[ -0.69723179291582163445262963600762, -0.86312308374866176586915777567782, 0.90787603708517215995857949120792, 0.53134050208298150303610149633743, -0.51450718729075094864284295767401]
代入Eqn100
K>> Eqn100

Eqn100 =

a3 + ((4*T11)/3 + T33 - b3/3 + (2*c3)/3 - 2/3)/(2*T11^2 + T33^2 + 2*b3^2 + 1)

K>> a3 + ((4*T11)/3 + T33 - b3/3 + (2*c3)/3 - 2/3)/(2*T11^2 + T33^2 + 2*b3^2 + 1)

ans =

-0.57062821968242266185604900221681

回复此楼