在评估矩阵的符号行列式时,我在 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 这里有一个大错误吗?