FUKA
|
Here lies the binary neutron star initial data solver and diagnostic codes. The initial data is constructed using maximal slicing with a flat spacial metric with which it is possible to achieve maximal spins up to approximately [-0.6, 0.6]
. The v2 code is a considerable improvement over the v1 code as it uses superimposed NS solutions to initialize the binary instead of building the binary from equal mass, non-boosted TOV solutions followed by iterative changes to the binary components. When comparing to v1, generating equal-mass non-rotating is automated (i.e. the user no longer needs to perform the setup and solve separately), however, generating the ID is slightly slower due to needing to perform a stage of fixed orbital frequency prior to solving the hydro fields consistently. Most importantly, when generating arbitrary spin and unequal mass, the cost savings is roughly (N)x faster where N is, in the case of v1, the number of iterative solutions needed to achieve a given mass ratio and spin configurations.
Note: When referring to Mtot
below, we will be referring to the sum of the ADM masses of the TOV solution as measured at infinite separation (i.e. isolated TOV solutions) - Mtot := (MADM_MINUS + MADM_PLUS)
where plus and minus simply refer to their location on the x-axis.
CMakeLists.txt
is the file needed by CMake to compile the codescompile
is a symbolic link to the script stored in $HOME_KADATH/build_release
to ease compilingsrc
directory contains the relevant source files:solve.cpp
: the one and done solve codereader.cpp
: the reader can provide diagnostics from ID solutions that are computed from the IDsolve
for the first timempirun ./bin/Release/solve initial_bns.info
solve
for the first timempirun ./bin/Release/solve initial_bns.info
BNS_ECC_RED.togashi.30.2.0.0.2.8.q1.0.0.09.<info/dat>
We can deconstruct the name to make it understandable:
BNS_ECC_RED.
denotes a converged BNS solution after the eccentricity reduction stage is completed which uses 3.5th order PN estimates for the orbital frequency and radial infall velocity. This is meant to distinguish the solution from earlier stages which will be discussed later.This also distinguishes it from checkpoints that can be turned on which are saved to file during each iteration of the solvertogashi
: the leading name of the eosfile30.2
: separation distance in geometric units!0.0.
: In this case, both NSs are not spinning2.8.q1
: the total mass 2.8
and mass ratio 1
0.0.09
: nshells = 0
for ns1
, nshells = 0
for ns2
where the resolution in each domain is 9
collocation points in the radial and theta direction with 8
points in the phi directioninfo
: The info file contains all the steering parameters, values used in creating the domain decomposition, stored fields, stages, settings, and controlsdat
: The dat file is a binary file that contains the numerical space and the variable fieldsThe default configuration for BNS has a total mass of 2.8M and non-spinning. Aside from the diagnostics observed during the solver stage, we can use the reader to verify the ID. This can be done by running:
./bin/Release/reader BNS_ECC_RED.togashi.28.0.0.2.8.q1.0.0.09.info
Which results in the following:
The first two blocks contain information related to the component NSs - the second has been abbreviated since it contains identical information. These details are covered in the NS README. The only additional parameter is the Center_COM
. This is the coordinate center of each object when shifted by the "center-of-mass" of the binary or, more specifically, the location of the axis of rotation for the binary that can approximate a quasi-stationary solution.
The third block contains information specifically related to the binary
res
: resolution in each numerical domain [r,theta,phi]
q
: mass ratio q <= 1
Separation
: coordinate separation d [d/Mtot] (d in km)
Orbital Omega
: orbital frequency of the binaryE_b / mu
: Binding energy normalized by the reduced massCOM<x/y>
: shift in the coordinates to find a helical killing vectorA-COM<x/y/z>
: A numerical calculation of the COM based on the analytical prescription from Osokine+ (REF)Note: The Diff
noted by the ADM mass is the symmetric difference between the ADM and Komar mass.
Using your favorite text editor, you can open up the initial_bns.info
. We will go through the file, but we'll discuss only the details relevant to the BNS case. For details on all the parameters you can read more in the Configurator README.
Notes:
initial_bns.info
. Using old initial data unless for very small changes in chi
is inefficient.In FUKAv2.2 a minimal Config file was introduced such that only the basic fixing parameters most relevant to users are shown. This minimal Config file can be bypassed by running: > solve full
to obtain the full Config file. Although useful for development, there is little advantage to using the full Config.
The above includes parameters that can be fixed by the user as well as parameters that are automated in the background and should not be changed. The parameters for each NS are simply copied from the isolated solution which is discussed in detail in the NS README - the same fixing applies also in the BNS. The relevant parameters to discuss are
res
The resolution shown for the individual compact objects is the highest resolution the isolated dataset will be ran at. This can be important for TOV solutions as the total baryonic mass is sensitive to the resolution. res 11
is the minimum recommended for production runsThe fixing parameters most relevant to the binary are
distance
: this is in geometric units! It is important to pick something reasonable. A general rule for a binary with a few orbits is distance = 10 * Mtot
, however this strongly depends on q
and the spins of component objectsouter_shells
: This allows for additional shells to be placed near the compactified domain. This can be helpful for more accurate quasi-equilibrium ID at lower resolution, but otherwise can be ignored and left to 0
q
: this parameter is computed. Changing it by hand does nothingres
: global resolution of the binary!adot
: (optional) This is the radial infall velocity parameter when performing eccentricity reduction. This will be discussed more in the eccentricity reduction section belowecc_omega
: (optional) This is the fixed orbital velocity parameter used when performing eccentricity reduction. This will be discussed more in the eccentricity reduction section below
Within the full Config or the output solutions, Fields
document the fields that are used in the solver and stored in the dat
file. Changing this has no impact.
checkpoint
: this will result in checkpoints being saved to file during each solving iteration - mainly helpful for high resolution binary ID where walltimes or server failures are a concern prior to a converged solution being obtainedcorot_binary
: the objects are no longer fixed based on chi
and instead provide a corotating ID solutionfixed_lapse
: toggling this control enables a fixed lapse on the horizon - not recommendedsequences
: this toggle is enabled by default and essentially tells the driver routine to start from scratch. If this is enabled when attempting to use a previous solution, the previous fields and numerical space (i.e. the dat
file) is ignoreduse_pn
: toggle whether to always use 3.5PN estimates. It is important to ensure this is off if the user wants to specify their own adot
and global_omega
parameters by hand (e.g. for iterative eccentricity reduction)resolve
: force all implicit compact object solutions to be resolved regardless of an existing previous solutioncentralized_cos
: stores all implicit COs into $HOME_KADATH/COs
solver_max_iterations
: set the number of iterations not to exceedsolver_precision
: Determines what is the maximum precision allowed by the solver that determines whether or not the solution has convergedinitial_resolution
: the resolution to use when solving the binary initial before regridding to a higher resolution
Now that you've generated the simplest case and we have a better understanding of the config file, we can try something more interesting
distance 30.2
res 11
ns1
set:madm 1.18
chi 0
res 11
ns2
set:madm 2.42
chi 0.52
res 11
mpirun ./bin/Release/solve initial_bns.info
This time around we see the iterative chi
increase being done for the primary NS as well as a regrid of the solution to the higher resolution before being imported into the initial binary setup. Overall, the main changes observed are related to the isolated NS solvers (see the NS README for details), but the binary solver itself is consistent when compared to the equal mass case.
In the event you changed the resolution to 11pts, the solver will solve the binary at the initial_resolution
until a solution in hydrostatic equilibrium has been obtained. At this point the solution will be regridded to a higher resolution at which point the solver will go through all the activated stages (e.g. TOTAL
, TOTAL_BC
, and ECC_RED
).
This results in the converged dataset of BNS_ECC_RED.togashi.28.0.52.0.3.6.q0.487603.0.0.11.info/dat
, however, the other implicit solutions have been saved as well:
BNS_TOTAL_FIXED_OMEGA.
: is the initial solution after the import of the two isolated solutions have been solved in the binary space for a fixed COM and orbital frequency. The hydro fields are simply rescaled to enforce the specified baryonic mass, but the fluid is not in hydrostatic equilibriumBNS_TOTAL.
: this is the solution in complete hydrostatic equilibrium with the ADM linear momenta and the orbital frequency being fixed by the force-balance equations for each NSBNS_TOTAL_BC.
: is the iterative solution such that the orbital frequency is set to a constant and the remaining linear momenta are driven to zero by varying the COM. In this stage, the fluid deviates from hydrostatic equilibrium since the orbital frequency is a constantBNS_ECC_RED.
: The final solution is one where the orbital frequency and radial infall velocity is fixed to either 3.5th order PN estimates based on the COM obtained in the TOTAL_BC
stage or the values for adot
and ecc_omega
are used in the case of iterative eccentricity reductionWe can of course verify that the ID matches our expectation using
./bin/Release/reader BNS_ECC_RED.togashi.30.2.0.0.52.3.6.q0.487603.0.0.09.info
To generate the initial setup for a BNS ID, we first need to make some guesses based on the input ADM masses and the coordinate separation
COM
- An estimate of the center-of-mass: this is purely Newtonianglobal_omega
- An estimate of the orbital frequency: 3.5th PN estimateOnce these estimates are computed an interface code is ran which
global_omega
At this point, the binary numerical space and fields are constructed and the isolated solutions are interpolated onto the new grid using the idea of superimposed solutions. Specifically:
decay_limit := w
is chosen such that w = distance / 2
decay_rate = exp(-(r_NS / w)^4)
where r_NS
is the coordinate distance to the respective NSFor example, if we wanted to compute the initial guess for the lapse at a given point x
, it would be
lapse(x) = 1. + decay_rate_NS1 * (lapse_NS1(x) - 1.) + decay_rate_NS2 * (lapse_NS2(x) - 1.)
This is then repeated for all fields in all numerical domains, with the compactified domain set the asymptotic values of the fields, i.e. lapse = psi = 1
, log(h) = shift = 0
.
Once the initial guess has been setup for the BNS, the first solver stage solves the full XCTS system of equations consistently and all using a fixed orbital velocity. As a consequence the matter is simply rescaled to achieve the desired baryonic mass for each NS. The reason for this is two fold:
phi
needs to initialize more accurately to the binary configurationThe output file from this stage, BNS_TOTAL_FIXED_OMEGA.
, can readily be discarded once a solution from a later stage is obtained.
Now with a more stable background, the matter becomes a variable field and is allowed to solve consistently in order to obtain hydrostatic equilibrium in the binary configuration. In this stage, the force balance equations are introduced which provides a quasi-equilibrium estimate of global_omega
and one component of the ADM linear momenta is minimized.
Note: in the event one wants to later increase the resolution or make iterative changes (e.g. make small changes to the MADM, MB, CHI of one or both stars), it can only be done using the solution from this stage, BNS_TOTAL.*<info/dat>
, as the initial starting point. Reusing the solutions from the hydro-rescaling stage will more often than not cause the solution to diverge or lead to unphysical results in the numerical evolution. Therefore, these solutions can be useful to retain. to retain.**
In the event the binary resolution was set to something higher than the sequence_setting initial_resolution
, the automated increase resolution will take place here and resolve the binary in hydro-static equilibrium. This is very important as increasing the resolution and repeating a matter-rescaling stage will result in a very inconsistent description of the fluid as it will include numerical errors from the interpolated solution at lower resolution!
In this stage we return to simply rescaling the matter now that it has obtained a more consistent solution and we fix it to the current global_omega
, however, we now minimize the x and y linear momentum components using the force-balance equations to determine the most accurante COM for the given binary. This is very important prior to the eccentricity reduction stage as a precise estimate of the COM is required for the post-Newtonian estimates.
The output file from this stage, BNS_TOTAL_BC.
, can be discarded when no longer needed.
Using either 3.5PN estimates or, in the case of iterative eccentricity reduction, ecc_omega
and adot
from the config file, a final stage of matter rescaling is performed based on the changes introduced by ecc_omega
and adot
. This is the recommended solution to use for evolutions and it is stored with a filename title of BNS_ECC_RED.*<info/dat>