1

我正在尝试解决非线性方程组。问题是我的输入值作为有效解决方案返回给我。为了实现这一点，必须忽略一些方程。我都试过了sympy.solvers.nsolve()和scipy.optimize.fsolve()。两者都给出相同的答案。我已经发布了以 scipy 和 sympy 开头的两个代码。最后，我发布了示例结果。

from scipy.optimize import fsolve  # required library
import sympy as sp
import scipy
# known values hard coded for testing
a = 1
b = 3
c = a + b
d = 0
kp_NO = 0.00051621
kp_H2O = 0.0000000127
kp_CO2 = 0.00000001733
p = 5 # pressure
dec = 3 # decimal point precision

# Solving the system of equations
def equations(vars):
    (e, f, g, h, j, k, l) = vars
    f1 = e + j - a
    f2 = e + 2*g + 2*j + k + l - a - 2*c
    f3 = f + k - b
    f4 = 2*h + l - 2*d - (2*79/21)*c
    f5 = kp_NO - (l/(c*(79/21)*c))
    f6 = kp_H2O - (k/(b*sp.sqrt((c*p)/(e + f+ g + h + j + k + l))))
    f7 = kp_CO2 - (j/(a*sp.sqrt((c*p)/(e + f + g + h + j + k + l))))
    return[f1, f2, f3, f4, f5, f6, f7]
e, f, g, h, j, k, l = scipy.optimize.root_scalar(equations, (0, 0, 0, 0, 0, 0, 0))
# CO, H2, O2, N2, CO2, H2O, NO
print(e, f, g, h, j, k, l)

import sympy as sp # required library
import math as m
# known values hard coded for testing
a = 1
b = 3
c = a + b
d = 0
kp_NO = 0.0000000015346
kp_H2O = 1.308*10**-23
kp_CO2 = 9.499*10**-46
p = 5 # pressure
dec = 10 # decimal point precision

# Solving the system of equations
e, f, g, h, j, k, l = sp.symbols('e, f, g, h, j, k, l')
f1 = e + j - a
f2 = e + 2*g + 2*j + k + l - a - 2*c
f3 = f + k - b
f4 = 2*h + l - 2*d - (2*79/21)*c
f5 = kp_NO - (l / (c * (79 / 21) * c))
f6 = kp_H2O - (k / (b * sp.sqrt((c * p) / (e + f + g + h + j + k + l))))
f7 = kp_CO2 - (j / (a * sp.sqrt((c * p) / (e + f + g + h + j + k + l))))
# f5 = m.log(kp_NO, p) + h/2 + g/2 - l
# f6 = m.log(kp_H2O, p) + g/2 + f - k
# f7 = m.log(kp_CO2, p) + e + g/2 - j
results = sp.solvers.nsolve([f1, f2, f3, f4, f5, f6, f7], [e, f, g, h, j, k, l], [0.00004, 0.00004, 0.49, 3.76, 0.25, 0.75, 0.01], manual=True)

e = round(results[0], dec)
f = round(results[1], dec)
g = round(results[2], dec)
h = round(results[3], dec)
j = round(results[4], dec)
k = round(results[5], dec)
l = round(results[6], dec)
# CO, H2, O2, N2, CO2, H2O, NO
print(e, f, g, h, j, k, l)

1.00000000000000 3.00000000000000 3.9999999538 15.0476190014 0.0 0.0 9.24e-8
4

1 回答 1

1

如果每个方程都在求解器的容差范围内，您也不会从最初的猜测中得到任何变化。随着方程的缩放，我怀疑情况就是这样。这是另一种方法：

让我们使用有理数而不是小数：

>>> from sympy import nsimplify
>>> eqs = [nsimplify(i, rational=True) for i in [f1, f2, f3, f4, f5, f6, f7]]
>>> v = [e, f, g, h, j, k, l]

除了最后一个方程之外的所有方程都很容易求解，但g; 有1个解决方案：

 >>> lpart = solve(eqs[:-1], set(v) - {g}, dict=True)[0]

当从最后一个方程的分子中去除根式时，它会产生一个可以找到根的三次方

>>> from sympy import real_roots
>>> from sympy.solvers.solvers import unrad
>>> u, cov = unrad(eqs[-1].subs(lpart).as_numer_denom()[0]); assert cov == []
>>> gsol = list(real_roots(u))

看起来前两个实根是最后一个方程的解：

>>> [abs(eqs[-1].simplify().subs(g, i).n()) for i in gsol]
[0.e-145, 0.e-145, 7.84800000000000e-23]

您可以检查它们是否也是其他解决方案。（真正的根不是根式的紧解。）

您可以将常量保留为符号并重复上述操作，并且仅在需要时才用值求解最后一个方程 - 如果需要。

于 2021-12-06T02:17:18.193 回答