python - FiPy 方程组发散。建议？

Question

我试图简化我的 PDE 集：

在这个新集合中，p 和 m 是我的变量，其中所有其他值都是常数。fc²m|m|/2DA²p 中的 m 和 p 值被视为每个单元格的平均值。此外，我想固定入口质量流量 (m(0, t)) 和出口压力 (p(L,t))，因此我将其定义为约束。

我有两个问题：我的代码似乎发散了，对于每个时间步，它要么给我一个非常高的 m 在管道中的值，而下一个给我一个非常负的值，依此类推。此外，似乎简单地定义入口质量流量和出口压力给出了一个不真实的解决方案。

有什么想法和建议吗？再次感谢您，亲爱的社区。

下面的代码：

import fipy as fi
import matplotlib.pyplot as plt
import numpy as np
import scipy as sci
from matplotlib import animation
from fipy.variables.cellVariable import CellVariable

#1. Domain
L = 101
dx = .1
mesh = fi.Grid1D(nx = L, dx=dx)

#2. Parameters values (Arbitrary) 
Lambda = 0.0001 # Friction factor
D = .4 # Pipe diameter
T = 350 # Gas Temperature
c = 380
A = 3.14 * (D **2) / 4

#3. Variables
## Pressure and mass flow
p = fi.CellVariable(name="pressure", 
                    hasOld=True, 
                    mesh=mesh, 
                    value=0.)
p.setValue(10000000.) # This is the initial condition

m = fi.CellVariable(name="mass flow", 
                    hasOld=True, 
                    mesh=mesh, 
                    value=0.)
m.setValue(300.) # This is the initial condition

#4. Boundary conditions
p.constrain (2000000., where = mesh.facesRight)
p.constrain (15., where = mesh.facesLeft)
m.constrain (500., where = mesh.facesLeft)
m.constrain (500., where = mesh.facesRight)

#5. PDE
eq1 = fi.TransientTerm(var=p) == fi.ConvectionTerm(coeff = [-(c**2) / A], var=m)
eq2 =  fi.TransientTerm(var=m) == fi.ConvectionTerm(coeff = [-A], var=p) - (Lambda*(c**2)*m.cellVolumeAverage*abs(m.cellVolumeAverage)/(2*D*A*p.cellVolumeAverage))

eqn = (eq1 & eq2)

timeStepDuration = .001
steps = 50

#This make a chart of the results
for step in range(steps):
    plt.plot(np.linspace(0, L, nx), m.value)
    #plt.xlim(0,1.1)
    #plt.ylim(0, 20)
    plt.show()
    
    eqn.sweep(dt=timeStepDuration)
    p.updateOld()
    m.updateOld()

标签： pythonpdefipy