nsnstohmns cannot be reproduced using LORENE2

[Ken Hui]

I have been using the Lorentz release of EinsteinToolkit to reproduce the simulation nsnstohmns in the ET gallery.
Compiling with LORENE1, the simulation works well but I encountered some errors from the Carpet thorn with LORENE2:

INFO (NSTracker): Found star at (4.21875,12.0938,0)
    29304   824.175 |  126.7906196 |    0.0001241    0.0019978 |    0.0007481 |    1.1329054 |    0.0008567 |      1182
INFO (CarpetRegrid2): Enforcing grid structure properties, iteration 0
INFO (CarpetRegrid2): Enforcing grid structure properties, iteration 1
WARNING level 1 from host chpc-cn024.rc.cuhk.edu.hk process 0
  in thorn CarpetLib, file /opt/Cactus/configs/sim/build/CarpetLib/dh.cc:161:
  ->
/opt/Cactus/configs/sim/build/CarpetLib/dh.cc:973:
   [ml=0 rl=6 c=0] The following grid structure consistency check failed:
   Synchronisation and boundary prolongation: All points must have been received
   needrecv.empty()

ml=0 rl=6
baseextent=([756,504,756]:[6408,11784,6408]:[4,4,4]/[189,126,189]:[1602,2946,1602]/[1414,2821,1414]/5640296116)

ml=0 rl=6 c=0
extent=([756,5748,756]:[932,5916,872]:[4,4,4]/[189,1437,189]:[233,1479,218]/[45,43,30]/58050)
outer_boundaries=[[1,0,1],[0,0,0]]
processor=0
...
...
...
memoryof(gh)=15396
memoryof(dh)=138360
memoryof(dh.light_boxes)=11304
memoryof(dh.local_boxes)=96952
memoryof(dh.level_boxes)=7328
memoryof(dh.fast_boxes)=16156
#gfs=181
memoryof(gfs)=193316
WARNING level 0 from host chpc-cn024.rc.cuhk.edu.hk process 0
  in thorn CarpetLib, file /opt/Cactus/configs/sim/build/CarpetLib/dh.cc:2105:
  -> The grid structure is inconsistent.  It is impossible to continue.
3 total processes killed (some possibly by mpirun during cleanup)
Stopping:

In previous release(s), only LORENE2 respect the manual input of eos_file_path parameters in Meudon_Bin_NS. But now it appears that the problem has been fixed in by:
Cactus/arrangements/ExternalLibraries/LORENE/dist/eos_table_filepath.patch

I wonder if there is no reason to use LORENE2 over LORENE1 anymore, and if someone could provide a hint to fix the above errors arose in Carpet. I am grateful for any help and comments.

[Ken Hui]

(Updates)
Okay, it turns out that if I amend some of the parameters in Lorene/C++/Include/unites.h in LORENE2 according to those in LORENE1, I can reproduce the nsnstohmns simulation.

But I would still like to know if there is a way to understand how it causes the error in Carpet, and how it can be fixed. Thank you.

[Roland Haas]

The error you see is very late and really related to the grid structure (and really should never happen).

We had some instances of small changes in initial data changing results of the BNS gallery example at very late times (close to merger), though never really in a production setup. This is discussed in ticket #2119 . It is very low resolution so some fragility is expected though this change is a bit unexpected.

Can you check that the results up to the point of failure are “sane”?

Can you provide a “diff” of the changes you made to unites.h, please?

‌

[Ken Hui]

I followed the particular ticket you mentioned to plot the Psi4 and maximum density, but I don’t think they give any information related to the failure of the simulation.

lorentz_lorene1: Lorentz release compiled with LORENE
lorentz_lorene2: … LORENE2 (the one giving the Carpet error)
lorentz_lorene2_amend: … LORENE2, but with some changes in unites.h

The values of the density differ even at the first step between lorentz_lorene1 and lorentz_lorene2, which is caused by differences in some physical constants, so I amended some numbers in unites.h of LORENE2 to match that of LORENE.

