1

在评估矩阵的符号行列式时,我在 evalf 函数中发现了一个严重但无法解释的错误,请参见下面的 MWE

import sympy
z1,z2,z3,mv1,mv2,mv3,R,T,D12,D23,D13,D14,D24,D34,x1,x2,x3 = sympy.symbols('z1,z2,z3,mv1,mv2,mv3,R,T,D12,D23,D13,D14,D24,D34,x1,x2,x3')
matrix = sympy.Matrix([ [  x2/D12 +  x3/D13 + 1/D14 ,         -x2/D12          ,           -x3/D13        ,  mv1/R/T  ] ,
                      [             -x1/D12         ,  x1/D12 + x3/D23 + 1/D24 ,           -x3/D23        ,  mv2/R/T  ] ,
                      [             -x1/D13         ,         -x2/D23          ,  x1/D13 + x2/D23 + 1/D34 ,  mv3/R/T  ] ,
                      [               z1*x1         ,           z2*x2          ,            z3*x3         ,     0     ]  ])

det1 = sympy.det(matrix.evalf(subs={ z1:-1,z2:-1,z3:1,mv1:50E-6,mv2:90E-6,mv3:150E-6,R:8.314,T:298,D12:1E-9,D23:9E-9,D13:3E-9,D14:0.165171466666667,D24:0.0917619259259259,D34:0.0550571555555556,x1:0.3,x2:0.2,x3:0.5}))
det2 = sympy.det(matrix).evalf(subs={z1:-1,z2:-1,z3:1,mv1:50E-6,mv2:90E-6,mv3:150E-6,R:8.314,T:298,D12:1E-9,D23:9E-9,D13:3E-9,D14:0.165171466666667,D24:0.0917619259259259,D34:0.0550571555555556,x1:0.3,x2:0.2,x3:0.5})
det3 = sympy.det(matrix).subs({      z1:-1,z2:-1,z3:1,mv1:50E-6,mv2:90E-6,mv3:150E-6,R:8.314,T:298,D12:1E-9,D23:9E-9,D13:3E-9,D14:0.165171466666667,D24:0.0917619259259259,D34:0.0550571555555556,x1:0.3,x2:0.2,x3:0.5})
print(det1)
print(det2)
print(det3)

产生以下答案:

det1:3.05615930451006e-7

det2:-439498375.118201

det3: 0

显示 sympy.det() 的符号输出上的 evalf() 有一个严重的错误——所有三个答案都应该接近于零(在数字误差的 ~1E-7 范围内)。方法 1 和 2 工作正常,但显然 det2 方法有一个巨大的错误。对 evalf 使用更高的精度并没有帮助。

感谢任何人的帮助:有人知道我是否在这里犯了错误,或者 sympy 这里有一个大错误吗?

4

0 回答 0