Commit 48820b6e authored by John Douglas Veitch's avatar John Douglas Veitch

Add option for fixed EOS

parent 510f7123
......@@ -737,6 +737,7 @@ LALInferenceModel *LALInferenceInitCBCModel(LALInferenceRunState *state) {
(--tidalT) Enables reparmeterized tidal corrections, only with LALSimulation.\n\
(--4PolyEOS) Enables 4-piece polytropic EOS parmeterization, only with LALSimulation.\n\
(--4SpectralDecomp) Enables 4-coeff. spectral decomposition EOS parmeterization, only with LALSimulation.\n\
(--eos EOS) Fix the neutron star EOS. Use \"--eos help\" for allowed names\n\
(--spinOrder PNorder) Specify twice the PN order (e.g. 5 <==> 2.5PN) of spin effects to use, only for LALSimulation (default: -1 <==> Use all spin effects).\n\
(--tidalOrder PNorder) Specify twice the PN order (e.g. 10 <==> 5PN) of tidal effects to use, only for LALSimulation (default: -1 <==> Use all tidal effects).\n\
(--numreldata FileName) Location of NR data file for NR waveforms (with NR_hdf5 approx).\n\
......@@ -1350,26 +1351,15 @@ LALInferenceModel *LALInferenceInitCBCModel(LALInferenceRunState *state) {
}
// For EOS, must pick to either use tidal, tidalT, or 4-piece polytrope parameters; otherwise throw error message
if(LALInferenceGetProcParamVal(commandLine,"--tidalT")&&LALInferenceGetProcParamVal(commandLine,"--tidal")){
XLALPrintError("Error: cannot use more than one of --tidalT and --tidal and --4PolyEOS and --4SpectralDecomp.\n");
XLAL_ERROR_NULL(XLAL_EINVAL);
} else if(LALInferenceGetProcParamVal(commandLine,"--tidalT")&&LALInferenceGetProcParamVal(commandLine,"--4PolyEOS")){
XLALPrintError("Error: cannot use more than one of --tidalT and --tidal and --4PolyEOS and --4SpectralDecomp.\n");
XLAL_ERROR_NULL(XLAL_EINVAL);
} else if(LALInferenceGetProcParamVal(commandLine,"--tidalT")&&LALInferenceGetProcParamVal(commandLine,"--4SpectralDecomp")){
XLALPrintError("Error: cannot use more than one of --tidalT and --tidal and --4PolyEOS and --4SpectralDecomp.\n");
XLAL_ERROR_NULL(XLAL_EINVAL);
} else if(LALInferenceGetProcParamVal(commandLine,"--tidal")&&LALInferenceGetProcParamVal(commandLine,"--4PolyEOS")){
XLALPrintError("Error: cannot use more than one of --tidalT and --tidal and --4PolyEOS and --4SpectralDecomp.\n");
XLAL_ERROR_NULL(XLAL_EINVAL);
} else if(LALInferenceGetProcParamVal(commandLine,"--tidal")&&LALInferenceGetProcParamVal(commandLine,"--4SpectralDecomp")){
XLALPrintError("Error: cannot use more than one of --tidalT and --tidal and --4PolyEOS and --4SpectralDecomp.\n");
XLAL_ERROR_NULL(XLAL_EINVAL);
} else if(LALInferenceGetProcParamVal(commandLine,"--4PolyEOS")&&LALInferenceGetProcParamVal(commandLine,"--4SpectralDecomp")){
XLALPrintError("Error: cannot use more than one of --tidalT and --tidal and --4PolyEOS and --4SpectralDecomp.\n");
XLAL_ERROR_NULL(XLAL_EINVAL);
} else if(LALInferenceGetProcParamVal(commandLine,"--tidalT")){
if((~~LALInferenceGetProcParamVal(commandLine,"--tidalT") + ~~LALInferenceGetProcParamVal(commandLine,"--tidal")
+~~LALInferenceGetProcParamVal(commandLine,"--4PolyEOS") + ~~LALInferenceGetProcParamVal(commandLine,"--4SpectralDecomp")
+~~LALInferenceGetProcParamVal(commandLine,"--eos")) > 1 )
{
XLALPrintError("Error: cannot use more than one of --tidalT, --tidal, --4PolyEOS, --4SpectralDecomp and --eos.\n");
XLAL_ERROR_NULL(XLAL_EINVAL);
}
if(LALInferenceGetProcParamVal(commandLine,"--tidalT")){
LALInferenceRegisterUniformVariableREAL8(state, model->params, "lambdaT", zero, lambdaTMin, lambdaTMax, LALINFERENCE_PARAM_LINEAR);
LALInferenceRegisterUniformVariableREAL8(state, model->params, "dLambdaT", zero, dLambdaTMin, dLambdaTMax, LALINFERENCE_PARAM_LINEAR);
......@@ -1389,6 +1379,15 @@ LALInferenceModel *LALInferenceInitCBCModel(LALInferenceRunState *state) {
LALInferenceRegisterUniformVariableREAL8(state, model->params, "SDgamma2", zero, SDgamma2Min, SDgamma2Max, LALINFERENCE_PARAM_LINEAR);
LALInferenceRegisterUniformVariableREAL8(state, model->params, "SDgamma3", zero, SDgamma3Min, SDgamma3Max, LALINFERENCE_PARAM_LINEAR);
}
else if((ppt=LALInferenceGetProcParamVal(commandLine,"--eos")))
{
LALSimNeutronStarEOS *eos=NULL;
INT4 errnum=XLAL_SUCCESS;
XLAL_TRY(eos=XLALSimNeutronStarEOSByName(ppt->value), errnum);
if(errnum!=XLAL_SUCCESS)
XLAL_ERROR_NULL(errnum);
LALInferenceAddVariable(model->params, "ns_eos", eos, LALINFERENCE_void_t, LALINFERENCE_PARAM_FIXED);
}
LALSimInspiralSpinOrder spinO = LAL_SIM_INSPIRAL_SPIN_ORDER_ALL;
ppt=LALInferenceGetProcParamVal(commandLine, "--spinOrder");
......
......@@ -804,6 +804,35 @@ void LALInferenceTemplateXLALSimInspiralChooseWaveform(LALInferenceModel *model)
sym_mass_ratio_eta = m1*m2/((m1+m2)*(m1+m2));
LALInferenceLambdaTsEta2Lambdas(lambdaT,dLambdaT,sym_mass_ratio_eta,&lambda1,&lambda2);
}
if(LALInferenceCheckVariable(model->params, "ns_eos"))
{
LALSimNeutronStarFamily *eos_fam = XLALCreateSimNeutronStarFamily((LALSimNeutronStarEOS *)LALInferenceGetVariable(model->params, "ns_eos"));
REAL8 r1=0, r2=0, k2_1=0, k2_2=0;
REAL8 mass_max = XLALSimNeutronStarMaximumMass(eos_fam) / LAL_MSUN_SI;
REAL8 mass_min = XLALSimNeutronStarFamMinimumMass(eos_fam) / LAL_MSUN_SI;
if(m1<mass_max && m1>mass_min)
{
/* Compute l1, l2 from mass and EOS */
r1 = SimNeutronStarRadius(mass1*LAL_MSUN_SI, eos_fam);
k2_1 = SimNeutronStarLoveNumberK2(mass1*LAL_MSUN_SI, eos_fam);
lambda1 = (2./3.)*k2_1 * pow(r1/(mass1*LAL_MRSUN_SI), 5.0);
}
if(m2<mass_max && m2>mass_min)
{
r2 = SimNeutronStarRadius(mass2*LAL_MSUN_SI, eos_fam);
k2_2 = SimNeutronStarLoveNumberK2(mass2*LAL_MSUN_SI, eos_fam);
lambda2 = (2./3.)*k2_2 * pow(r2/(mass2*LAL_MRSUN_SI), 5.0);
}
/* Add derived quantities for output */
LALInferenceAddVariable(model->params, "radius1", &r1, LALINFERENCE_REAL8_t, LALINFERENCE_PARAM_OUTPUT);
LALInferenceAddVariable(model->params, "radius2", &r2, LALINFERENCE_REAL8_t, LALINFERENCE_PARAM_OUTPUT);
LALInferenceAddVariable(model->params, "lambda1", &lambda1, LALINFERENCE_REAL8_t, LALINFERENCE_PARAM_OUTPUT);
LALInferenceAddVariable(model->params, "lambda2", &lambda2, LALINFERENCE_REAL8_t, LALINFERENCE_PARAM_OUTPUT);
/* ClLean up */
if(eos_fam) XLALDestroySimNeutronStarFamily(eos_fam);
}
/* ==== 4-PIECE POLYTROPE EOS PARAMETERS ==== */
REAL8 logp1 = 0.;
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment