问如何将面板变量转换为单元格变量，将其插入方程中？

EN

from fipy import CellVariable, FaceVariable, Grid2D, DiffusionTerm, Viewer, ExponentialConvectionTerm
from fipy.tools import numerix
from fipy.meshes.gmshMesh import Gmsh2D
import numpy as np

L = 1.0
N = 100
dL = L / N
viscosity = 1
U = 1.

meshsize = 1./N
geo = '''
cl__1 = 1;
Point(1) = {1, 1, 0, %f};
Point(2) = {2, 1, 0, %f};
Point(3) = {2, 2, 0, %f};
Point(4) = {1, 2, 0, %f};
Line(1) = {4, 3};
Line(2) = {3, 2};
Line(3) = {2, 1};
Line(4) = {1, 4};
Line Loop(6) = {1, 2, 3, 4};
Plane Surface(6) = {6};
Transfinite Line {4, 1, 2, 3} = %d Using Progression 1;
Transfinite Surface {6};
Recombine Surface {6};
'''

geo = geo %(meshsize ,meshsize, meshsize, meshsize , 1/meshsize )

mesh = Gmsh2D(geo , background=None)

#0.8 for pressure and 0.5 for velocity are typical relaxation values for SIMPLE
pressureRelaxation = 0.1
velocityRelaxation = 0.1
if __name__ == '__main__':
    sweeps = 100
else:
    sweeps = 5

x = mesh.cellCenters.value[0]
y = mesh.cellCenters.value[1]

sourcevx = CellVariable(mesh=mesh , name = 'sourcevx')
sourcevy = CellVariable(mesh=mesh , name = 'sourcevy')
source_p = FaceVariable(mesh=mesh , rank=1)
source_pcorr = CellVariable(mesh=mesh , name = 'source_pcor')

pressure = CellVariable(mesh=mesh, name='pressure')
pressureCorrection = CellVariable(mesh=mesh)
xVelocity = CellVariable(mesh=mesh, name='X velocity')
yVelocity = CellVariable(mesh=mesh, name='Y velocity')

velocity = FaceVariable(mesh=mesh, rank=1)

sourcevx.setValue(-x**2*numerix.sin(x*y) - 2*x*numerix.cos(x**2 + y**2) - \
            y**2*numerix.sin(x*y))
sourcevy.setValue(-x*y*numerix.sin(x*y) - 2*y*numerix.cos(x**2 + y**2) \
           + numerix.cos(x*y) - y**3*numerix.sin(x*y)/x - \
           3*y**2*numerix.cos(x*y)/x**2 + 6*y*numerix.sin(x*y)/x**3 \
           + 6*numerix.cos(x*y)/x**4)
xVelocityEq = DiffusionTerm(coeff=viscosity) - pressure.grad.dot([1.,0.]) \
               == sourcevx
yVelocityEq = DiffusionTerm(coeff=viscosity) - pressure.grad.dot([0.,1.]) \
               == sourcevy



ap = CellVariable(mesh=mesh, value=1.)
coeff = 1./ ap.arithmeticFaceValue*mesh._faceAreas * mesh._cellDistances

x , y = mesh.faceCenters

source_p.setValue(coeff * (-4*x**2*numerix.sin(x**2 + y**2) \
            - 4*y**2*numerix.sin(x**2 + y**2) + 4*numerix.cos(x**2 + y**2))\
                  - 2*y*numerix.cos(x*y) )

pressureCorrectionEq = DiffusionTerm(coeff=coeff) - velocity.divergence \
                       == source_p