1. The Enzo Project
  2. Untitled project
  3. enzo-dev
  4. Pull requests

Pull requests

#176 Merged
Repository
jsoishi
Branch
week-of-code
Repository
enzo
Branch
week-of-code

Fixing a bug in CT MHD outflow boundary conditions

Author
  1. J. S. Oishi
Reviewers
Description

This fixes the outflow boundary condition for the magnetic field in MHDCT. Previously, DivB != 0 on any boundary set to outflow. This rectifies this. A good test is Orzsag-Tang with outflow instead of periodic boundary conditions.

As David points out in the comments below, the problem was that the boundary was actually filling an active zone on the grid, rather than the ghost cells. The confusion arises because of the Add[3], which is necessary to handle the fact that enzo has one more active zone only in the direction of the magnetic field, in order to handle the face centering. So, for Bx along an x boundary, Add = {1, 0, 0}, for By along a y boundary, Add = {0, 1, 0}, and for Bz along a z boundary, Add = {0,0,1}. Prior to this fix, the limits of the loop as well as the selection of the cell to copy over were incorrect. For reference, in enzo, outflow boundary conditions enforce a Neumann boundary condition (in this case, zero gradient) by simply copying the last active cell into all ghost zone cells.

Comments (5)

  1. dcollins4096

    Prior tests were not sufficiently complex to show the indexing error. I had only ever used this with 1d shock tubes, which washes out the indexing error.

  2. Matt Turk

    Hi J. S. Oishi could you update this with some explanatory comments? Specifically, since the presence/absence of Add seems to be a big part of this, having some type of note describing what it means and why it is or is not used in the various loops would be very helpful.

  3. dcollins4096

    For cell centered fields, the first ghost zone on the right boundary is NumberOfGhostZones + NumberOfZones. For the face centered vector field along a given direction, it's one beyond that. The variable "Add" keeps track of this, so for instance when SetMagneticBoundary is called on B_{x,face}, Add={1,0,0} so the boundaries along the x face is extended, but not along y or z. I've added a comment in the code to clarify.