在 SymPy 中求解一组非线性方程以获得离散值

Question

我正在尝试求解 SymPy 中离散值的非线性方程组。该系统由 7 个方程组成。

这是我的代码：

import math
from CoolProp.CoolProp import PropsSI
import sympy as sym

alpha_Fluid, Re_l, v_Fluid, mdot_Fluid, Nu_Fluid, Nu_l_lam, Nu_l_turb = sym.symbols('alpha_Fluid, Re_l, v_Fluid, mdot_Fluid, Nu_Fluid, Nu_l_lam, Nu_l_turb')

################################ Discrete Parameters ################################
fluid = 'Nitrogen';
p_Fluid = 1.1e5; 
T_Fluid = 80; 
T_ZP = 300; 
Pr_Fluid = PropsSI("PRANDTL",'P|gas', p_Fluid, 'T', (T_Fluid+T_ZP)/2, fluid);
eta_Fluid = PropsSI("VISCOSITY",'P|gas', p_Fluid, 'T', (T_Fluid+T_ZP)/2, fluid);
lambda_Fluid = PropsSI("CONDUCTIVITY",'P|gas', p_Fluid, 'T', (T_Fluid+T_ZP)/2, fluid);
rho_Fluid = PropsSI("DMASS",'P|gas', p_Fluid, 'T', (T_Fluid+T_ZP)/2, fluid);
Bi = 0.1;
r_ZP = 7e-3; 
lambda_ZP = 230; 
L_ZP_Biot = r_ZP; 
L_ZP_KONV = (math.pi/2)*2*r_ZP; 
A_Inflow = 4e-3; 
A_Anstr = 2*math.pi*r_ZP*A_Inflow; 


Eq1 = sym.Eq((alpha_Fluid*L_ZP_Biot)/lambda_ZP, Bi)
Eq2 = sym.Eq((alpha_Fluid*L_ZP_KONV)/lambda_Fluid, Nu_Fluid) 
Eq3 = sym.Eq(0.664*(Re_l)**(1/2)*Pr_Fluid**(1/3),Nu_l_lam)
Eq4 = sym.Eq((0.037*Re_l**(0.8)*Pr_Fluid)/(1+2.443*Re_l**(-0.1)*(Pr_Fluid**(2/3)-1)), Nu_l_turb)
Eq5 = sym.Eq(0.3 + (Nu_l_lam**2 + Nu_l_turb**2)**(1/2), Nu_Fluid)
Eq6 = sym.Eq((rho_Fluid*v_Fluid*L_ZP_KONV)/eta_Fluid, Re_l)
Eq7 = sym.Eq((mdot_Fluid*rho_Fluid)/A_Anstr, v_Fluid)

result = sym.solve([Eq1, Eq2, Eq3, Eq4, Eq5, Eq6, Eq7],mdot_Fluid)
print(result)

我期望变量“mdot_Fluid”有一个离散值。

但是我得到的解决方案如下：{mdot_Fluid: 8.99341199537194e-5*v_Fluid}。

我的第一个猜测是，整个系统是不确定的，这很奇怪，因为我有 7 个方程和 7 个未知数。

Answer 1

这可以解为 6 个线性方程和第 7 个非线性方程。求解除方程 4 之外的所有方程（不包括 propsSI 变量 - 使用 Dummy() ），然后将解代入方程 4。当您知道 propsSI 的值时，将它们代入方程 4 并使用

nsolve

进行初始猜测

v_Fluid

得到答案。

import math
import sympy as sym
from sympy import *
math=sym  # don't use math.pi use sympy.pi
alpha_Fluid, Re_l, v_Fluid, mdot_Fluid, Nu_Fluid, Nu_l_lam, Nu_l_turb = sym.symbols('alpha_Fluid, Re_l, v_Fluid, mdot_Fluid, Nu_Fluid, Nu_l_lam, Nu_l_turb')

################################ Discrete Parameters ################################
fluid = 'Nitrogen';
p_Fluid = 1.1e5; 
T_Fluid = 80; 
T_ZP = 300; 
Pr_Fluid = Dummy('Pr_Fluid')#(T_Fluid+T_ZP)/2
eta_Fluid = Dummy('eta_Fluid')#(T_Fluid+T_ZP)/2
lambda_Fluid =Dummy('lambda_Fluid')#(T_Fluid+T_ZP)/2
rho_Fluid = Dummy('rho_Fluid')#(T_Fluid+T_ZP)/2
Bi = 0.1;
r_ZP = 7e-3; 
lambda_ZP = 230; 
L_ZP_Biot = r_ZP; 
L_ZP_KONV = (math.pi/2)*2*r_ZP; 
A_Inflow = 4e-3; 
A_Anstr = 2*math.pi*r_ZP*A_Inflow; 


Eq1 = sym.Eq((alpha_Fluid*L_ZP_Biot)/lambda_ZP, Bi)
Eq2 = sym.Eq((alpha_Fluid*L_ZP_KONV)/lambda_Fluid, Nu_Fluid) 
Eq3 = sym.Eq(0.664*(Re_l)**(1/2)*Pr_Fluid**(1/3),Nu_l_lam)
Eq4 = sym.Eq((0.037*Re_l**(0.8)*Pr_Fluid)/(1+2.443*Re_l**(-0.1)*(Pr_Fluid**(2/3)-1)), Nu_l_turb)
Eq5 = sym.Eq(0.3 + (Nu_l_lam**2 + Nu_l_turb**2)**(1/2), Nu_Fluid)
Eq6 = sym.Eq((rho_Fluid*v_Fluid*L_ZP_KONV)/eta_Fluid, Re_l)
Eq7 = sym.Eq((mdot_Fluid*rho_Fluid)/A_Anstr, v_Fluid)
req = nsimplify(Tuple(*[Eq1, Eq2, Eq3, Eq5, Eq6, Eq7]))
props = [Pr_Fluid,eta_Fluid,lambda_Fluid,rho_Fluid]
result = sym.solve(req, dict=True, exclude=props)[0]
nl = Eq4.subs(result)
reps = dict(zip(props, [.1,.2,.3,.4])) # or whatever values are pertinent
nsolve(nl.subs(reps), guess)

或者，对

nl

两侧求平方，将

v_Fluid**.1

替换为正值

x

，并在已知

real_roots

的值后对结果使用

props

。下面我展示了正在使用的 .1,.2,.3,.4 的值

>>> eq=nsimplify(nl.lhs**2-nl.rhs**2).as_numer_denom()[0].subs(root(v_Fluid,10),var('x',positive=True)).n()
>>> [i.n(3) for i in real_roots(ex.subs(dict(zip(reps,[.1,.2,.3,.4]))))]
[-5.10, 2.61, 2.63, 4.77]

您必须将这些值替换为

nl

才能查看哪些值（如果有）实际上是解决方案。