给定速度场下的对流流动

Question

我用FiPy来模拟化学物质在给定速度场中的对流/扩散，我在把这个场设为对流项系数时遇到了问题

我的速度场是用两个二维数组来表示的，比如说Ugrid表示水平速度，Vgrid表示垂直速度，表示每个单元格表面的值。出于这个原因，我认为将这个字段设置为conventionterm的coeff参数的正确方法是将它赋给FaceVariable

但是，我不知道如何将这两个数组作为FaceVariable的value传递。显然，我在使用value = [np.ravel(Ugrid),np.ravel(Vgrid)]将字段设置为CellVariable方面没有问题，对流模拟似乎也有意义，但我不认为这是正确的，正如我在上面简要提到的那样

有什么建议吗

Answer 1

FiPy的一个重要考虑因素是网格索引不是基于网格索引的，因此从网格数组到单元值的映射不应假定网格是以任何特定方式排序的。这就需要使用某种形式的插值。让我们首先构造一个具有网格索引的速度场，其中点位于网格点（而不是单元中心）上

from fipy import Grid2D, CellVariable, FaceVariable
from fipy import ConvectionTerm, DiffusionTerm, TransientTerm
import numpy as np
from scipy import interpolate

dx = 0.5
nx = 7
dy = 0.5
ny = 7

xy = np.zeros((nx, ny, 2))
xy[:, :,  0] = np.mgrid[0:nx, 0:ny][0] * dx
xy[:, :, 1] = np.mgrid[0:nx, 0:ny][1] * dy

u = xy[..., 0] + xy[..., 1]
v = xy[..., 0] * xy[..., 1]

我们现在有一个(u, v)速度场，每个形状nx, ny和相应的坐标，xy，形状(nx, ny, 2)。假设我们想把它插值到一个具有相同域但不同网格的网格上

m = Grid2D(nx=3, ny=3, dx=1.0, dy=1.0)

网格不一定需要与速度场的网格点对齐。然后我们可以用

xy_interp = np.array(m.faceCenters).swapaxes(0, 1)
u_interp = interpolate.griddata(xy.reshape(-1, 2), u.flatten(), xy_interp, method='cubic')
v_interp = interpolate.griddata(xy.reshape(-1, 2), v.flatten(), xy_interp, method='cubic')

其中xy_interp是网格的面中心。请注意，使用griddata要求xy_interp在xy内，否则会给出nan值。一旦我们有了插值，我们就可以建立FiPy速度场

velocity = FaceVariable(mesh=m, rank=1)
velocity[0, :] = u_interp
velocity[1, :] = v_interp

注意，ConvectionTerm的系数可以是FaceVariable或CellVariable。一旦我们有了速度，我们就可以建立和求解方程了

var = CellVariable(mesh=m)
eqn = (TransientTerm() + ConvectionTerm(velocity) == DiffusionTerm(1.0))
eqn.solve(var, dt=1.0)

这对我来说没有错误

Link to complete script for this problem