GRHydro C2P inconsistency with and without magnetic field evolution

[Lorenzo Sala]

I have found an inconsistency in simulating static TOVs with tabulated EoS with GRHydro. These TOV are not magnetized, thus evolving with or without magnetic fields evolution should not change the results. Unfortunately, evolving magnetic fields leads to crash.
I’m importing initial data from Lorene/Nrotstar using my thorn ID_Nrotstar, setting up T and Y_e in beta-equilibrium, with uniform initial specific entropy at 1 kb/baryon; this creates a gradient in temperature, from a central ~17MeV to an atmosphere value of 0.01MeV.
Attached you can find plots of maximum rho vs time: if magnetic field is initially set up to 0 and magnetic field evolution is activated, as you can see in the attached .par, rho gets stuck and the simulation crashes after a while due to NaN in metric. If I don’t evolve magnetic field (thus, using a different C2P), I don’t have this issue. Attached you can find the .par files, as you can see they only differ in two lines concerning Bvec evolution.
Thank you for your attention,
Lorenzo Sala

[Roland Haas]

Would you mind trying out what happens if you keep the MHD parfile but actually change the code (schedule.ccl) to call the non-MHD con2prim ie change

if (CCTK_Equals(GRHydro_eos_type,"General")) 
{

  if(CCTK_Equals(Bvec_evolution_method,"GRHydro"))
  {
    schedule Conservative2PrimitiveM IN HydroBase_Con2Prim AS Con2Prim IF GRHydro::execute_MoL_PostStep
    {
      LANG: Fortran
    } "Convert back to primitive variables (general) - MHD version"

to

if (CCTK_Equals(GRHydro_eos_type,"General")) 
{

  if(CCTK_Equals(Bvec_evolution_method,"GRHydro"))
  {
    schedule Conservative2Primitive IN HydroBase_Con2Prim AS Con2Prim IF GRHydro::execute_MoL_PostStep
    {
      LANG: Fortran
    } "Convert back to primitive variables (general) - non-MHD version" 

please?

There are a small number of other changes that happen when you change from non-MHD to the MHD code, in particular when using an tabulated (nuclear) EOS eg the initial prim2con treatment.

You may also want to

1. provide some 3d dump of the initial primitive variables and parfiles that can use those as initial data (or the version of ID_NrotStar that you used if you feel comfortable sharing the code)
2. send the description of your issue to the [e-mail] mailing list which has a larger audience than the bug tracking mailing list
3. you could also try some of the alternative con2prim methods in the paper by Siegel and Moesta (https://arxiv.org/abs/1712.07538)

[Lorenzo Sala]

Thank you very much for suggesting the paper by Siegel and Moesta, I’m already using it for some code I’m writing, I wanted to point out this issue in GRHydro. Unfortunately, using the non-MHD-C2P instead of the MHD-C2P (modifying schedule.ccl) leads to "EOS Problem in C2P hot!!!" (Con2PrimHot.F90, line 255). I’m going to attach as soon as possible a 3d checkpoint taken after the initial set-up, that can be used with the two parfiles I had already attached.
Thank you for your attention

[Roland Haas]

I had a chat with @pmoesta who is using GRHydro for MHD simulations. Philipp had a look at your parfile and pointed out two issues that I had missed:

1. you need to set GRHydro::transport_constraints = yes otherwise the default is a free B field evolution which is unstable and quite likely very much untested.
2. you should set HydroBase::initial_Bvec = "zero" rather than None (unless some code of your own actually sets Bvec to zero). Otherwise Bvec is not initialized and would retain the poison value that Carpet sets it two when allocating memory for it (or even worse it would be left unitalized which could be any random value if poisoning is not used).

[Lorenzo Sala]

Thank you so much for pointing this out, this has fixed everything! Fixing those two parameters, rho and temperature do not get stuck and the plot perfectly overlaps the non-MHD one.