Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALSimulation 6.2.0.1-8a6b96f
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
LALSimIMRPhenomHM.h
Go to the documentation of this file.
1#ifndef _LALSIM_IMR_PHENOMHM_H
2#define _LALSIM_IMR_PHENOMHM_H
3
4#include <lal/LALDatatypes.h>
5#include <lal/Sequence.h>
6#include <lal/LALDict.h>
7#include <lal/LALConstants.h>
8#include <lal/XLALError.h>
9#include <lal/FrequencySeries.h>
10#include <math.h>
11
12#ifdef __cplusplus
13extern "C" {
14#endif
15
16#ifdef __GNUC__
17#define UNUSED __attribute__ ((unused))
18#else
19#define UNUSED
20#endif
21
22/**
23 * Dimensionless frequency (Mf) at which define the end of the waveform
24 */
25#define PHENOMHM_DEFAULT_MF_MAX 0.5
26
27/**
28 * eta is the symmetric mass-ratio.
29 * This number corresponds to a mass-ratio of 20
30 * The underlying PhenomD model was calibrated to mass-ratio 18
31 * simulations. We choose mass-ratio 20 as a conservative
32 * limit on where the model should be accurate.
33 */
34#define MAX_ALLOWED_ETA 0.045351
35
36/**
37 * Maximum number of (l,m) mode paris PhenomHM models.
38 * Only count positive 'm' modes.
39 * Used to set default mode array
40 */
41#define NMODES_MAX 6
42
43/**
44 * Highest ell multipole PhenomHM models + 1.
45 * Used to set array sizes
46 */
47#define L_MAX_PLUS_1 5
48
49/* Activates amplitude part of the model */
50#define AmpFlagTrue 1
51#define AmpFlagFalse 0
52
54 LALDict *extraParams);
55
56static int IMRPhenomHM_check_mode_array(LALValue *ModeArray);
57
58/**
59 * useful powers in GW waveforms: 1/6, 1/3, 2/3, 4/3, 5/3, 2, 7/3, 8/3, -1, -1/6, -7/6, -1/3, -2/3, -5/3
60 * calculated using only one invocation of 'pow', the rest are just multiplications and divisions
61 */
62typedef struct tagPhenomHMUsefulPowers
63{
77
78/**
79 * Useful powers of Mf: 1/6, 1/3, 2/3, 4/3, 5/3, 2, 7/3, 8/3, -7/6, -5/6, -1/2, -1/6, 1/2
80 * calculated using only one invocation of 'pow' and one of 'sqrt'.
81 * The rest are just multiplications and divisions. Also including Mf itself in here.
82 */
83typedef struct tagPhenomHMUsefulMfPowers
84{
100
101/**
102 * must be called before the first usage of *p
103 */
105
106/**
107 * must be called before the first usage of *p
108 */
110
111/**
112 * Structure holding Higher Mode Phase pre-computations
113 */
114typedef struct tagHMPhasePreComp
115{
116 double ai;
117 double bi;
118 double am;
119 double bm;
120 double ar;
121 double br;
122 double fi;
123 double fr;
124 double PhDBconst;
125 double PhDCconst;
126 double PhDBAterm;
128
129/**
130 * Structure storing pre-determined quantities
131 * that describe the frequency array
132 * and tells you over which indices will contain non-zero values.
133 */
134typedef struct tagPhenomHMFrequencyBoundsStorage
135{
140 UINT4 freq_is_uniform; /**< If = 1 then assume uniform spaced, If = 0 then assume arbitrarily spaced. */
141 size_t npts; /**< number of points in waveform array */
142 size_t ind_min; /**< start index containing non-zero values */
143 size_t ind_max; /**< end index containing non-zero values */
145
148 REAL8Sequence *freqs,
149 REAL8 Mtot,
150 REAL8 deltaF,
151 REAL8 f_ref_in);
152
154 REAL8Sequence *freqs,
155 REAL8 deltaF);
156
157/**
158 * Structure storing pre-determined quantities
159 * complying to the conventions of the PhenomHM model.
160 * convensions such as m1>=m2
161 */
162typedef struct tagPhenomHMStorage
163{
164 REAL8 m1; /**< mass of larger body in solar masses */
165 REAL8 m2; /**< mass of lighter body in solar masses */
166 REAL8 m1_SI; /**< mass of larger body in kg */
167 REAL8 m2_SI; /**< mass of lighter body in kg */
168 REAL8 Mtot; /**< total mass in solar masses */
169 REAL8 eta; /**< symmetric mass-ratio */
170 REAL8 chi1x; /**< x-component of the dimensionless spin of object 1 w.r.t. Lhat = (0,0,1) */
171 REAL8 chi1y; /**< y-component of the dimensionless spin of object 1 w.r.t. Lhat = (0,0,1) */
172 REAL8 chi1z; /**< z-component of the dimensionless spin of object 1 w.r.t. Lhat = (0,0,1) */
173 REAL8 chi2x; /**< x-component of the dimensionless spin of object 2 w.r.t. Lhat = (0,0,1) */
174 REAL8 chi2y; /**< y-component of the dimensionless spin of object 2 w.r.t. Lhat = (0,0,1) */
175 REAL8 chi2z; /**< z-component of the dimensionless spin of object 2 w.r.t. Lhat = (0,0,1) */
176 REAL8 chip; /**< precession parameter */
182 REAL8 Mf_ref; /**< reference frequnecy in geometric units */
184 UINT4 freq_is_uniform; /**< If = 1 then assume uniform spaced, If = 0 then assume arbitrarily spaced. */
185 size_t npts; /**< number of points in waveform array */
186 size_t ind_min; /**< start index containing non-zero values */
187 size_t ind_max; /**< end index containing non-zero values */
194 REAL8 Rholm[L_MAX_PLUS_1][L_MAX_PLUS_1]; /**< ratio of (2,2) mode to (l,m) mode ringdown frequency */
195 REAL8 Taulm[L_MAX_PLUS_1][L_MAX_PLUS_1]; /**< ratio of (l,m) mode to (2,2) mode damping time */
197
199 PhenomHMStorage *p, /**< [out] PhenomHMStorage struct */
200 const REAL8 m1_SI, /**< mass of companion 1 (kg) */
201 const REAL8 m2_SI, /**< mass of companion 2 (kg) */
202 const REAL8 chi1x, /**< x-component of the dimensionless spin of object 1 w.r.t. Lhat = (0,0,1) */
203 const REAL8 chi1y, /**< y-component of the dimensionless spin of object 1 w.r.t. Lhat = (0,0,1) */
204 const REAL8 chi1z, /**< z-component of the dimensionless spin of object 1 w.r.t. Lhat = (0,0,1) */
205 const REAL8 chi2x, /**< x-component of the dimensionless spin of object 2 w.r.t. Lhat = (0,0,1) */
206 const REAL8 chi2y, /**< y-component of the dimensionless spin of object 2 w.r.t. Lhat = (0,0,1) */
207 const REAL8 chi2z, /**< z-component of the dimensionless spin of object 2 w.r.t. Lhat = (0,0,1) */
208 REAL8Sequence *freqs, /**< Input frequency sequence (Hz) */
209 const REAL8 deltaF, /**< frequency spacing (Hz) */
210 const REAL8 f_ref, /**< reference GW frequency (hz) */
211 const REAL8 phiRef /**< orbital phase at f_ref */
212);
213
214double IMRPhenomHMTrd(
215 REAL8 Mf,
216 REAL8 Mf_RD_22,
217 REAL8 Mf_RD_lm,
218 const INT4 AmpFlag,
219 const INT4 ell,
220 const INT4 mm,
221 PhenomHMStorage *pHM);
222
223double IMRPhenomHMTi(
224 REAL8 Mf,
225 const INT4 mm);
226
228 double *Am,
229 double *Bm,
230 const INT4 mm,
231 REAL8 fi,
232 REAL8 fr,
233 REAL8 Mf_RD_22,
234 REAL8 Mf_RD_lm,
235 const INT4 AmpFlag,
236 const INT4 ell,
237 PhenomHMStorage *pHM);
238
240 REAL8 *a,
241 REAL8 *b,
242 REAL8 flm,
243 REAL8 fi,
244 REAL8 fr,
245 REAL8 Ai,
246 REAL8 Bi,
247 REAL8 Am,
248 REAL8 Bm,
249 REAL8 Ar,
250 REAL8 Br);
251
253 REAL8 *a,
254 REAL8 *b,
255 REAL8 *fi,
256 REAL8 *fr,
257 REAL8 *f1,
258 const REAL8 flm,
259 const INT4 ell,
260 const INT4 mm,
261 PhenomHMStorage *pHM,
262 const int AmpFlag);
263
265 REAL8 Mflm,
266 const INT4 ell,
267 const INT4 mm,
268 PhenomHMStorage *pHM,
269 const int AmpFlag);
270
273 const INT4 ell,
274 const INT4 mm,
275 PhenomHMStorage *pHM,
276 LALDict *extraParams);
277
279 REAL8 fM,
280 INT4 l,
281 INT4 m,
282 REAL8 M1,
283 REAL8 M2,
284 REAL8 X1z,
285 REAL8 X2z);
286
288 COMPLEX16FrequencySeries **hptilde,
289 COMPLEX16FrequencySeries **hctilde,
290 REAL8Sequence *freqs,
291 REAL8 m1_SI,
292 REAL8 m2_SI,
293 REAL8 chi1z,
294 REAL8 chi2z,
295 const REAL8 distance,
296 const REAL8 inclination,
297 const REAL8 phiRef,
298 const REAL8 deltaF,
299 REAL8 f_ref,
300 LALDict *extraParams);
301
304 REAL8Sequence *amps,
305 REAL8Sequence *phases,
306 REAL8Sequence *freqs_geom,
307 PhenomHMStorage *pHM,
308 UINT4 ell,
309 INT4 mm,
310 REAL8 phi0,
311 LALDict *extraParams);
312
314 REAL8Sequence *amps,
315 REAL8Sequence *freqs_geom,
316 PhenomHMStorage *pHM,
317 UINT4 ell,
318 INT4 mm,
319 LALDict *extraParams);
320
322 REAL8Sequence *phases,
323 REAL8Sequence *freqs_geom,
324 PhenomHMStorage *pHM,
325 UINT4 ell,
326 INT4 mm,
327 LALDict *extraParams);
328
330 REAL8 *fringdown,
331 REAL8 *fdamp,
332 UINT4 ell,
333 INT4 mm,
334 REAL8 finalmass,
335 REAL8 finalspin);
336
337#ifdef __cplusplus
338}
339#endif
340
341#endif /* _LALSIM_IMR_PHENOMHM_H */
static double fdamp(double eta, double chi1, double chi2, double finspin)
fdamp is the complex part of the ringdown frequency 1508.07250 figure 9
int IMRPhenomHMPhasePreComp(HMPhasePreComp *q, const INT4 ell, const INT4 mm, PhenomHMStorage *pHM, LALDict *extraParams)
int IMRPhenomHMCore(COMPLEX16FrequencySeries **hptilde, COMPLEX16FrequencySeries **hctilde, REAL8Sequence *freqs, REAL8 m1_SI, REAL8 m2_SI, REAL8 chi1z, REAL8 chi2z, const REAL8 distance, const REAL8 inclination, const REAL8 phiRef, const REAL8 deltaF, REAL8 f_ref, LALDict *extraParams)
int PhenomHM_init_useful_powers(PhenomHMUsefulPowers *p, REAL8 number)
must be called before the first usage of *p
double IMRPhenomHMTi(REAL8 Mf, const INT4 mm)
mathematica function Ti domain mapping function - inspiral
COMPLEX16 IMRPhenomHMOnePointFiveSpinPN(REAL8 fM, INT4 l, INT4 m, REAL8 M1, REAL8 M2, REAL8 X1z, REAL8 X2z)
Define function for FD PN amplitudes.
static int IMRPhenomHM_check_mode_array(LALValue *ModeArray)
static int init_PhenomHM_Storage(PhenomHMStorage *p, const REAL8 m1_SI, const REAL8 m2_SI, const REAL8 chi1x, const REAL8 chi1y, const REAL8 chi1z, const REAL8 chi2x, const REAL8 chi2y, const REAL8 chi2z, REAL8Sequence *freqs, const REAL8 deltaF, const REAL8 f_ref, const REAL8 phiRef)
#define L_MAX_PLUS_1
Highest ell multipole PhenomHM models + 1.
int IMRPhenomHMFreqDomainMapParams(REAL8 *a, REAL8 *b, REAL8 *fi, REAL8 *fr, REAL8 *f1, const REAL8 flm, const INT4 ell, const INT4 mm, PhenomHMStorage *pHM, const int AmpFlag)
helper function for IMRPhenomHMFreqDomainMap
int IMRPhenomHMAmplitude(REAL8Sequence *amps, REAL8Sequence *freqs_geom, PhenomHMStorage *pHM, UINT4 ell, INT4 mm, LALDict *extraParams)
int IMRPhenomHMPhase(REAL8Sequence *phases, REAL8Sequence *freqs_geom, PhenomHMStorage *pHM, UINT4 ell, INT4 mm, LALDict *extraParams)
int init_IMRPhenomHMGet_FrequencyBounds_storage(PhenomHMFrequencyBoundsStorage *p, REAL8Sequence *freqs, REAL8 Mtot, REAL8 deltaF, REAL8 f_ref_in)
derive frequency variables for PhenomHM based on input.
int IMRPhenomHMEvaluateOnehlmMode(COMPLEX16FrequencySeries **hlm, REAL8Sequence *amps, REAL8Sequence *phases, REAL8Sequence *freqs_geom, PhenomHMStorage *pHM, UINT4 ell, INT4 mm, REAL8 phi0, LALDict *extraParams)
double IMRPhenomHMFreqDomainMap(REAL8 Mflm, const INT4 ell, const INT4 mm, PhenomHMStorage *pHM, const int AmpFlag)
IMRPhenomHMFreqDomainMap Input waveform frequency in Geometric units (Mflm) and computes what frequen...
int IMRPhenomHMSlopeAmAndBm(double *Am, double *Bm, const INT4 mm, REAL8 fi, REAL8 fr, REAL8 Mf_RD_22, REAL8 Mf_RD_lm, const INT4 AmpFlag, const INT4 ell, PhenomHMStorage *pHM)
helper function for IMRPhenomHMFreqDomainMap
int IMRPhenomHMMapParams(REAL8 *a, REAL8 *b, REAL8 flm, REAL8 fi, REAL8 fr, REAL8 Ai, REAL8 Bi, REAL8 Am, REAL8 Bm, REAL8 Ar, REAL8 Br)
helper function for IMRPhenomHMFreqDomainMap
int IMRPhenomHMGetRingdownFrequency(REAL8 *fringdown, REAL8 *fdamp, UINT4 ell, INT4 mm, REAL8 finalmass, REAL8 finalspin)
returns the real and imag parts of the complex ringdown frequency for the (l,m) mode.
LALDict * IMRPhenomHM_setup_mode_array(LALDict *extraParams)
read in a LALDict.
int PhenomHM_init_useful_mf_powers(PhenomHMUsefulMfPowers *p, REAL8 number)
must be called before the first usage of *p
UINT4 IMRPhenomHM_is_freq_uniform(REAL8Sequence *freqs, REAL8 deltaF)
helper function to easily check if the input frequency sequence is uniformly space or a user defined ...
double IMRPhenomHMTrd(REAL8 Mf, REAL8 Mf_RD_22, REAL8 Mf_RD_lm, const INT4 AmpFlag, const INT4 ell, const INT4 mm, PhenomHMStorage *pHM)
domain mapping function - ringdown
int l
Definition: bh_qnmode.c:135
double complex COMPLEX16
double REAL8
uint32_t UINT4
int32_t INT4
static const INT4 m
static const INT4 q
static const INT4 a
Structure holding Higher Mode Phase pre-computations.
Structure storing pre-determined quantities that describe the frequency array and tells you over whic...
size_t npts
number of points in waveform array
size_t ind_min
start index containing non-zero values
UINT4 freq_is_uniform
If = 1 then assume uniform spaced, If = 0 then assume arbitrarily spaced.
size_t ind_max
end index containing non-zero values
Structure storing pre-determined quantities complying to the conventions of the PhenomHM model.
REAL8 m1_SI
mass of larger body in kg
REAL8 chi1x
x-component of the dimensionless spin of object 1 w.r.t.
REAL8 chi2y
y-component of the dimensionless spin of object 2 w.r.t.
REAL8 chip
precession parameter
REAL8 eta
symmetric mass-ratio
size_t npts
number of points in waveform array
size_t ind_max
end index containing non-zero values
REAL8Sequence * freqs
REAL8 chi2z
z-component of the dimensionless spin of object 2 w.r.t.
REAL8 m2
mass of lighter body in solar masses
REAL8 chi2x
x-component of the dimensionless spin of object 2 w.r.t.
REAL8 chi1z
z-component of the dimensionless spin of object 1 w.r.t.
REAL8 Mtot
total mass in solar masses
REAL8 m1
mass of larger body in solar masses
REAL8 m2_SI
mass of lighter body in kg
size_t ind_min
start index containing non-zero values
REAL8 Mf_ref
reference frequnecy in geometric units
REAL8 chi1y
y-component of the dimensionless spin of object 1 w.r.t.
UINT4 freq_is_uniform
If = 1 then assume uniform spaced, If = 0 then assume arbitrarily spaced.
Useful powers of Mf: 1/6, 1/3, 2/3, 4/3, 5/3, 2, 7/3, 8/3, -7/6, -5/6, -1/2, -1/6,...
useful powers in GW waveforms: 1/6, 1/3, 2/3, 4/3, 5/3, 2, 7/3, 8/3, -1, -1/6, -7/6,...