Spherical shell support

[nilsdeppe]

Below are the notes from the discussion at ICERM about spherical shells on Aug. 9.

• Spherical shells are not allowed to be used in portions of the grid where there may be hydrodynamics or non-smooth solutions. For binary neutron star mergers we only have shells in the wavezone. For black hole-neutron star mergers we can likely have a spherical shell around the BH initially, but will need to switch to cubes later. Nils D. generally expects we will do Z4c for any neutron star simulations anyway so we don't have to deal with FD on cubed spheres.
• Non-conforming blocks, i.e. blocks can have more than one neighbor per logical direction. Larry wants to have a Region which is a collection of blocks that have only 1 neighbor in each logical direction. The blocks in neighboring Regions can have more than one block neighbor per logical direction. Nils D. is not convinced that we should actually introduce this as a class/object into the code, but it's definitely helpful for discussions to have this concept.
• We need the serial parts of the code to work correctly:
  1. Mesh needs to work with 2d bases
  2. Spherical harmonic filtering will be necessary and is non-trivial.
  3. All the linear operators like interpolation, differentiation, filtering, and integration need to work on spherical shells. (Depends on 1 only for deciding when to branch, depends on 2. for the filtering but the other operators can be done.)
  4. Element needs to be made aware of that it might have different neighbors and their topologies. Nils D. thinks that we should probably not change the low-level class like Element and Domain unless strictly necessary since these may propagate a lot of work.
  5. The domain creator needs to be able to create a domain with both cubes and spheres. (Depends on 1 and 3 in terms of what to return, but the decision of how to decompose the domain can be figured out before, including maps, options, etc.)
• The ComputeTimeDerivative and ApplyBoundaryCorrections need to be able to correctly handle cube-shell interfaces.
  • The ComputeTimeDerivative action will:
    • spherical shells send it's 2d surface packaged data to all cube neighbors.
    • cubes send their interpolated packaged data to the shell. This may be zero grid points depending on resolution jumps and spherical harmonic locations. However, we will still need to send a message for bookkeeping reasons.
  • The ApplyBoundaryCorrections action will:
    • spherical shells: will receive points from all of the cube neighbors and need to assemble them into the right locations in the vectors. The Riemann problem/upwinding is solved locally on the shell at the expense of not having conservation.
    • cubes: will receive the full spherical harmonic grid from the shell neighbor. They then interpolate to the face mesh of the element, and then solve the Riemann problem there. This gives up strict conservation.
  • We will need to not use FixedHashMap for neighbor information on spherical shells because they may (will) have many more neighbors than the cube-cube neighbors do.
• What do we need for elliptic problems and do we even need shells for elliptic problems? We should weigh cost and interest of adding spherical shells to the elliptic code with the cost and interest of getting the elliptic code to run on GPUs.