22 #include <gsl/gsl_const.h>
23 #include <gsl/gsl_errno.h>
24 #include <gsl/gsl_math.h>
25 #include <gsl/gsl_odeiv.h>
27 #include <lal/LALSimInspiral.h>
28 #include <lal/LALConstants.h>
29 #include <lal/LALStdlib.h>
30 #include <lal/TimeSeries.h>
31 #include <lal/Units.h>
36 #define UNUSED __attribute__ ((unused))
49 tagexpnCoeffsEccTaylorT4 {
78 tagexpnFuncEccTaylorT4
112 OTS = sqrt(1. - et * et);
114 CFdvdt = ak->
eta * pow(v, 9.0) / (2. * ak->
mt * ak->
ST);
115 POLdvdt0 = (192. + 584. * pow(et, 2.0) + 74 * pow(et, 4.0)) / (15. * pow(OTS, 7.0));
117 ans = CFdvdt * POLdvdt0;
136 OTS = sqrt(1. - et * et);
138 CFdvdt = ak->
eta * pow(v, 9.0) / (2. * ak->
mt * ak->
ST);
139 POLdvdt0 = (192. + 584. * pow(et, 2.0) + 74 * pow(et, 4.0)) / (15. * pow(OTS, 7.0));
140 POLdvdt2 = ((-11888. + pow(et, 2.0) * (87720. - 159600. * ak->
eta) + pow(et, 4.0) * (171038. - 141708. * ak->
eta)
141 + pow(et, 6.0) * (11717. - 8288. * ak->
eta) - 14784. * ak->
eta) * pow(v, 2.0)) / (420. * pow(OTS, 9.0));
143 ans = CFdvdt * (POLdvdt0 + POLdvdt2);
164 OTS = sqrt(1. - et * et);
166 CFdvdt = ak->
eta * pow(v, 9.0) / (2. * ak->
mt * ak->
ST);
167 POLdvdt0 = (192. + 584. * pow(et, 2.0) + 74 * pow(et, 4.0)) / (15. * pow(OTS, 7.0));
168 POLdvdt2 = ((-11888. + pow(et, 2.0) * (87720. - 159600. * ak->
eta) + pow(et, 4.0) * (171038. - 141708. * ak->
eta)
169 + pow(et, 6.0) * (11717. - 8288. * ak->
eta) - 14784. * ak->
eta) * pow(v, 2.0)) / (420. * pow(OTS, 9.0));
170 FACdvdt3 = (256./5.) *
LAL_PI * pow(v, 3.0);
171 RFTdvdt3 = (1. + 7.260831042 * pow(et, 2.0) + 5.844370473 * pow(et, 4.0) + 0.8452020270 * pow(et, 6.0)
172 + 0.07580633432 * pow(et, 8.0) + 0.002034045037 * pow(et, 10.0)) / (1. - 4.900627291 * pow(et, 2.0)
173 + 9.512155497 * pow(et, 4.0) - 9.051368575 * pow(et, 6.0) + 4.096465525 * pow(et, 8.0)
174 - 0.5933309609 * pow(et, 10.0) - 0.05427399445 * pow(et, 12.0) - 0.009020225634 * pow(et, 14.0));
176 ans = CFdvdt * (POLdvdt0 + POLdvdt2 + (FACdvdt3 * RFTdvdt3));
198 OTS = sqrt(1. - et * et);
200 CFdvdt = ak->
eta * pow(v, 9.0) / (2. * ak->
mt * ak->
ST);
201 POLdvdt0 = (192. + 584. * pow(et, 2.0) + 74 * pow(et, 4.0)) / (15. * pow(OTS, 7.0));
202 POLdvdt2 = ((-11888. + pow(et, 2.0) * (87720. - 159600. * ak->
eta) + pow(et, 4.0) * (171038. - 141708. * ak->
eta)
203 + pow(et, 6.0) * (11717. - 8288. * ak->
eta) - 14784. * ak->
eta) * pow(v, 2.0)) / (420. * pow(OTS, 9.0));
204 FACdvdt3 = (256./5.) *
LAL_PI * pow(v, 3.0);
205 RFTdvdt3 = (1. + 7.260831042 * pow(et, 2.0) + 5.844370473 * pow(et, 4.0) + 0.8452020270 * pow(et, 6.0)
206 + 0.07580633432 * pow(et, 8.0) + 0.002034045037 * pow(et, 10.0)) / (1. - 4.900627291 * pow(et, 2.0)
207 + 9.512155497 * pow(et, 4.0) - 9.051368575 * pow(et, 6.0) + 4.096465525 * pow(et, 8.0)
208 - 0.5933309609 * pow(et, 10.0) - 0.05427399445 * pow(et, 12.0) - 0.009020225634 * pow(et, 14.0));
209 POLdvdt4 = ((-360224. + 4514976. * ak->
eta + 1903104. * pow(ak->
eta, 2.0)
210 + pow(et, 8.0) * (3523113. - 3259980. * ak->
eta + 1964256. * pow(ak->
eta, 2.0))
211 + pow(et, 2.0) * (-92846560. + 15464736. * ak->
eta + 61282032. * pow(ak->
eta, 2.0))
212 + pow(et, 6.0) * (83424402. - 123108426. * ak->
eta + 64828848. * pow(ak->
eta, 2.0))
213 + pow(et, 4.0) * (783768. - 207204264. * ak->
eta + 166506060. * pow(ak->
eta, 2.0))
214 - 3024. * (96. + 4268. * pow(et, 2.0) + 4386. * pow(et, 4.0)
215 + 175. * pow(et, 6.0))*(-5. + 2. * ak->
eta) * OTS) * pow(v, 4.0)) / (45360. * pow(OTS, 11.0));
217 ans = CFdvdt * (POLdvdt0 + POLdvdt2 + (FACdvdt3 * RFTdvdt3) + POLdvdt4);
249 OTS = sqrt(1. - et * et);
251 CFdedt = - et * ak->
eta * pow(v, 8.0) / (ak->
mt * ak->
ST);
252 POLdedt0 = (304. + 121. * pow(et, 2.0)) / (15. * pow(OTS, 5.0));
254 ans = CFdedt * POLdedt0;
273 OTS = sqrt(1. - et * et);
275 CFdedt = - et * ak->
eta * pow(v, 8.0) / (ak->
mt * ak->
ST);
276 POLdedt0 = (304. + 121. * pow(et, 2.0)) / (15. * pow(OTS, 5.0));
277 POLdedt2 = - ((67608. + 228704. * ak->
eta + pow(et, 4.0) * (-125361. + 93184. * ak->
eta)
278 + pow(et, 2.0) * (-718008. + 651252. * ak->
eta)) * pow(v, 2.0)) / (2520. * pow(OTS, 7.0));
280 ans = CFdedt * (POLdedt0 + POLdedt2);
301 OTS = sqrt(1. - et * et);
303 CFdedt = - et * ak->
eta * pow(v, 8.0) / (ak->
mt * ak->
ST);
304 POLdedt0 = (304. + 121. * pow(et, 2.0)) / (15. * pow(OTS, 5.0));
305 POLdedt2 = - ((67608. + 228704. * ak->
eta + pow(et, 4.0) * (-125361. + 93184. * ak->
eta)
306 + pow(et, 2.0) * (-718008. + 651252. * ak->
eta)) * pow(v, 2.0)) / (2520. * pow(OTS, 7.0));
307 FACdedt3 = (394./3.) *
LAL_PI * pow(v, 3.0);
308 RFTdedt3 = 0.1949238579 * OTS * (OTS * (1. + 7.260831042 * pow(et, 2.0) + 5.844370473 * pow(et, 4.0)
309 + 0.8452020270 * pow(et, 6.0) + 0.07580633432 * pow(et, 8.0)
310 + 0.002034045037 * pow(et, 10.0)) / (1. - 4.900627291 * pow(et, 2.0) + 9.512155497 * pow(et, 4.0)
311 - 9.051368575 * pow(et, 6.0) + 4.096465525 * pow(et, 8.0) - 0.5933309609 * pow(et, 10.0)
312 - 0.05427399445 * pow(et, 12.0) - 0.009020225634 * pow(et, 14.0)) - (1. + 1.893242666 * pow(et, 2.0)
313 - 2.708117333 * pow(et, 4.0) + 0.6192474531 * pow(et, 6.0) + 0.0500847462 * pow(et, 8.0)
314 - 0.01059040781 * pow(et, 10.0)) / (1. - 4.638007334 * pow(et, 2.0) + 8.716680569 * pow(et, 4.0)
315 - 8.451197591 * pow(et, 6.0) + 4.435922348 * pow(et, 8.0) - 1.199023304 * pow(et, 10.0)
316 + 0.1398678608 * pow(et, 12.0) - 0.004254544193 * pow(et, 14.0))) / pow(et, 2.0);
318 ans = CFdedt * (POLdedt0 + POLdedt2 + (FACdedt3 * RFTdedt3));
340 OTS = sqrt(1. - et * et);
342 CFdedt = - et * ak->
eta * pow(v, 8.0) / (ak->
mt * ak->
ST);
343 POLdedt0 = (304. + 121. * pow(et, 2.0)) / (15. * pow(OTS, 5.0));
344 POLdedt2 = - ((67608. + 228704. * ak->
eta + pow(et, 4.0) * (-125361. + 93184. * ak->
eta)
345 + pow(et, 2.0) * (-718008. + 651252. * ak->
eta)) * pow(v, 2.0)) / (2520. * pow(OTS, 7.0));
346 FACdedt3 = (394./3.) *
LAL_PI * pow(v, 3.0);
347 RFTdedt3 = 0.1949238579 * OTS * (OTS * (1. + 7.260831042 * pow(et, 2.0) + 5.844370473 * pow(et, 4.0)
348 + 0.8452020270 * pow(et, 6.0) + 0.07580633432 * pow(et, 8.0)
349 + 0.002034045037 * pow(et, 10.0)) / (1. - 4.900627291 * pow(et, 2.0) + 9.512155497 * pow(et, 4.0)
350 - 9.051368575 * pow(et, 6.0) + 4.096465525 * pow(et, 8.0) - 0.5933309609 * pow(et, 10.0)
351 - 0.05427399445 * pow(et, 12.0) - 0.009020225634 * pow(et, 14.0)) - (1. + 1.893242666 * pow(et, 2.0)
352 - 2.708117333 * pow(et, 4.0) + 0.6192474531 * pow(et, 6.0) + 0.0500847462 * pow(et, 8.0)
353 - 0.01059040781 * pow(et, 10.0)) / (1. - 4.638007334 * pow(et, 2.0) + 8.716680569 * pow(et, 4.0)
354 - 8.451197591 * pow(et, 6.0) + 4.435922348 * pow(et, 8.0) - 1.199023304 * pow(et, 10.0)
355 + 0.1398678608 * pow(et, 12.0) - 0.004254544193 * pow(et, 14.0))) / pow(et, 2.0);
356 POLdedt4 = ((-15238352. + 12823920. * ak->
eta + 4548096. * pow(ak->
eta, 2.0)
357 + pow(et, 6.0) * (3786543. - 4344852. * ak->
eta + 2758560. * pow(ak->
eta, 2.0))
358 + pow(et, 4.0) * (46566110. - 78343602. * ak->
eta + 42810096 * pow(ak->
eta, 2.0))
359 + pow(et, 2.0) * (-37367868 - 41949252 * ak->
eta + 48711348 * pow(ak->
eta, 2.0)) - 1008. * (2672.
360 + 6963. * pow(et, 2.0) + 565. * pow(et, 4.0)) * (-5. + 2. * ak->
eta) * OTS) * pow(v, 4.0)) / (30240. * pow(OTS, 9.0));
362 ans = CFdedt * (POLdedt0 + POLdedt2 + (FACdedt3 * RFTdedt3) + POLdedt4);
393 CFdldt = pow(v, 3.0) / (ak->
mt * ak->
ST);
396 ans = CFdldt * POLdldt0;
415 OTS = sqrt(1. - et * et);
417 CFdldt = pow(v, 3.0) / (ak->
mt * ak->
ST);
419 POLdldt2 = - (3. * pow(v, 2.0)) / pow(OTS, 2.0);
421 ans = CFdldt * (POLdldt0 + POLdldt2);
441 OTS = sqrt(1. - et * et);
443 CFdldt = pow(v, 3.0) / (ak->
mt * ak->
ST);
445 POLdldt2 = - (3. * pow(v, 2.0)) / pow(OTS, 2.0);
446 POLdldt4 = ((-18. + 28. * ak->
eta + pow(et, 2.0)*(-51. + 26. * ak->
eta)) * pow(v, 4.0)) / (4. * pow(OTS, 4.0));
448 ans = CFdldt * (POLdldt0 + POLdldt2 + POLdldt4);
478 ydot[0] =
p->funcv(
y[0],
y[1],&
p->ak);
479 ydot[1] =
p->funcet(
y[0],
y[1],&
p->ak);
480 ydot[2] =
p->funcl(
y[0],
y[1],&
p->ak);
481 ydot[3] =
y[0]*
y[0]*
y[0]*
p->ak.av;
505 ak->
m = ak->
m1 + ak->
m2;
506 ak->
mu = m1 * m2 / ak->
m;
507 ak->
nu = ak->
mu / ak->
m;
509 ak->
eta = m1 * m2 / pow(ak->
m, 2.0);
514 ak->
av = 1./(ak->
ST * ak->
mt);
524 XLALPrintError(
"XLAL Error - %s: PN approximant not supported for PN order %d\n", __func__,O);
544 XLALPrintError(
"XLAL Error - %s: PN approximant not yet implemented for PN order %d\n", __func__,O);
548 XLALPrintError(
"XLAL Error - %s: PN approximant not yet implemented for PN order %d\n", __func__,O);
552 XLALPrintError(
"XLAL Error - %s: PN approximant not yet implemented for PN order %d\n", __func__,O);
556 XLALPrintError(
"XLAL Error - %s: PN approximant not supported for PN order %d\n", __func__,O);
560 XLALPrintError(
"XLAL Error - %s: Unknown PN order in switch\n", __func__);
588 len = (*V)->data->length;
592 REAL8 ND,
SIGN, SIGN2,
alpha, beta1, zplus, zminus, SIGN3, z,
s,
w, E0,
u;
593 REAL8 f, f1, f2, f3, f4, u1,
u2,
u3,
u4, xi, sol;
595 for(k = 0; k < len; k++){
597 l = TIME->
data->
data[k]; et = (*ECC)->data->data[k];
600 if(
l > 0.0 ||
l == 0.0){ ND = floor(
l/(2.*
LAL_PI));
602 else{ ND = ceil(
l/(2.*
LAL_PI));
616 alpha = (1. - et)/(4.*et + 1./2.);
617 beta1 = (
l/2.)/(4.*et + 1./2.);
618 zplus = pow(beta1 + sqrt(pow(
alpha, 3.0) + pow(beta1, 2.0)), (1./3.));
619 zminus = pow(beta1 - sqrt(pow(
alpha, 3.0) + pow(beta1, 2.0)), (1./3.));
621 if(beta1 > 0.0 ||beta1 == 0.0){ SIGN3 = 1.; }
624 if(SIGN3 > 0.0){ z = zplus; }
631 w = (
s - (0.078*pow(
s, 5.0))/(1. + et));
634 E0 = (
l + et*(3.*
w - 4.*pow(
w, 3.0)));
638 f =
u - et*sin(
u) -
l;
646 u2 = -f/(f1 + (1./2.)*f2*u1);
647 u3 = -f/(f1 + (1./2.)*f2*
u2 + (1./6.)*f3*(pow(
u2, 2.0)));
648 u4 = -f/(f1 + (1./2.)*f2*
u3 + (1./6.)*f3*(pow(
u3, 2.0)) + (1./24.)*f4*(pow(
u3, 3.0)));
651 if(SIGN2 > 0.0){ sol = xi; }
652 else{ sol = (2.*
LAL_PI - xi); }
654 if(
SIGN < 0.0){ sol = -sol; }
688 len = (*V)->data->length;
692 REAL8 OTS, betaN, vmuN;
695 for(k = 0; k < len; k++){
698 u = U->
data->
data[k]; et = (*ECC)->data->data[k];
699 v = (*V)->data->data[k]; xp = (*TIME)->data->data[k];
703 OTS = sqrt(1. - et*et);
704 betaN = (1. - OTS)/et;
705 vmuN = 2*atan(betaN*sin(
u)/(1. - betaN*cos(
u)));
707 L2PN = pow(v, 4.0)*(
dt*(60. - 24.*ak->
eta)*vmuN
708 + et*(15.*ak->
eta - pow(ak->
eta, 2.0))*OTS*sin(
u))/(8.*
dt*OTS);
741 len = (*V)->data->length;
745 REAL8 beta1PN, beta2PN;
746 REAL8 vmu1PN, vmu2PN;
748 for(k = 0; k < len; k++){
750 v = (*V)->data->data[k]; et = (*ECC)->data->data[k];
751 u = (*U)->data->data[k];
753 dt = 1. - et * cos(
u);
754 OTS = sqrt(1. - et * et);
756 beta1PN = (1. - OTS) / et
758 + ((-4. + pow(et, 2.0) * (8. - 2. * ak->
eta) + ak->
eta + (4. - ak->
eta) * OTS) * pow(v, 2.0)) / (et * OTS);
762 + (-528. - 220. * ak->
eta + 4. * pow(ak->
eta, 2.0) + pow(et, 4.0) * (-3840. + 2086. * ak->
eta
763 - 178. * pow(ak->
eta, 2.0)) + pow(et, 2.0) * (5232. - 1659. * ak->
eta + 177. * pow(ak->
eta, 2.0))
764 + (528. + 220. * ak->
eta - 4. * pow(ak->
eta, 2.0) + pow(et, 2.0) * (288. + 83. * ak->
eta
765 - 41. * pow(ak->
eta, 2.0))) * OTS) * pow(v, 4.0) / (96. * et * pow(OTS, 3.0));
768 vmu1PN = 2. * atan(beta1PN * sin(
u) / (1. - beta1PN * cos(
u)));
769 vmu2PN = 2. * atan(beta2PN * sin(
u) / (1. - beta2PN * cos(
u)));
771 w = et * sin(
u) + vmu2PN
773 + 3. * (et * sin(
u) + vmu1PN) * pow(v, 2.0) / pow(OTS, 2.0)
775 + et * sin(
u) * (4. * pow(
dt, 2.0) * (
dt * (108. + pow(et, 2.0) * (102. - 52. * ak->
eta) - 56. * ak->
eta)
776 - 15. * ak->
eta + pow(ak->
eta, 2.0) + pow(et, 2.0) * (30. * ak->
eta - 2. * pow(ak->
eta, 2.0))
777 + pow(et, 4.0) * (-15. * ak->
eta + pow(ak->
eta, 2.0))) + (4. * ak->
eta - 12. * pow(ak->
eta, 2.0)
778 +
dt * (8. + pow(et, 2.0) * (-8. - 144. * ak->
eta) + 144 * ak->
eta) + pow(et, 4.0) * (4. * ak->
eta
779 - 12. * pow(ak->
eta, 2.0)) + pow(et, 2.0) * (-8. * ak->
eta + 24. * pow(ak->
eta, 2.0))
780 + pow(
dt, 2.0) * (-8. - 148. * ak->
eta + 12. * pow(ak->
eta, 2.0) + pow(et, 2.0) * (- ak->
eta
781 + 3. * pow(ak->
eta, 2.0)))) * OTS) * pow(v, 4.0) / (32. * pow(
dt, 3.0) * pow(OTS, 4.0));
826 const UINT4 blocklen = 1024;
827 const REAL8 visco = 1./sqrt(6.);
841 UINT4 j, len, idxRef = 0;
845 const gsl_odeiv_step_type *T = gsl_odeiv_step_rk4;
847 gsl_odeiv_system sys;
860 if ( !v || !et || !
l || !lambda )
864 y[1] = (*et)->data->data[0] = e_min;
865 y[2] = (*l)->data->data[0] = 0.;
866 y[3] = (*lambda)->data->data[0] = 0.;
870 s = gsl_odeiv_step_alloc(T, 4);
873 gsl_odeiv_step_apply(
s, j*
deltaT,
deltaT,
y, yerr, NULL, NULL, &sys);
875 if (
y[0] > visco ) {
876 XLALPrintInfo(
"XLAL Info - %s: PN inspiral terminated at ISCO\n", __func__);
879 if ( j >= (*v)->data->length ) {
889 (*v)->data->data[j] =
y[0];
890 (*et)->data->data[j] =
y[1];
891 (*l)->data->data[j] =
y[2];
892 (*lambda)->data->data[j] =
y[3];
894 gsl_odeiv_step_free(
s);
913 LEN = (*v)->data->length;
935 memset(
u2->data->data, 0,
u2->data->length
936 *
sizeof(*
u2->data->data));
939 for(k = 0; k < LEN; k++){
940 l1->
data->
data[k] = (*l)->data->data[k];
952 for(k = 0; k < LEN; k++){
953 (*u)->data->data[k] =
u2->data->data[k];
966 memset(
w->data->data, 0,
w->data->length
967 *
sizeof(*
w->data->data));
974 for(k = 0; k < LEN; k++){
975 (*phi)->data->data[k] = (*lambda)->data->data[k] +
w->data->data[k];
981 len = (*phi)->data->length;
984 phiRef -= (*phi)->data->data[len-1];
986 else if( fRef ==
f_min )
987 phiRef -= (*phi)->data->data[0];
996 }
while ((*v)->data->data[j] <= VRef);
997 phiRef -= (*phi)->data->data[idxRef];
999 for (j = 0; j < len; ++j)
1000 (*phi)->data->data[j] += phiRef;
1002 return (
int)(*v)->data->length;
1039 if( fRef != 0. && fRef <
f_min )
1041 XLALPrintError(
"XLAL Error - %s: fRef = %f must be > fStart = %f\n",
1042 __func__, fRef,
f_min);
1047 XLALPrintError(
"XLAL Error - %s: fRef = %f must be < Schwar. ISCO=%f\n",
1048 __func__, fRef, fISCO);
1063 m1, m2,
f_min, fRef, e_min, phaseO);
1069 m1, m2,
r,
i, amplitudeO, phaseO);
1150 #include <lal/PrintFTSeries.h>
1151 #include <lal/PrintFTSeries.h>
1166 XLALSimInspiralEccentricTDPNRestricted(&hplus, &hcross, &tc, phic,
deltaT, m1, m2,
f_min, fRef,
r,
i, e_min, O);
void LALCheckMemoryLeaks(void)
static REAL8 XLALSimInspiralMeanAnomalyEvolution4_0PN(REAL8 v, REAL8 UNUSED et, expnCoeffsEccTaylorT4 *ak)
static int XLALSimInspiralKeplerEquationSolver(REAL8TimeSeries **V, REAL8TimeSeries **ECC, REAL8TimeSeries *TIME, REAL8TimeSeries *U)
static int XLALSimInspiralEccTaylorT4PNEvolveOrbitIntegrand(double UNUSED t, const double y[], double ydot[], void *params)
static REAL8 XLALSimInspiralMeanAnomalyEvolution4_4PN(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
static REAL8 XLALSimInspiralOrbitalVelocityEvolution4_2PN(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
static REAL8 XLALSimInspiralOrbitalVelocityEvolution4_0PN(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
static REAL8 XLALSimInspiralOrbitalVelocityEvolution4_4PN(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
static REAL8 XLALSimInspiralOrbitalEccentricityEvolution4_3PN(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
static REAL8 XLALSimInspiralOrbitalEccentricityEvolution4_2PN(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
static int XLALSimInspiralPhaseContributionW(expnCoeffsEccTaylorT4 *ak, REAL8TimeSeries **V, REAL8TimeSeries **ECC, REAL8TimeSeries **U, REAL8TimeSeries *W)
REAL8() SimInspiralEvolutionEquations4(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
static int XLALSimInspiralKeplerEquationLHS_PN(expnCoeffsEccTaylorT4 *ak, REAL8TimeSeries **V, REAL8TimeSeries *U, REAL8TimeSeries **ECC, REAL8TimeSeries **TIME, REAL8TimeSeries *L2)
static REAL8 XLALSimInspiralOrbitalEccentricityEvolution4_4PN(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
static int XLALSimInspiralEccTaylorT4Setup(expnCoeffsEccTaylorT4 *ak, expnFuncEccTaylorT4 *f, REAL8 m1, REAL8 m2, int O)
static REAL8 XLALSimInspiralOrbitalEccentricityEvolution4_0PN(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
static REAL8 XLALSimInspiralMeanAnomalyEvolution4_2PN(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
static REAL8 XLALSimInspiralOrbitalVelocityEvolution4_3PN(REAL8 v, REAL8 et, expnCoeffsEccTaylorT4 *ak)
int main(int argc, char *argv[])
int XLALSimInspiralPNPolarizationWaveformsEccentric(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, REAL8TimeSeries *V, REAL8TimeSeries *Ecc, REAL8TimeSeries *U, REAL8TimeSeries *Phi, REAL8 m1, REAL8 m2, REAL8 r, REAL8 i, int ampO, int ph_O)
Given time series for a binary's orbital dynamical variables, computes the radial and angular orbital...
int XLALSimInspiralEccentricTDPNRestricted(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, REAL8 phiRef, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, REAL8 fRef, REAL8 r, REAL8 i, REAL8 e_min, int O)
Driver routine to compute the restricted post-Newtonian inspiral waveform.
int XLALSimInspiralEccentricTDPNEvolveOrbit(REAL8TimeSeries **v, REAL8TimeSeries **et, REAL8TimeSeries **l, REAL8TimeSeries **lambda, REAL8TimeSeries **u, REAL8TimeSeries **phi, REAL8 phiRef, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, REAL8 fRef, REAL8 e_min, int O)
Evolves a post-Newtonian orbit using the eccentric Taylor T4 method.
int XLALSimInspiralEccentricTDPN(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, REAL8 phiRef, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, REAL8 fRef, REAL8 r, REAL8 i, REAL8 e_min, int O)
Driver routine to compute the post-Newtonian inspiral waveform.
int XLALSimInspiralEccentricTDPNGenerator(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, REAL8 phiRef, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, REAL8 fRef, REAL8 r, REAL8 i, REAL8 e_min, int amplitudeO, int phaseO)
Driver routine to compute the post-Newtonian inspiral waveform.
void LALDPrintTimeSeries(REAL8TimeSeries *series, const CHAR *filename)
REAL8TimeSeries * XLALResizeREAL8TimeSeries(REAL8TimeSeries *series, int first, size_t length)
REAL8TimeSeries * XLALCreateREAL8TimeSeries(const CHAR *name, const LIGOTimeGPS *epoch, REAL8 f0, REAL8 deltaT, const LALUnit *sampleUnits, size_t length)
void XLALDestroyREAL8TimeSeries(REAL8TimeSeries *series)
const LALUnit lalDimensionlessUnit
int int int XLALPrintInfo(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
int XLALPrintError(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
LIGOTimeGPS * XLALGPSAdd(LIGOTimeGPS *epoch, REAL8 dt)
SimInspiralEvolutionEquations4 * meananom4
SimInspiralEvolutionEquations4 * orbvel4
SimInspiralEvolutionEquations4 * orbecc4