usnistgov / fipy

FiPy is a Finite Volume PDE solver written in Python
http://pages.nist.gov/fipy/en/latest
Other
504 stars 148 forks source link

Coefficients in vector equations can't be scalar #919

Open guyer opened 1 year ago

guyer commented 1 year ago

A DiffusionTerm cannot have a scalar diffusion coefficient if the solution variable is a vector. E.g.,

import fipy as fp

mesh = fp.Grid2D(nx=10, ny=10)
vec = fp.CellVariable(mesh=mesh, name="vec", rank=1, hasOld=True)

D = 1

eq = fp.TransientTerm(var=vec) == fp.DiffusionTerm(coeff=D, var=vec)
eq.solve(dt=1)

raises

IndexError: boolean index did not match indexed array along dimension 0; dimension is 180 but corresponding boolean dimension is 720

This can be remedied by changing to

eq = fp.TransientTerm(var=vec) == fp.DiffusionTerm(coeff=[[[D, 0],
                                                           [0, D]]], var=vec)

Note: Similarly, it's necessary to write

fp.ConvectionTerm(coeff=[[[1, 0],
                          [0, 1]], 
                         [[2, 0],
                          [0, 2]]], var=vec)

instead of

fp.ConvectionTerm(coeff=[[1],
                         [2]], var=vec)

Note: While it is acceptable to write

fp.TransientTerm(coeff=[[1, 0],
                        [0, 1]], var=vec)

it is not required. It should also be possible to pass scalar D.