存储 ImplicitSourceTerm 的值

Question

我正在使用 fipy，我希望在某些选定面上模拟具有自由通量 BC 的流动。按照 other examples，我尝试了 2 种不同的技术：

from fipy import *
from fipy.meshes.nonUniformGrid1D import NonUniformGrid1D as Grid1D
from fipy.tools import numerix as np

dx = [1.,  1. , 0.78584047] 
coeff = 2.
dt = min(0.1, 0.9 * min(dx) ** 2 / (2 * coeff))
steps = 10

mesh = Grid1D(nx = len(dx), dx = dx)

#phi values 
phi = CellVariable(mesh = mesh, value=[10., 0., 0.], hasOld = True)
phi.constrain(0., mesh.facesRight) #dirichlet BC
phi2 = CellVariable(mesh = mesh, value=[10., 0., 0.], hasOld = True)
phi2.constrain(0., mesh.facesRight) #dirichlet BC

#outflow value
phi2Out = CellVariable(mesh = mesh, value=[0., 0., 0.], hasOld = True)

intCoef = FaceVariable(mesh, value=coeff)
extCoef = FaceVariable(mesh, value=coeff)
intCoef.constrain(0., where = mesh.facesRight)
extCoef.constrain(0., where = ~mesh.facesRight)

eq1= TransientTerm() == DiffusionTerm(coeff= intCoef * phi.faceValue) \
  + ImplicitSourceTerm(coeff= (extCoef * phi.faceGrad).divergence)
  
eq2 =TransientTerm() == DiffusionTerm(coeff= intCoef * phi2.faceValue) \
  + phi2 * (extCoef * phi2.faceGrad).divergence
  
eq3 =TransientTerm() == - phi2 * (extCoef * phi2.faceGrad).divergence

phi.updateOld()
phi2.updateOld()
phi2Out.updateOld()
    
for _ in range(steps):
    res =  1e+10
    loop = 0
    while res > 1e-4 and loop < 1000:
        res =  eq1.sweep(var= phi,dt= dt)
        res2 =  eq2.sweep(var=phi2 ,dt= dt)
        res3 =  eq3.sweep(var=phi2Out ,dt= dt)
        res = max(res, res2, res3)
        loop += 1
    if res > 1e-4:
        print('no convergence! res: ', res)
        break
        
    phi.updateOld()
    phi2.updateOld()
    phi2Out.updateOld()
    
    print('\nphi: ', phi * mesh.cellVolumes, sum(phi* mesh.cellVolumes))
    print('phi2: ', phi2* mesh.cellVolumes,sum(phi2* mesh.cellVolumes))
    print('phi2_outFlow: ', phi2Out* mesh.cellVolumes, sum(phi2Out* mesh.cellVolumes))
    print('phi2_total: ', sum(phi2* mesh.cellVolumes + phi2Out* mesh.cellVolumes))

10步后的输出：

phi: [2.21786728 1.82700906 0.63381385] 4.678690190704105

phi2: [2.21786872 1.82701185 0.6338187] 4.678699267591382

phi2_outFlow: [0. 0. 5.32132192] 5.321321918864122

phi2_total: 10.000021186455504

理想情况下，我想提取流出值以检查质量平衡并了解每个单元格的每个项的确切值（其他 sink/source 项将在稍后添加，以及其他自由流动边界）。

我的问题是：

为什么 phi 和 phi2 的值略有不同？
如何在保持 'ImplicitSourceTerm' 的同时提取每个单元格的流出项（当将使用更复杂的网格时），这样效率更高？

感谢您的帮助！

PS：到目前为止，我一直在使用 python3，因此，虽然我对涉及其他求解器的解决方案很感兴趣，但我特别想寻找可以使用 [=44] 实现的东西=].

Answer 1

why are the values of phi and phi2 slightly different?

phi 和 phi2 不同，因为 eq2 收敛速度不如 eq1。这是因为 eq1 比 eq2 更隐含。如果您更改残差的公差，例如 res > 1e-10，您会发现这两个解决方案更加一致。

how could I extract the outflow term for each cell (when a more complex grid will be used) while keeping 'ImplicitSourceTerm', which is more efficient?

您仍然可以评估通量 phi2 * extCoef * phi2.faceGrad，即使您使用 ImplicitSourceTerm.

一般来说，提取每个 Term 在物理上所做的事情并不容易（请参阅问题 #461). You can use the FIPY_DISPLAY_MATRIX 环境变量以了解每个 Term 如何对解决方案矩阵做出贡献，但这可能会也可能不会给你很多关于正在发生的事情的物理直觉。

存储 ImplicitSourceTerm 的值

storing value of ImplicitSourceTerm

fipy