Hi all,
I was debugging some simulation stuff and was comparing the DivV field to
that calculated within Enzo for a 2D simulation, and found them to be
different. I noticed that in the definition of DivV the dimensionality of
the data doesn't play into the calculation, and instead it assumes that the
contribution from the z dimension will be zero. Anyways, if I comment out
the z contribution, it matches the calculation from within Enzo.
def _DivV(field, data):
# We need to set up stencils
if data.pf["HydroMethod"] == 2:
sl_left = slice(None,-2,None)
sl_right = slice(1,-1,None)
div_fac = 1.0
else:
sl_left = slice(None,-2,None)
sl_right = slice(2,None,None)
div_fac = 2.0
ds = div_fac * data['dx'].flat[0]
f = data["x-velocity"][sl_right,1:-1,1:-1]/ds
f -= data["x-velocity"][sl_left ,1:-1,1:-1]/ds
ds = div_fac * data['dy'].flat[0]
f += data["y-velocity"][1:-1,sl_right,1:-1]/ds
f -= data["y-velocity"][1:-1,sl_left ,1:-1]/ds
If I comment out the following lines, I get matching results to within Enzo.
ds = div_fac * data['dz'].flat[0]
f += data["z-velocity"][1:-1,1:-1,sl_right]/ds
f -= data["z-velocity"][1:-1,1:-1,sl_left ]/ds
Since DivV needs the ValidateSpatial, I'm assuming that there is something
going wrong in the calculation of the ghost zones here. Does anyone have an
idea of what might be going on?
I could see us using a fix that does something like:
if len(data.pf.dimensions) > 2:
ds = div_fac * data['dz'].flat[0]
f += data["z-velocity"][1:-1,1:-1,sl_right]/ds
f -= data["z-velocity"][1:-1,1:-1,sl_left ]/ds
but I'm not sure if that would be the right way to handle it.
Any thoughts?
Sam
