一组发散的方程组。建议?

2024-10-04 03:24:38 发布

您现在位置:Python中文网/ 问答频道 /正文

我试图简化我的PDE集:Set of equations.

在这个新的集合中,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()

Tags: theimportvaluevaraspltthiswhere