--- lorene2_unites.h.ori    2022-03-11 02:08:27.660497092 +0800
+++ lorene2_unites.h.amend  2022-03-11 02:03:28.053007609 +0800
@@ -78,13 +78,13 @@
    * and black holes.\ingroup (unites)
    */
   namespace Unites {
-    const double g_si = 6.6738E-11 ;    ///< Newton gravitational constant [SI]
+    const double g_si = 6.6726E-11 ;    ///< Newton gravitational constant [SI]
     const double c_si = 2.99792458E+8 ;     ///< Velocity of light [m/s]
     const double kB_si = 1.3806488E-23 ; ///< Boltzmann constant [J/K]
     const double rhonuc_si = 1.66E+17 ;     ///< Nuclear density [kg/m3] (arbitrary)
     const double km_si = 1.E+3 ;    ///< One kilometer [m]
-    const double msol_si = 1.9885E+30 ;     ///< Solar mass [kg]
-    const double mev_si = 1.602176565E-13 ;   ///< One MeV [J]
+    const double msol_si = 1.989E+30 ;  ///< Solar mass [kg]
+    const double mev_si = 1.6021892E-13 ;   ///< One MeV [J]

     const double r_unit = 1.e4 ;  ///< Lorene's unit of length = 10 km
     const double v_unit = c_si ; ///< Lorene's unit of velocity = c 

‌

[Roland Haas]

Thanks for the update. The changes in solar mass and and gravitational constant I can understand (and expect). Though how can MeV in Joule change? It’s essentially the electron charge which I would have expected to be very well known.

I agree that the plots are not directly useful. The interesting bit would be that failure happens once the central density rises ie star start to interact.

Another useful plot, since the grid structure is what fails would be a plot of the grid structure. Eg using the plot-bboxes tool in https://bitbucket.org/eschnett/carpet/src/master/Carpet/src/util/plot-bboxes It expects as input one of Carpet’s carpet-grid-structure files:

> gnuplot
$ load "plot-bboxes"

then paste (since it reads from /dev/tty) one of the sections from the grid structure file, eg:

iteration 0
maps 1
0 mglevels 1
0 0 reflevels 3
0 0 0 components 2
0 0 0 0   0 ([0,0,0]:[136,136,68]:[4,4,4]/[0,0,0]:[34,34,17]/[35,35,18]/22050) [[1,1,1],[1,1,0]]
0 0 0 1   1 ([0,0,72]:[136,136,136]:[4,4,4]/[0,0,18]:[34,34,34]/[35,35,17]/20825) [[1,1,0],[1,1,1]]
0 0 1 components 2
0 0 1 0   0 ([36,36,36]:[100,100,68]:[2,2,2]/[18,18,18]:[50,50,34]/[33,33,17]/18513) [[0,0,0],[0,0,0]]
0 0 1 1   1 ([36,36,70]:[100,100,100]:[2,2,2]/[18,18,35]:[50,50,50]/[33,33,16]/17424) [[0,0,0],[0,0,0]]
0 0 2 components 2
0 0 2 0   0 ([52,52,52]:[84,84,68]:[1,1,1]/[52,52,52]:[84,84,68]/[33,33,17]/18513) [[0,0,0],[0,0,0]]
0 0 2 1   1 ([52,52,69]:[84,84,84]:[1,1,1]/[52,52,69]:[84,84,84]/[33,33,16]/17424) [[0,0,0],[0,0,0]]

which should produce this plot:

This may also be possible with a dedicated Cactus analysis package like SimulationTools or Kuibit or PyCactus/PostCactus listed on https://docs.einsteintoolkit.org/et-docs/Analysis_and_post-processing .

‌

[Ken Hui]

The definition of those 3 physical constants do not match between LORENE1 and LORENE2. The initial data of the gallery simulation nsnstohmns is (I supposed) taken from https://lorene.obspm.fr/data/binNS.html, which should have been generated with an older LORENE version. That’s why I tried to change the physical constants in LORENE2, to recover their old values.

--- lorene1_unites.h.ori    2022-03-11 13:14:02.810753230 +0800
+++ lorene2_unites.h.ori    2022-03-11 02:08:27.660497092 +0800
@@ -26,8 +26,20 @@


 /*
- * $Id: unites.h,v 1.4 2004/12/01 12:28:32 p_grandclement Exp $
+ * $Id: unites.h,v 1.8 2015/03/17 14:20:00 j_novak Exp $
  * $Log: unites.h,v $
+ * Revision 1.8  2015/03/17 14:20:00  j_novak
+ * New class Hot_eos to deal with temperature-dependent EOSs.
+ *
+ * Revision 1.7  2014/10/13 08:52:37  j_novak
+ * Lorene classes and functions now belong to the namespace Lorene.
+ *
+ * Revision 1.6  2014/06/30 14:34:57  j_novak
+ * Update of the values of some constants (G, M_sol and MeV).
+ *
+ * Revision 1.5  2014/05/13 10:06:12  j_novak
+ * Change of magnetic units, to make the Lorene unit system coherent. Magnetic field is now expressed in Lorene units. Improvement on the comments on units.
+ *
  * Revision 1.4  2004/12/01 12:28:32  p_grandclement
  * Include math.h in unite.h
  *
@@ -47,50 +59,55 @@
  * Initial revision
  *
  *
- * $Header: /cvsroot/Lorene/C++/Include/unites.h,v 1.4 2004/12/01 12:28:32 p_grandclement Exp $
+ * $Header: /cvsroot/Lorene/C++/Include/unites.h,v 1.8 2015/03/17 14:20:00 j_novak Exp $

Indeed I have tried to update the values of G, M_sol and MeV in LORENE1 according to the 2014 revision, recompile ET and run the nsnstohmns simulations. The Carpet errors came up again at the exact same time step.

For the plot of the grids, I am having problem visualising it in gnuplot. Allow me to attach the raw output here for the moment, I will also try to see if I can plot them using other methods later.
https://drive.google.com/drive/folders/1iutL_5avKqVTQxDCvqb8dIY5QwB6A_PN?usp=sharing

‌

[Roland Haas]

So based on your tests there seems to be good evidence that indeed the change in constants in unites.h is responsible for the failure (since you can make LORENE1 runs fail but only changing unites.h).

Unfortunately, I doubt there is much the ET can do about this since the changes in the constants are all on the <1% level so one would not have expected such a large effect. Other than a warning to use “compatible” versions of LORENE when reading in data files, which in itself is not terribly helpful. If LORENE encodes a version number in its output files, then the ET can try and check the number and warn if the version number in the file disagrees with the runtime version. It will generate a lot of false positive “version does not match, be careful” warnings so most likely will be ignored by most researchers.

For the gallery example, it may quite likely be possible to make things work by using a slightly higher resolution (the one the gallery is very low, and really not science-worthy).