翻译2- d的PDE和条件成准确Fipy码

Question

enter image description here

图像包括控制方程，初始和边界condition.It descripe板和流体之间的热传递的问题。 我不知道如何使用fipy编码包含var的二维问题和边界条件。 这是我的尝试。

from fipy import * 
import numpy as np 
#constant 
Pe=2400 
le_L=1/20000 
L_l=20000 
alphas=1 
alphaf=1 
a=1/Pe+le_L 
b=1/Pe+L_l 
Bi=0.4 
c=Bi/Pe*L_l 
#generate 
mesh=Grid2D(dx=1,dy=1) 
Ts=CellVariable(mesh=mesh,name='Ts',value=900) 
Tf=CellVariable(mesh=mesh,name='Tf',value=300) 
#condition 
Ts.faceGrad.constrain([0.],mesh.facesLeft) 
Ts.faceGrad.constrain([0.],mesh.facesRight) 
Ts.faceGrad.constrain([-1.*Bi*(Tf.value-Ts.value)],mesh.facesBottom) 
Ts.faceGrad.constrain([0.],mesh.facesTop) 

Tf.constrain(300,mesh.facesLeft) 
Tf.grad.constrain(0,mesh.facesRight) 



a=CellVariable(mesh=mesh,rank=1) 
a[:]=1 

#eq 
eq1=TransientTerm(var=Ts)==DiffusionTerm(coeff=[[a,b]],var=Ts) 
eq2=TransientTerm(var=Tf)==DiffusionTerm(coeff=[[a,0]],var=Tf)-  
ExponentialConvectionTerm(a,var=Tf)+ImplicitSourceTerm(c,var=Tf)-   
ImplicitSourceTerm(c,var=Ts) 

eq=eq1&eq2 

#solve 
dt=0.1 
steps=100 
viewer=Viewer(vars=(Ts,Tf),datamax=1000,datamin=0) 
for i in range(steps): 
eq.solve(dt=dt) 
viewer.plot() 

我发现它失败了。我不知道哪里出了问题，我会欢迎任何帮助;非常感谢！ 顺便说一句，我想要得到的最终图像就像enter image description here ......很多thx！

Answer 1

[编辑，以解决一般的边界条件]

下运行，并似乎给自然的结果，你正在寻找：

from fipy import * 
import numpy as np 
#constant 
Pe=2400. 
le_L=1./20000. 
L_l=20000. 
alphasx=alphasy=1. 
alphaf=1. 
Bi=0.4 
c=Bi/Pe*L_l 

Dsxx = alphasx 
Dsyy = alphasy * L_l**2 
Ds = 1./Pe * le_L * (1./alphaf) * Variable([[alphasx, 0.], 
              [0., alphasy * L_l**2]]) 

Df = Variable([[1./Pe * le_L, 0], 
       [0., 0.]]) 

#generate 
mesh=Grid2D(Lx=1.,Ly=1.,nx=100, ny=100) 
Ts=CellVariable(mesh=mesh,name='Ts',value=900.) 
Tf=CellVariable(mesh=mesh,name='Tf',value=900.) 
#condition 
bottom_mask = (mesh.facesBottom * mesh.faceNormals).divergence 
dPR = mesh._cellDistances[mesh.facesBottom.value][0] 
Af = mesh._faceAreas[mesh.facesBottom.value][0] 
bottom_coeff = bottom_mask * Ds[1,1] * Af/(1 + dPR) 

Tf.constrain(300,mesh.facesLeft) 



#eq 
eq1=(TransientTerm(var=Ts)==DiffusionTerm(coeff=Ds,var=Ts) 
    + ImplicitSourceTerm(coeff=bottom_coeff * -Bi, var=Tf) 
    - ImplicitSourceTerm(coeff=bottom_coeff * -Bi, var=Ts)) 
eq2=(TransientTerm(var=Tf)==DiffusionTerm(coeff=Df,var=Tf) 
    -ExponentialConvectionTerm(coeff=[[1.], [0]],var=Tf) 
    +ImplicitSourceTerm(c,var=Tf) 
    -ImplicitSourceTerm(c,var=Ts)) 

eq=eq1&eq2 

#solve 
dt=0.01 
steps=100 
viewer=Viewer(vars=(Ts,Tf),datamax=1000,datamin=0) 
for i in range(steps): 
    eq.solve(dt=dt) 
    viewer.plot() 

• 我改变了一些系数，以同意你提供的数学。
• 我固定的扩散系数有由FiPy预期各向异性扩散的形状
• 我改变了很多整数到浮动的，因为整数不FiPy很好地工作
• 我提供了一个域名，解决了（你的目只有一个单元格，使空间变化不可能）
• 我减少了时间步
• 我介绍了我们知道如何处理general boundary conditions最好的方法。这很丑陋，但我认为它是正确的。

您可能还需要引入扫描来解释方程和边界条件之间的非线性依赖关系。