Communication Steps After ConservedToPrimitive Call

[ltlancas]

The Problem

Hey everyone! I am running some radiative SN blast wave test problems at low-ish resolution. In order to better take care of the numerous bad cells that come up in these sorts of calculations, near the inside edge of the shell, I am correcting bad pressures by taking the average of nearby cell pressures rather than just applying a floor.

Currently, this is done over both ghost and active cells. Since, the average is taken only over cells in the MeshBlock this means that this averaging could give different answers for ghost cells of a MeshBlock compared to the active cells they correspond to. Meaning that subsequent fluxes that are calculated will be different. This is just in non-amr, plain hydro, nothing fancy (other than the cooling I suppose).

This post isn't specifically concerned with this but I would also like to apply the cell-averaging in a way that maintains the symmetry, but I think this would require another copy of the hydro variables to be kept. Any thoughts on the easiest way to do that would also be appreciated.

My Initial Solution

To get around this I figured I would restrict the ConservedToPrimitive() call to the active cells and then perform another Boundary communication. Having not worked with the TaskList before, my first attempt didn't work. I first added more TaskIDs here to codify another SEND_HYD, RECV_HYD, and SETB_HYD for dependencies. I called these SEND_HYD2, RECV_HYD2, and SETB_HYD2. I then added

AddTask(SEND_HYD2,CONS2PRIM);
AddTask(RECV_HYD2,CONS2PRIM);
AddTask(SETB_HYD2,RECV_HYD2);

here after copying this in the AddTask() function to handle the new TaskIDs. I tried running this and got a SegFault 🙁

Any Help Appreciated

I know I probably did a bunch of things wrong so any help is appreciated!

[c-white]

A quick question before I have time to look in detail: Did you look at the CLEAR_ALLBND task, which resets communication at the end of the step in preparation for the next step? Could it be that you are trying to use communication channels twice in between necessary resets.

[ltlancas]

I did not look at that and it could be that... I was assuming that if I just put an AddTask() in the right spot in the TaskList constructor and changed the dependencies that it would work out 😬

[ltlancas]

I tried adding another task to Clear All Boundaries but got another SegFault :/

I added:

  AddTask(CLEAR_ALLBND2,CONS2PRIM);
  AddTask(SEND_HYD2,CLEAR_ALLBND2);
  AddTask(RECV_HYD2,CLEAR_ALLBND2);
  AddTask(SETB_HYD2,RECV_HYD2);

  // everything else
  TaskID before_bval = (CONS2PRIM|SETB_HYD2);

here and changed the dependency on the AddTask(PHY_BVAL,before_bval). I also added in definitions for the new task IDs in task)list.hpp as:

const TaskID SEND_HYD2(91);
const TaskID RECV_HYD2(92);
const TaskID SETB_HYD2(93);
const TaskID CLEAR_ALLBND2(94);

and added

   else if (id == CLEAR_ALLBND2) {
    task_list_[ntasks].TaskFunc=
        static_cast<TaskStatus (TaskList::*)(MeshBlock*,int)>
        (&TimeIntegratorTaskList::ClearAllBoundary);
    task_list_[ntasks].lb_time = false;
    task_list_[ntasks].task_name.append("ClearAllBoundary");
 }

else if (id == SEND_HYD2) {
    task_list_[ntasks].TaskFunc=
        static_cast<TaskStatus (TaskList::*)(MeshBlock*,int)>
        (&TimeIntegratorTaskList::SendHydro);
    task_list_[ntasks].lb_time = true;
    task_list_[ntasks].task_name.append("SendHydro");
  } 

else if (id == RECV_HYD2) {
    task_list_[ntasks].TaskFunc=
        static_cast<TaskStatus (TaskList::*)(MeshBlock*,int)>
        (&TimeIntegratorTaskList::ReceiveHydro);
    task_list_[ntasks].lb_time = false;
    task_list_[ntasks].task_name.append("ReceiveHydro");
  } 

and

else if (id == SETB_HYD2) {
    task_list_[ntasks].TaskFunc=
        static_cast<TaskStatus (TaskList::*)(MeshBlock*,int)>
        (&TimeIntegratorTaskList::SetBoundariesHydro);
    task_list_[ntasks].lb_time = true;
    task_list_[ntasks].task_name.append("SetBoundariesHydro");
  }

to the AddTask() function.

[ltlancas]

We decided that the best way around this was to simply do the neighbor averaging in a flag-and-apply manner and then increase the number of ghost zones by one more than is needed for reconstruction. This way, all ghost zones that matter in each MeshBlock have all the information they need to do a neighbor average in exactly the same way it is done in the corresponding active cells. The flag-and-apply strategy also avoids errors (and asymmetries) brought on by the order in which the corrections are done. See the bottom of this page.