25 #define UNUSED __attribute__ ((unused))
33 #include <gsl/gsl_const.h>
34 #include <gsl/gsl_errno.h>
35 #include <gsl/gsl_math.h>
36 #include <gsl/gsl_odeiv.h>
44 static const REAL8 ChlmNewt_ampli[35] = {2.1137745587232057, 6.341323676169617, 0.1412325034218127, 1.7864655618418102, 4.9229202032627635, 0.023872650234580958, 0.2250735048768909, 1.7053495827316825, 4.763908164911493, 0.001122165903318321, 0.06333806741197714, 0.2945348827200268, 1.755276012972272, 5.0817902739730565, 0.00014954736544380072, 0.005296595280255638, 0.10342548105284892, 0.3713362832603404, 1.8983258440274462, 5.727111757630886, 5.184622059790144e-06, 0.0012191691413436815, 0.011783593824950922, 0.14639129995388936, 0.4653654097044924, 2.125638973894669, 6.685178460621457, 5.54485779151375621e-7, 0.0000763473331250837455, 0.00353250998285463003,
45 0.0204988821800401766, 0.19579402814926015, 0.584571015778149663, 2.44207899966355693, 7.9955401278763745};
48 static const REAL8 ChlmNewt_phase[35] = {4.71238898038469, 3.141592653589793, 4.71238898038469, 0.0, 1.5707963267948966, 1.5707963267948966, 0.0, 4.71238898038469, 3.141592653589793, 1.5707963267948966, 3.141592653589793, 4.71238898038469, 0.0, 1.5707963267948966, 4.71238898038469, 3.141592653589793, 1.5707963267948966, 0.0, 4.71238898038469, 3.141592653589793, 4.71238898038469, 0.0, 1.5707963267948966, 3.141592653589793, 4.71238898038469, 0.0, 1.5707963267948966, 4.7123889803846898577, 0.0, 1.5707963267948966192, 3.1415926535897932385,
49 4.7123889803846898577, 0.0, 1.5707963267948966192, 3.1415926535897932385};
52 static const double CNlm[35] = {
53 0.17777777777777778, 6.4,
54 0.0007936507936507937, 0.5079365079365079, 8.678571428571429,
55 2.2675736961451248e-05, 0.008062484252960444, 1.0414285714285714, 14.447971781305114,
56 5.010421677088344e-08, 0.0006384836014465644, 0.03106534090909091, 1.9614216236438458, 25.688197074915823,
57 8.898517797026696e-10, 4.464920289836115e-06, 0.003830520129221428, 0.08778390483440988, 3.584607999290539, 46.98226573426573,
58 1.0695333890657086e-12, 2.3656367312710969e-07, 4.972309783123969e-05, 0.013643024456211269, 0.21542115380351798, 6.472046810332524, 87.1329124642076,
59 1.2233225038333268e-14, 9.27700678929842e-10, 4.468579053461083e-06, 0.0002675102834551229, 0.03813282951314997, 0.4894825318738884, 11.627213264559293, 162.79300083906728
71 const REAL8 z3 = 2.0*nu*(4.0-3.0*nu);
76 const REAL8 prstar2 =
SQ(prstar);
77 const REAL8 prstar3 = prstar2*prstar;
78 const REAL8 prstar4 = prstar2*prstar2;
80 *Heff = sqrt(A*(1.0 + pphi2*
u2) + prstar2 +
z3*A*
u2*prstar4);
81 *
H = sqrt( 1.0 + 2.0*nu*(*Heff - 1) )/nu;
83 if (dHeff_dr != NULL) *dHeff_dr = 0.5*(dA + (pphi2 +
z3*prstar4)*(dA*
u2 - 2*A*
u3))/(*Heff);
84 if (dHeff_dprstar != NULL) *dHeff_dprstar = (prstar +
z3*2.0*A*
u2*prstar3)/(*Heff);
85 if (dHeff_dpphi != NULL) *dHeff_dpphi = A*pphi*
u2/(*Heff);
110 REAL8 *dHeff_dprstar,
112 REAL8 *d2Heff_dprstar20
116 const REAL8 z3 = 2.0*nu*(4.0-3.0*nu);
118 const REAL8 prstar2 =
SQ(prstar);
119 const REAL8 UNUSED prstar3 = prstar2*prstar;
120 const REAL8 prstar4 = prstar2*prstar2;
121 const REAL8 uc = 1./rc;
122 const REAL8 uc2 = uc*uc;
123 const REAL8 uc3 = uc2*uc;
127 eob_dyn_s_GS(
r, rc, drc_dr, aK2, prstar, pphi, nu, chi1, chi2, X1, X2, c3, ggm);
128 const REAL8 GS = ggm[2];
129 const REAL8 GSs = ggm[3];
130 const REAL8 dGS_dprstar = ggm[4];
131 const REAL8 dGSs_dprstar = ggm[5];
132 const REAL8 dGS_dr = ggm[6];
133 const REAL8 dGSs_dr = ggm[7];
134 const REAL8 dGSs_dpphi = ggm[9];
135 const REAL8 d2GS_dprstar20 = ggm[12];
136 const REAL8 d2GSs_dprstar20 = ggm[13];
139 *Heff_orb = sqrt( prstar2+A*(1. + pphi2*uc2 +
z3*prstar4*uc2) );
140 *Heff = *Heff_orb + (GS*S + GSs*Sstar)*pphi;
141 *
H = sqrt( 1. + 2.*nu*(*Heff - 1.) )/nu;
142 if (dHeff_dr != NULL) *dHeff_dr = pphi*(dGS_dr*S + dGSs_dr*Sstar) + 1./(2.*(*Heff_orb))*( dA*(1. + pphi2*uc2 +
z3*prstar4*uc2) - 2.*A*uc3*drc_dr*(pphi2 +
z3*prstar4) );
143 if (dHeff_dprstar != NULL) *dHeff_dprstar = pphi*(dGS_dprstar*S + dGSs_dprstar*Sstar) + (prstar/(*Heff_orb))*(1. + 2.*A*uc2*
z3*prstar2);
144 if (d2Heff_dprstar20 != NULL) *d2Heff_dprstar20 = pphi*(d2GS_dprstar20*S + d2GSs_dprstar20*Sstar) + (1./(*Heff_orb))*(1. + 2.*A*uc2*
z3*prstar2);
145 if (dHeff_dpphi != NULL) *dHeff_dpphi = GS*S + (GSs + pphi*dGSs_dpphi)*Sstar + pphi*A*uc2/(*Heff_orb);
154 REAL8 rorb,
A,
dA,
rc,
drc_dr,
ak2,
S,
Ss,
nu,
chi1,
chi2,
X1,
X2,
c3;
180 eob_dyn_s_GS(
rorb,
rc,
drc_dr,
ak2, 0.,
x,
nu,
chi1,
chi2,
X1,
X2,
c3, ggm0);
181 REAL8 dGS_dr = ggm0[6];
182 REAL8 dGSs_dr = ggm0[7];
190 REAL8 Horbeff0 = sqrt(
A*(1. + x2*uc2));
191 REAL8 dHeff_dr =
x*(dGS_dr*
S + dGSs_dr*
Ss) + 1./(2.*Horbeff0)*(
dA*(1. + x2*uc2) - 2.*
A*uc3*
drc_dr*x2);
216 const gsl_root_fsolver_type *T;
220 double x_lo = 0.5*(double) pph, x_hi = 1.5*(
double) pph;
223 struct DHeff0_tmp_params p = {
rorb,
A,
dA,
rc,
drc_dr,
ak2,
S,
Ss,
nu,
chi1,
chi2,
X1,
X2,
c3};
227 T = gsl_root_fsolver_bisection;
228 s = gsl_root_fsolver_alloc (T);
229 gsl_root_fsolver_set (
s, &F, x_lo, x_hi);
233 status = gsl_root_fsolver_iterate (
s);
234 r = gsl_root_fsolver_root (
s);
235 x_lo = gsl_root_fsolver_x_lower (
s);
236 x_hi = gsl_root_fsolver_x_upper (
s);
240 gsl_root_fsolver_free (
s);
249 return pow(omg_orb0, -2./3.);
295 REAL8 A,
B,dA,rc,drc_dr,
G,dG_dr,uc,uc2,dAuc2_dr,j02,j0,
H,Heff,Heff_orb,dHeff_dj0,omg_orb;
296 REAL8 pl_hold,a_coeff,b_coeff,c_coeff,Delta,sol_p,sol_m;
303 dyn->
eob_dyn_s_get_rc(
r, nu, a1,
a2, aK2, C_Q1, C_Q2, C_Oct1, C_Oct2, C_Hex1, C_Hex2, usetidal, &rc, &drc_dr, &pl_hold);
304 eob_dyn_s_GS(
r, rc, drc_dr, aK2, 0.0, 0.0, nu, chi1, chi2, X1, X2, c3, ggm);
305 G = ggm[2]*S + ggm[3]*Sstar;
306 dG_dr = ggm[6]*S + ggm[7]*Sstar;
320 dAuc2_dr = uc2*(dA-2*A*uc*drc_dr);
327 a_coeff =
SQ(dAuc2_dr) - 4*A*uc2*
SQ(dG_dr);
328 b_coeff = 2*dA*dAuc2_dr - 4*A*
SQ(dG_dr);
331 Delta =
SQ(b_coeff) - 4*a_coeff*c_coeff;
333 if (S==0 && Sstar==0)
336 sol_p = (-b_coeff + sqrt(Delta))/(2*a_coeff);
337 sol_m = (-b_coeff - sqrt(Delta))/(2*a_coeff);
348 j02 = -b_coeff/a_coeff;
354 Heff_orb = sqrt(A*(1+j02*uc2));
355 Heff = Heff_orb + j0*
G;
356 H = sqrt(1+2*nu*(Heff-1))/nu;
359 dHeff_dj0 =
G + A*j0*uc2/Heff_orb;
360 omg_orb = dHeff_dj0/nu/
H;
371 const gsl_root_fsolver_type *T;
374 double x_lo = 0.5*r0_kepl, x_hi = 1.5*r0_kepl;
381 T = gsl_root_fsolver_bisection;
382 s = gsl_root_fsolver_alloc (T);
383 gsl_root_fsolver_set (
s, &F, x_lo, x_hi);
386 status = gsl_root_fsolver_iterate (
s);
387 r0 = (
REAL8) gsl_root_fsolver_root (
s);
388 x_lo = gsl_root_fsolver_x_lower (
s);
389 x_hi = gsl_root_fsolver_x_upper (
s);
393 gsl_root_fsolver_free (
s);
414 const REAL8 z3 = 2.0*nu*(4.0-3.0*nu);
417 #define TEOB_IC_N (6)
426 REAL8 H0eff, H0, psi, r_omega, v_phi, jhat,
x, dprstardt;
437 j2[
i] = r3*dA[
i]/(2.*A-
r[
i]*dA[
i]);
440 djdr[
i] = -j3/r3*( 2.0 - 3.0*A/(
r[
i]*dA[
i]) - A*d2A/(dA[
i]*dA[
i]) );
443 H0eff = sqrt(A*(1.0 + j2[
i]/
r2));
444 E0[
i] = sqrt(1.0 + 2.0*nu*(H0eff - 1.0) );
446 Omega_j[
i] = A*j[
i]/(nu*
r2*H0*H0eff);
447 psi = 2.*(1.0 + 2.0*nu*(H0eff - 1.0))/(
r2*dA[
i]);
448 r_omega =
r[
i]*cbrt(psi);
449 v_phi = Omega_j[
i]*r_omega;
451 jhat = j[
i]/(r_omega*v_phi);
456 Ctmp[
i] = sqrt(
B/A)*nu*H0*H0eff;
457 prstar[
i] = Ctmp[
i]*Fphi[
i]/djdr[
i];
460 pr[
i] = prstar[
i]*sqrt(
B/A);
468 dprstardt = dprstardr[
i] * Fphi[
i]/djdr[
i];
469 pph[
i] = j[
i]*sqrt(1. + 2.*Ctmp[
i]/dA[
i]*dprstardt -
z3*gsl_pow_int(prstar[
i],4)/j2[
i]);
504 const REAL8 S = S1 + S2;
505 const REAL8 Ss = X2*a1 + X1*
a2;
506 const REAL8 z3 = 2.0*nu*(4.0-3.0*nu);
509 #define TEOB_IC_S_N (6)
517 REAL8 pphorb, uc, uc2, psic, r_omg, v_phi, jhat,
x, Omg;
518 REAL8 UNUSED H0eff, H0, Horbeff0, Heff0, one_H0, dHeff_dprstarbyprstar, dHeff_dpph, Heff,
H, Horbeff;
519 REAL8 ggm0[14], GS_0, GSs_0, dGS_dr_0, dGSs_dr_0, dGSs_dpph_0, dGS_dprstarbyprstar_0, dGSs_dprstarbyprstar_0, GS, GSs, dGS_dr, dGSs_dr;
521 REAL8 Gtilde, dGtilde_dr, duc_dr;
531 pphorb =
r[
i]/sqrt(
r[
i]-3.);
532 dyn->
eob_dyn_s_get_rc(
r[
i], nu, a1,
a2, aK2, C_Q1, C_Q2, C_Oct1, C_Oct2, C_Hex1,
535 dA[
i], rc[
i], drc_dr[
i], aK2, S, Ss);
547 sqrtAbyB = sqrt(A[
i]/
B[
i]);
552 Horbeff0 = sqrt(A[
i]*(1. +
SQ(pph[
i])*uc2));
555 eob_dyn_s_GS(
r[
i], rc[
i], drc_dr[
i], aK2, 0, pph[
i], nu, chi1, chi2, X1, X2, c3, ggm0);
560 dGSs_dpph_0 = ggm0[9];
561 dGS_dprstarbyprstar_0 = ggm0[10];
562 dGSs_dprstarbyprstar_0 = ggm0[11];
565 Heff0 = (GS_0*S + GSs_0*Ss)*pph[
i] + Horbeff0;
568 H0 = sqrt( 1. + 2.*nu*(Heff0 - 1.));
572 dHeff_dprstarbyprstar = pph[
i]*(dGS_dprstarbyprstar_0*S + dGSs_dprstarbyprstar_0*Ss) + 1./Horbeff0;
574 C0 = sqrtAbyB*one_H0*dHeff_dprstarbyprstar;
577 dHeff_dpph = GS_0*S + (GSs_0 + pph[
i]*dGSs_dpph_0)*Ss + pph[
i]*A[
i]*uc2/Horbeff0;
578 Omg = one_H0*dHeff_dpph;
581 Gtilde = GS_0*S + GSs_0*Ss;
582 dGtilde_dr = dGS_dr_0*S + dGSs_dr_0*Ss;
583 duc_dr = -uc2*drc_dr[
i];
584 psic = (duc_dr + dGtilde_dr*rc[
i]*sqrt(A[
i]/(
SQ(pph[
i])) + A[
i]*uc2)/A[
i])/(-0.5*dA[
i]);
585 r_omg = pow((pow(rc[
i]*rc[
i]*rc[
i]*psic,-1./2)+Gtilde)*one_H0,-2./3.);
588 jhat = pph[
i]/(r_omg*v_phi);
591 prstar[
i] = Fphi[
i]/(dpph_dr[
i]*
C0);
592 pr[
i] = prstar[
i]/sqrtAbyB;
600 #if (POSTPOSTCIRCULAR)
613 sqrtAbyB = sqrt(A[
i]/
B[
i]);
617 dpi1bydj = dprstardr[
i]/djdr[
i];
618 dpi1dt = dpi1bydj*Fphi[
i];
619 prstar4 = prstar[
i]*prstar[
i]*prstar[
i]*prstar[
i];
622 Horbeff = sqrt(A[
i]*(1. +
SQ(pph[
i])*uc2));
624 eob_dyn_s_GS(
r[
i], rc[
i], drc_dr[
i], aK2, 0, pph[
i], nu, chi1, chi2, X1, X2, c3, ggm0);
631 Heff = (GS*S + GSs*Ss)*pph[
i] + Horbeff;
634 H = sqrt( 1. + 2.*nu*(Heff - 1.));
637 a = -sqrtAbyB*uc2/(2.*
H*Horbeff)*(dA[
i] - 2.*A[
i]*uc*drc_dr[
i]);
638 b = -sqrtAbyB/
H*(dGS_dr*S + dGSs_dr*Ss);
639 c = -dpi1dt - sqrtAbyB/(2.*
H*Horbeff)*(dA[
i] +
z3*prstar4*uc2*(dA[
i] - 2.*A[
i]*uc*drc_dr[
i]));
642 pph[
i] = 0.5*(-b + sqrt(b*b-4*
a*
c))/
a;
669 const REAL8 z3 = 2.0*nu*(4.0-3.0*nu);
682 REAL8 H, Heff, dHeff_dr,dHeff_dprstar;
683 eob_ham(nu,
r,pphi,prstar,A,dA, &
H,&Heff,&dHeff_dr,&dHeff_dprstar,NULL);
691 const REAL8 prstar2 = prstar*prstar;
692 const REAL8 prstar3 = prstar2*prstar;
693 const REAL8 prstar4 = prstar3*prstar;
694 const REAL8 sqrtAbyB = sqrt(A/
B);
695 const REAL8 divHE = 1./(Heff*E);
696 const REAL8 Omega = A*pphi*
u2*divHE;
705 dy[
TEOB_EVOLVE_PRSTAR] = - 0.5*sqrtAbyB*( pphi2*
u2*(dA-2.0*A*
u) + dA + 2.0*nu*(4.0-3.0*nu)*(dA*
u2 - 2.0*A*
u3)*prstar4 )*divHE;
708 const REAL8 sqrtW = sqrt(A*(1. + pphi2*
u2));
709 const REAL8 psi = 2.*(1.0 + 2.0*nu*(sqrtW - 1.0))/(
SQ(
r)*dA);
711 const REAL8 r_omega =
r*cbrt(psi);
712 const REAL8 v_phi = r_omega*Omega;
713 const REAL8 x = v_phi * v_phi;
714 const REAL8 jhat = pphi/(r_omega*v_phi);
715 const REAL8 tmpE = 1./Heff+nu/(E*E);
718 const REAL8 ddotr_dr = sqrtAbyB*( (prstar +
z3*2.*A*
u2*prstar3)*(0.5*(dA/A-dB/
B)-dHeff_dr*tmpE)+ 2.0*
z3*(dA*
u2 - 2.*A*
u3)*prstar3)*divHE;
719 const REAL8 ddotr_dprstar = sqrtAbyB*( 1.+
z3*6.*A*
u2*prstar2-(prstar +
z3*2.*A*
u2*prstar3)*dHeff_dprstar*tmpE)*divHE;
722 const REAL8 ddotr = dprstar_dt*ddotr_dprstar + dr_dt*ddotr_dr;
789 const REAL8 pphi2 = pphi*pphi;
796 REAL8 rc, drc_dr, d2rc_dr;
797 dyn->
eob_dyn_s_get_rc(
r, nu, a1,
a2, aK2, C_Q1, C_Q2, C_Oct1, C_Oct2, C_Hex1, C_Hex2, usetidal, &rc, &drc_dr, &d2rc_dr);
798 const REAL8 uc = 1./rc;
799 const REAL8 uc2 = uc*uc;
800 const REAL8 UNUSED uc3 = uc2*uc;
803 REAL8 Heff_orb, Heff,
H, dHeff_dr, dHeff_dprstar, d2Heff_dprstar20, dHeff_dpphi;
804 eob_ham_s(nu,
r, rc, drc_dr, pphi, prstar, S, Sstar, chi1, chi2, X1, X2, aK2, c3, A, dA,
805 &
H, &Heff, &Heff_orb, &dHeff_dr, &dHeff_dprstar, &dHeff_dpphi, &d2Heff_dprstar20);
810 const REAL8 ooH = 1./E;
812 const REAL8 sqrtAbyB = sqrt(A/
B);
813 const REAL8 dp_rstar_dt_0 = - sqrtAbyB*dHeff_dr*ooH;
814 const REAL8 ddotr_dp_rstar = sqrtAbyB*d2Heff_dprstar20*ooH;
815 const REAL8 Omg = dHeff_dpphi*ooH;
816 const REAL8 ddotr = dp_rstar_dt_0*ddotr_dp_rstar;
831 eob_dyn_s_GS(
r, rc, drc_dr, aK2, 0., pphi, nu, chi1, chi2, X1, X2, c3, ggm0);
833 const REAL8 GS_0 = ggm0[2];
834 const REAL8 GSs_0 = ggm0[3];
835 const REAL8 dGS_dr_0 = ggm0[6];
836 const REAL8 dGSs_dr_0 = ggm0[7];
837 const REAL8 Heff_orb_0 = sqrt(A*(1.0 + pphi2*uc2));
838 const REAL8 Heff_0 = Heff_orb_0 + (GS_0*S + GSs_0*Sstar)*pphi;
839 const REAL8 H0 = sqrt(1.0 + 2.0*nu*(Heff_0 - 1.0) );
840 const REAL8 ooH0 = 1./H0;
841 const REAL8 Gtilde = GS_0*S + GSs_0*Sstar;
842 const REAL8 dGtilde_dr = dGS_dr_0*S + dGSs_dr_0*Sstar;
843 const REAL8 duc_dr = -uc2*drc_dr;
844 const REAL8 psic = (duc_dr + dGtilde_dr*rc*sqrt(A/pphi2 + A*uc2)/A)/(-0.5*dA);
845 const REAL8 r_omg = pow( ((1./sqrt(rc*rc*rc*psic))+Gtilde)*ooH0, -2./3. );
846 const REAL8 v_phi = r_omg*Omg;
847 const REAL8 x = v_phi*v_phi;
848 const REAL8 jhat = pphi/(r_omg*v_phi);
906 REAL8 c10,c20,c30,c02,c12,c04;
907 REAL8 cs10,cs20,cs30,cs40,cs02,cs12,cs04;
917 c20 = 51./8.*nu + 41./256.*nu2;
920 c12 = 12.*nu - 49./128.*nu2;
921 c04 = -5./16.*nu + 169./256.*nu2;
923 cs10 = 3./4. + nu/2.;
924 cs20 = 27./16. + 29./4.*nu + 3./8.*nu2;
925 cs02 = 5./4. + 3./2.*nu;
926 cs12 = 4. + 11.*nu - 7./8.*nu2;
927 cs04 = 5./48. + 25./12.*nu + 3./8.*nu2;
928 cs30 = nu*cN3LO + 135./32.;
939 REAL8 prstar2 = prstar*prstar;
940 REAL8 prstar4 = prstar2*prstar2;
943 REAL8 dGS0_duc = 2.*
u2/drc_dr + 4.*
u*uc;
945 REAL8 GSs0 = 3./2.*uc3;
946 REAL8 dGSs0_duc = 9./2.*uc2;
947 REAL8 dGSs0_dprstar = 0.0;
948 REAL8 dGSs0_dpph = 0.0;
950 REAL8 hGS = 1./(1. + c10*uc + c20*uc2 + c30*uc3 + c02*prstar2 + c12*uc*prstar2 + c04*prstar4);
951 REAL8 hGSs = 1./(1. + cs10*uc + cs20*uc2 + cs30*uc3 + cs40*uc4 + cs02*prstar2 + cs12*uc*prstar2 + cs04*prstar4);
955 REAL8 GSs = GSs0*hGSs;
958 REAL8 dhGS_dprstar = -2.*prstar*hGS*hGS *( c02 + c12*uc + 2.*c04*prstar2);
959 REAL8 dhGSs_dprstar = -2.*prstar*hGSs*hGSs*(cs02 + cs12*uc + 2.*cs04*prstar2);
961 REAL8 dGS_dprstar = GS0 *dhGS_dprstar;
962 REAL8 dGSs_dprstar = GSs0*dhGSs_dprstar + dGSs0_dprstar*hGSs;
965 REAL8 dhGS_duc = -hGS*hGS*(c10 + 2.*c20*uc + 3.*c30*uc2 + c12*prstar2);
966 REAL8 dhGSs_duc = -hGSs*hGSs*(cs10 + 2.*cs20*uc + 3.*cs30*uc2 + 4.*cs40*uc3 + cs12*prstar2);
969 REAL8 dGS_duc = dGS0_duc*hGS + GS0*dhGS_duc;
970 REAL8 dGSs_duc = dGSs0_duc*hGSs + GSs0*dhGSs_duc;
973 REAL8 dGS_dr = -drc_dr*uc2*dGS_duc;
974 REAL8 dGSs_dr = -drc_dr*uc2*dGSs_duc;
978 REAL8 dGSs_dpph = dGSs0_dpph*hGSs;
981 const REAL8 dGS_dprstarbyprstar = -2.*GS0*hGS*hGS *( c02 + c12*uc + 2.*c04*prstar2);
982 const REAL8 dGSs_dprstarbyprstar = -2.*GSs0*hGSs*hGSs*(cs02 + cs12*uc + 2.*cs04*prstar2);
985 const REAL8 d2GS_dprstar20 = GS0*(-2.*hGS*hGS *( c02 + c12*uc + 2.*c04*prstar2));
986 const REAL8 d2GSs_dprstar20 = GSs0*(-2.*hGSs*hGSs*(cs02 + cs12*uc + 2.*cs04*prstar2));
998 ggm[10]=dGS_dprstarbyprstar;
999 ggm[11]=dGSs_dprstarbyprstar;
1000 ggm[12]=d2GS_dprstar20;
1001 ggm[13]=d2GSs_dprstar20;
1025 REAL8 UNUSED C_Oct1,
1026 REAL8 UNUSED C_Oct2,
1027 REAL8 UNUSED C_Hex1,
1028 REAL8 UNUSED C_Hex2,
1041 #if (RC_EXCLUDESPINSPINTIDES)
1059 REAL8 a02 = C_Q1*at1*at1 + 2.*at1*at2 + C_Q2*at2*at2;
1062 *drc_dr =
r/(*rc)*(1.-a02*
u3);
1063 *d2rc_dr2 = 1./(*rc)*(1.-(*drc_dr)*
r/(*rc)*(1.-a02*
u3)+2.*a02*
u3);
1075 REAL8 X12 = sqrt(1.-4.*nu);
1076 REAL8 c_ss_nlo = (- at2*at2*(1.25 + 1.25*X12 + 0.5*nu) - at1*at1*(1.25 - 1.25*X12 + 0.5*nu) + at1*at2*(-2.+nu));
1077 REAL8 rc2 =
r2 + aK2*(1. + 2.*
u) +
u*c_ss_nlo;
1079 REAL8 divrc = 1.0/(*rc);
1080 *drc_dr =
r*divrc*(1-(aK2 + 0.5*c_ss_nlo)*
u3);
1081 *d2rc_dr2 = divrc*(1.-(*drc_dr)*
r*divrc*(1.-(aK2+0.5*c_ss_nlo)*
u3)+ (2.*aK2 + c_ss_nlo)*
u3);
1095 REAL8 UNUSED C_Oct1,
1096 REAL8 UNUSED C_Oct2,
1097 REAL8 UNUSED C_Hex1,
1098 REAL8 UNUSED C_Hex2,
1109 REAL8 X12 = sqrt(1.-4.*nu);
1114 REAL8 a02 = C_Q1*at1*at1 + 2.*at1*at2 + C_Q2*at2*at2;
1116 REAL8 delta_a2 = X12*(at1*at1*(C_Q1+0.25) - at2*at2*(C_Q2+0.25))
1117 + at1*at1*(-17./4.+3.*C_Q1-0.5*nu)
1118 + at2*at2*(-17./4.+3.*C_Q2-0.5*nu)
1121 REAL8 rc2 =
r2 + a02*(1. + 2.*
u) + delta_a2*
u;
1123 REAL8 divrc = 1.0/(*rc);
1124 *drc_dr = divrc*(
r - (a02 + 0.5*delta_a2)*
u2);
1125 *d2rc_dr2 = divrc*(1 + (2.*a02 + delta_a2)*
u3 - (*drc_dr)*(*drc_dr));
1130 REAL8 c_ss_nlo = (- at2*at2*(1.25 + 1.25*X12 + 0.5*nu) - at1*at1*(1.25 - 1.25*X12 + 0.5*nu) + at1*at2*(-2.+nu));
1131 REAL8 rc2 =
r2 + aK2*(1. + 2.*
u) +
u*c_ss_nlo;
1133 REAL8 divrc = 1.0/(*rc);
1134 *drc_dr =
r*divrc*(1-(aK2 + 0.5*c_ss_nlo)*
u3);
1135 *d2rc_dr2 = divrc*(1.-(*drc_dr)*
r*divrc*(1.-(aK2+0.5*c_ss_nlo)*
u3)+ (2.*aK2 + c_ss_nlo)*
u3);
1149 REAL8 UNUSED C_Oct1,
1150 REAL8 UNUSED C_Oct2,
1151 REAL8 UNUSED C_Hex1,
1152 REAL8 UNUSED C_Hex2,
1165 REAL8 X12 = sqrt(1.-4.*nu);
1171 REAL8 a02 = C_Q1*at1*at1 + 2.*at1*at2 + C_Q2*at2*at2;
1173 REAL8 delta_a2 = X12*(at1*at1*(C_Q1+0.25) - at2*at2*(C_Q2+0.25))
1174 + at1*at1*(-17./4.+3.*C_Q1-0.5*nu)
1175 + at2*at2*(-17./4.+3.*C_Q2-0.5*nu)
1178 REAL8 delta_a2_nnlo =
1179 ( 387./28. - 207./28.*nu ) *a02
1180 + (-2171./212. - 269./28.*nu + 0.375*nu*nu) *(at1*at1+at2*at2)
1181 + (- 281./7 - 187./56.*nu - 0.75 *nu*nu) *at1*at2
1182 + 163./28. *X12*(C_Q1*at1*at1-C_Q2*at2*at2)
1183 + ( -29./112. - 2.625 *nu ) *X12*(at1*at1-at2*at2);
1185 REAL8 rc2 =
r2 + a02*(1. + 2.*
u) + delta_a2*
u + delta_a2_nnlo*
u2;
1187 REAL8 divrc = 1.0/(*rc);
1188 *drc_dr = divrc*(
r - (a02 + 0.5*delta_a2)*
u2 - delta_a2_nnlo*
u3);
1189 *d2rc_dr2 = divrc*(1 + (2.*a02 + delta_a2)*
u3
1190 + 3*delta_a2_nnlo*
u4 - (*drc_dr)*(*drc_dr));
1196 REAL8 a0 = at1 + at2;
1197 REAL8 a12 = at1 - at2;
1199 REAL8 c_ss_nlo = -1.125*a0*a0 -(0.125+0.5+nu)*a12*a12 + 1.25*X12*a0*a12;
1201 REAL8 c_ss_nnlo = - (189./32. + 417.32*nu) *a0 *a0
1202 + ( 11./32. - 127.32*nu + 0.375*nu*nu) *a12*a12
1203 + ( 87.16 - 2.625*nu )*X12 *a0 *a12;
1206 REAL8 rc2 =
r2 + aK2*(1. + 2.*
u) +
u*c_ss_nlo +
u2*c_ss_nnlo;
1208 REAL8 divrc = 1.0/(*rc);
1209 *drc_dr =
r*divrc*(1-(aK2 + 0.5*c_ss_nlo)*
u3 - 0.5*
u4*c_ss_nnlo);
1210 *d2rc_dr2 = 1./
r*(*drc_dr) +
r*divrc*((3.*aK2+c_ss_nlo)*
u4 + 2.*c_ss_nnlo*
u5);
1240 REAL8 X12 = sqrt(1.-4.*nu);
1246 REAL8 a02 = C_Q1*at1*at1 + 2.*at1*at2 + C_Q2*at2*at2;
1248 REAL8 delta_a2 = X12*(at1*at1*(C_Q1+0.25) - at2*at2*(C_Q2+0.25))
1249 + at1*at1*(-17./4.+3.*C_Q1-0.5*nu)
1250 + at2*at2*(-17./4.+3.*C_Q2-0.5*nu)
1253 REAL8 delta_a2_nnlo =
1254 ( 387./28. - 207./28.*nu ) *a02
1255 + (-2171./212. - 269./28.*nu + 0.375*nu*nu) *(at1*at1+at2*at2)
1256 + (- 281./7 - 187./56.*nu - 0.75 *nu*nu) *at1*at2
1257 + 163./28. *X12*(C_Q1*at1*at1-C_Q2*at2*at2)
1258 + ( -29./112. - 2.625 *nu ) *X12*(at1*at1-at2*at2);
1260 REAL8 delta_a4_lo = 0.75*(C_Hex1 - C_Q1*C_Q1)*at1*at1*at1*at1
1261 + 3.*(C_Oct1 - C_Q1) *at1*at1*at1*at2
1262 + 3.*(C_Q1*C_Q2 - 1) *at1*at1*at2*at2
1263 + 3.*(C_Oct2 - C_Q2) *at1*at2*at2*at2
1264 + 0.75*(C_Hex2 - C_Q2*C_Q2)*at2*at2*at2*at2;
1266 REAL8 rc2 =
r2 + a02*(1. + 2.*
u) + delta_a2*
u + (delta_a2_nnlo+delta_a4_lo)*
u2;
1268 REAL8 divrc = 1.0/(*rc);
1269 *drc_dr = divrc*(
r - (a02 + 0.5*delta_a2)*
u2 - (delta_a2_nnlo+delta_a4_lo)*
u3);
1270 *d2rc_dr2 = divrc*(1 + (2.*a02 + delta_a2)*
u3
1271 + 3*(delta_a2_nnlo+delta_a4_lo)*
u4 - (*drc_dr)*(*drc_dr));
1277 REAL8 a0 = at1 + at2;
1278 REAL8 a12 = at1 - at2;
1280 REAL8 c_ss_nlo = -1.125*a0*a0 -(0.125+0.5+nu)*a12*a12 + 1.25*X12*a0*a12;
1282 REAL8 c_ss_nnlo = - (189./32. + 417.32*nu ) *a0 *a0
1283 + ( 11./32. - 127.32*nu + 0.375*nu*nu) *a12*a12
1284 + ( 87.16 - 2.625*nu )*X12*a0 *a12;
1287 REAL8 rc2 =
r2 + aK2*(1. + 2.*
u) +
u*c_ss_nlo +
u2*c_ss_nnlo;
1289 REAL8 divrc = 1.0/(*rc);
1290 *drc_dr =
r*divrc*(1-(aK2 + 0.5*c_ss_nlo)*
u3 - 0.5*
u4*c_ss_nnlo);
1291 *d2rc_dr2 = 1./
r*(*drc_dr) +
r*divrc*((3.*aK2+c_ss_nlo)*
u4 + 2.*c_ss_nnlo*
u5);
1305 REAL8 UNUSED C_Oct1,
1306 REAL8 UNUSED C_Oct2,
1307 REAL8 UNUSED C_Hex1,
1308 REAL8 UNUSED C_Hex2,
1309 int UNUSED usetidal,
1311 REAL8 UNUSED *drc_dr,
1312 REAL8 UNUSED *d2rc_dr2)
1327 REAL8 UNUSED C_Oct1,
1328 REAL8 UNUSED C_Oct2,
1329 REAL8 UNUSED C_Hex1,
1330 REAL8 UNUSED C_Hex2,
1345 REAL8 a02 = 2.*at1*at2;
1348 *drc_dr =
r/(*rc)*(1.-a02*
u3);
1349 *d2rc_dr2 = 1./(*rc)*(1.-(*drc_dr)*
r/(*rc)*(1.-a02*
u3)+2.*a02*
u3);
1353 REAL8 X12 = sqrt(1.-4.*nu);
1354 REAL8 c_ss_nlo = (- at2*at2*(1.25 + 1.25*X12 + 0.5*nu) - at1*at1*(1.25 - 1.25*X12 + 0.5*nu) + at1*at2*(-2.+nu));
1355 REAL8 rc2 =
r2 + aK2*(1. + 2.*
u) +
u*c_ss_nlo;
1357 REAL8 divrc = 1.0/(*rc);
1358 *drc_dr =
r*divrc*(1-(aK2 + 0.5*c_ss_nlo)*
u3);
1359 *d2rc_dr2 = divrc*(1.-(*drc_dr)*
r*divrc*(1.-(aK2+0.5*c_ss_nlo)*
u3)+ (2.*aK2 + c_ss_nlo)*
u3);
1386 REAL8 pi4 = pi2*pi2;
1400 REAL8 a5c1 = -221./6. + 41./32.*pi2;
1401 REAL8 a5 = a5c0 + nu*a5c1;
1402 REAL8 a6 = 3097.3*nu2 - 1330.6*nu + 81.38;
1406 REAL8 a6tot = a6 + (-7004./105. - 144./5.*nu)*
logu;
1407 REAL8 a5tot2 = a5tot*a5tot;
1410 REAL8 N1 = (-3*(-512 - 32*nu2 + nu*(3520 + 32*a5tot + 8*a6tot - 123*pi2)))/(-768 + nu*(3584 + 24*a5tot - 123*pi2));
1411 REAL8 D1 = (nu*(-3392 - 48*a5tot - 24*a6tot + 96*nu + 123*pi2))/(-768 + nu*(3584 + 24*a5tot - 123*pi2));
1412 REAL8 D2 = (2*nu*(-3392 - 48*a5tot - 24*a6tot + 96*nu + 123*pi2))/(-768 + nu*(3584 + 24*a5tot - 123*pi2));
1413 REAL8 D3 = (-2*nu*(6016 + 48*a6tot + 3392*nu + 24*a5tot*(4 + nu) - 246*pi2 - 123*nu*pi2))/(-768 + nu*(3584 + 24*a5tot - 123*pi2));
1414 REAL8 D4 = -(nu*(-4608*a6tot*(-4 + nu) + a5tot*(36864 + nu*(72192 - 2952*pi2)) + nu*(2048*(5582 + 9*nu) - 834432*pi2 + 15129*pi4)))/(96.*(-768 + nu*(3584 + 24*a5tot - 123*pi2)));
1415 REAL8 D5 = (nu*(-24*a6tot*(1536 + nu*(-3776 + 123*pi2)) + nu*(-2304*a5tot2 + 96*a5tot*(-3392 + 123*pi2) - (-3776 + 123*pi2)*(-3008 - 96*nu + 123*pi2))))/(96.*(-768 + nu*(3584 + 24*a5tot - 123*pi2)));
1418 REAL8 dN1 = (160*nu*(-828672 - 32256*nu2 + 756*nu*(-768 + nu*(3584 + 24*a5 - 123*pi2)) + nu*(5006848 + 42024*a5 + 8064*a6 - 174045*pi2)))/(7.*gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),2)*
u);
1419 REAL8 dD1 = (160*nu*(-828672 - 32256*nu2 + 756*nu*(-768 + nu*(3584 + 24*a5 - 123*pi2)) + nu*(5006848 + 42024*a5 + 8064*a6 - 174045*pi2)))/(7.*gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),2)*
u);
1420 REAL8 dD2 = (320*nu*(-828672 - 32256*nu2 + 756*nu*(-768 + nu*(3584 + 24*a5 - 123*pi2)) + nu*(5006848 + 42024*a5 + 8064*a6 - 174045*pi2)))/(7.*gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),2)*
u);
1421 REAL8 dD3 = (640*nu*(-828672 - 32256*nu2 + 756*nu*(-768 + nu*(3584 + 24*a5 - 123*pi2)) + nu*(5006848 + 42024*a5 + 8064*a6 - 174045*pi2)))/(7.*gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),2)*
u);
1422 REAL8 dD4 = (-320*(-4 + nu)*nu*(-828672 - 32256*nu2 + 756*nu*(-768 + nu*(3584 + 24*a5 - 123*pi2)) + nu*(5006848 + 42024*a5 + 8064*a6 - 174045*pi2)))/(7.*gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),2)*
u);
1423 REAL8 dD5 = (nu*(-8400*nu*(-24*(a6 - (4*
logu*(1751 + 756*nu))/105.)*(1536 + nu*(-3776 + 123*pi2)) + nu*(-2304*gsl_pow_int(a5 + (64*
logu)/5.,2) + 96*(a5 + (64*
logu)/5.)*(-3392 + 123*pi2) - (-3776 + 123*pi2)*(-32*(94 + 3*nu) + 123*pi2))) - (1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)))*(4128768*
logu*nu + 5*(-2689536 + nu*(11170624 + 64512*a5 - 380685*pi2) - 756*nu*(1536 + nu*(-3776 + 123*pi2))))))/(2625.*gsl_pow_int(-768 + nu*(3584 + 24*(a5 + (64*
logu)/5.) - 123*pi2),2)*
u);
1432 REAL8 dDen = D1 +
u*(dD1 + 2*
D2) +
u2*(dD2 + 3*D3) +
u3*(dD3 + 4*D4) +
u4*(dD4 + 5*D5) + dD5*
u5;
1435 REAL8 prefactor = (*A)/(Num*Den);
1436 REAL8 dA_u = prefactor*(dNum*Den - dDen*Num);
1446 REAL8 d2N1 = (160*nu*(-3840 + 1536*
logu*nu + nu*(20992 + 120*a5 - 615*pi2))*(828672 + nu*(-42024*a5 - 8064*a6 + 3584*(-1397 + 9*nu) + 174045*pi2) + 756*nu*(768 + nu*(-3584 - 24*a5 + 123*pi2))))/(7.*gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),3)*
u2);
1447 REAL8 d2D1 = (160*nu*(-3840 + 1536*
logu*nu + nu*(20992 + 120*a5 - 615*pi2))*(828672 + nu*(-42024*a5 - 8064*a6 + 3584*(-1397 + 9*nu) + 174045*pi2) + 756*nu*(768 + nu*(-3584 - 24*a5 + 123*pi2))))/(7.*gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),3)*
u2);
1448 REAL8 d2D2 = (320*nu*(-3840 + 1536*
logu*nu + nu*(20992 + 120*a5 - 615*pi2))*(828672 + nu*(-42024*a5 - 8064*a6 + 3584*(-1397 + 9*nu) + 174045*pi2) + 756*nu*(768 + nu*(-3584 - 24*a5 + 123*pi2))))/(7.*gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),3)*
u2);
1449 REAL8 d2D3 = (640*nu*(-3840 + 1536*
logu*nu + nu*(20992 + 120*a5 - 615*pi2))*(828672 + nu*(-42024*a5 - 8064*a6 + 3584*(-1397 + 9*nu) + 174045*pi2) + 756*nu*(768 + nu*(-3584 - 24*a5 + 123*pi2))))/(7.*gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),3)*
u2);
1450 REAL8 d2D4 = (320*(-4 + nu)*nu*(-828672 + 756*nu*(-768 + nu*(3584 + 24*a5 - 123*pi2)) + nu*(5006848 + 42024*a5 + 8064*a6 - 32256*nu - 174045*pi2))*(-3840 + 1536*
logu*nu + nu*(20992 + 120*a5 - 615*pi2)))/(7.*gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),3)*
u2);
1451 REAL8 d2D5 = (nu*(gsl_pow_int(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)),2)*(4128768*
logu*nu - 7680*(1751 + 756*nu) + nu*(64*(808193 + 5040*a5 + 223020*nu) - 615*(3095 + 756*nu)*pi2)) + 3072*nu*(1536*
logu*nu + 5*(-768 + nu*(3584 + 24*a5 - 123*pi2)))*(4128768*
logu*nu - 7680*(1751 + 756*nu) + 5*nu*(64*(174541 + 1008*a5 + 44604*nu) - 123*(3095 + 756*nu)*pi2)) + 25804800*nu2*(-24*(a6 - (4*
logu*(1751 + 756*nu))/105.)*(1536 + nu*(-3776 + 123*pi2)) + nu*(-2304*gsl_pow_int(a5 + (64*
logu)/5.,2) + 96*(a5 + (64*
logu)/5.)*(-3392 + 123*pi2) - (-3776 + 123*pi2)*(-32*(94 + 3*nu) + 123*pi2))) + 42000*nu*(-768 + nu*(3584 + 24*(a5 + (64*
logu)/5.) - 123*pi2))*(-24*(a6 - (4*
logu*(1751 + 756*nu))/105.)*(1536 + nu*(-3776 + 123*pi2)) + nu*(-2304*gsl_pow_int(a5 + (64*
logu)/5.,2) + 96*(a5 + (64*
logu)/5.)*(-3392 + 123*pi2) - (-3776 + 123*pi2)*(-32*(94 + 3*nu) + 123*pi2)))))/(13125.*gsl_pow_int(-768 + nu*(3584 + 24*(a5 + (64*
logu)/5.) - 123*pi2),3)*
u2);
1454 REAL8 d2Num = 2.*dN1 + d2N1*
u;
1455 REAL8 d2Den = 2.*(
D2 + dD1) +
u*(6.*D3 + 4.*dD2 + d2D1) +
u2*(12.*D4 + 6.*dD3 + d2D2) +
u3*(20.*D5 + 8.*dD4 + d2D3) +
u4*(10.*dD5 + d2D4) +
u5*d2D5;
1458 REAL8 d2A_u = prefactor*(2.*dDen*dDen*(*A) - 2.*dNum*dDen + Den*d2Num - d2Den*Num);
1478 REAL8 A, dA_u, d2A_u, UNUSED dA, UNUSED d2A;
1479 A = dA_u = d2A_u = 0.;
1481 const REAL8 elsix = 1.833333333333333333333;
1482 const REAL8 eightthird = 2.6666666666666666667;
1510 REAL8 UNUSED nu2 = nu*nu;
1512 REAL8 UNUSED pi4 = pi2*pi2;
1524 REAL8 oom3u = 1./(1.-rLR*
u);
1529 A = -(kapT4*u10) - kapT2*u6*(1. + bar_alph2_1*
u + bar_alph2_2*
u2) - kapT3*u8*(1. + bar_alph3_1*
u + bar_alph3_2*
u2);
1530 dA_u = -10.*kapT4*u9 - kapT2*u6*(bar_alph2_1 + 2.*bar_alph2_2*
u) - kapT3*u8*(bar_alph3_1 + 2.*bar_alph3_2*
u)
1531 - 6.*kapT2*
u5*(1. + bar_alph2_1*
u + bar_alph2_2*
u2) - 8.*kapT3*u7*(1. + bar_alph3_1*
u + bar_alph3_2*
u2);
1534 d2A_u = -90.*kapT4*u8
1535 - kapT2*(2*bar_alph2_2*u6 + 12.*
u5*(bar_alph2_1 + 2*bar_alph2_2*
u)
1536 + 30.*
u4*(1 + bar_alph2_1*
u + bar_alph2_2*
u2))
1537 - kapT3*(2.*bar_alph3_2*u8 + 16.*u7*(bar_alph3_1 + 2*bar_alph3_2*
u) + 56.*u6*(1 + bar_alph3_1*
u + bar_alph3_2*
u2));
1544 const REAL8 n1 = 0.8400636422;
1545 const REAL8 d2 = 17.7324036;
1551 REAL8 f23 = (1. + n1*
u)*Den;
1552 REAL8 df23 = (n1 - 2*d2*
u - n1*d2*
u2)*(Den*Den);
1553 REAL8 A1SF = Acub*f23;
1554 REAL8 dA1SF = dAcub*f23 + Acub*df23;
1556 REAL8 dA2SF = 674./28.*
u;
1559 REAL8 f1 = A1SF *pow(oom3u,7./2.);
1560 REAL8 f2 = A2SF *pow(oom3u,
p);
1562 REAL8 df0 = 3*
u*(2.-rLR*
u)*(oom3u*oom3u);
1563 REAL8 df1 = 0.5*(7*rLR*A1SF + 2*(1.-rLR*
u)*dA1SF)*pow(oom3u,9./2.);
1564 REAL8 df2 = (rLR*
p*A2SF + (1.-rLR*
u)*dA2SF)*pow(oom3u,
p+1);
1567 REAL8 AT2 = - kapA2*u6*( f0 + XA*f1 + XA*XA*f2 ) - kapB2*u6*( f0 + XB*f1 + XB*XB*f2 );
1568 REAL8 AT3 = - kapT3*u8*(1. + bar_alph3_1*
u + bar_alph3_2*
u2);
1569 REAL8 AT4 = - kapT4*u10;
1571 REAL8 dAT2 = - kapA2*6.*
u5*( f0 + XA*f1 + XA*XA*f2 ) - kapB2*6.*
u5*( f0 + XB*f1 + XB*XB*f2 ) - kapA2*u6*( df0 + XA*df1 + XA*XA*df2 ) - kapB2*u6*( df0 + XB*df1 + XB*XB*df2 );
1572 REAL8 dAT3 = - kapT3*(8.*u7 + 9*bar_alph3_1*u8 + 10*bar_alph3_2*u9);
1573 REAL8 dAT4 = - kapT4*10.*u9;
1575 A = AT2 + AT3 + AT4;
1576 dA_u = dAT2 + dAT3 + dAT4;
1579 REAL8 d2f23 = 2*d2*(-1 + 3*d2*
u2 + n1*
u*(-3+d2*
u2))*(Den*Den*Den);
1580 REAL8 d2A1SF = d2Acub*f23 + 2*dAcub*df23 + Acub*d2f23;
1581 REAL8 d2A2SF = 674./28.;
1582 REAL8 d2f0 = 6*(oom3u*oom3u*oom3u);
1583 REAL8 d2f1 = 0.25*(63*(rLR*rLR)*A1SF + 4*(-1+rLR*
u)*(-7*rLR*dA1SF + (-1+rLR*
u)*d2A1SF))*pow(oom3u,11./2.);
1584 REAL8 d2f2 = ( rLR*
p*(1+
p)*rLR*A2SF +(-1+rLR*
u)*( -2.*
p*rLR*dA2SF +(-1.+rLR*
u)*d2A2SF ) )*pow(oom3u,
p+2);
1586 REAL8 d2AT2 = - kapA2*30*
u4*( f0 + XA*f1 + XA*XA*f2 ) - kapB2*30*
u4*( f0 + XB*f1 + XB*XB*f2 ) - 2*kapA2*6*
u5*( df0 + XA*df1 + XA*XA*df2 ) - 2*kapB2*6*
u5*( df0 + XB*df1 + XB*XB*df2 ) - kapA2*u6*( d2f0 + XA*d2f1 + XA*XA*d2f2 ) - kapB2*u6*( d2f0 + XB*d2f1 + XB*XB*d2f2 );
1587 REAL8 d2AT3 = - kapT3*(56*u6 + 72*bar_alph3_1*u7 + 90*bar_alph3_2*u8);
1588 REAL8 d2AT4 = - kapT4*90*u8;
1590 d2A_u = d2AT2 + d2AT3 + d2AT4;
1597 const REAL8 n1 = 0.8400636422;
1598 const REAL8 d2 = 17.7324036;
1604 REAL8 f23 = (1. + n1*
u)*Den;
1605 REAL8 df23 = (n1 - 2*d2*
u - n1*d2*
u2)*(Den*Den);
1606 REAL8 A1SF = Acub*f23;
1607 REAL8 dA1SF = dAcub*f23 + Acub*df23;
1609 REAL8 dA2SF = 674./28.*
u;
1612 REAL8 f1 = A1SF *pow(oom3u,7./2.);
1613 REAL8 f2 = A2SF *pow(oom3u,
p);
1615 REAL8 df0 = 3*
u*(2.-rLR*
u)*(oom3u*oom3u);
1616 REAL8 df1 = 0.5*(7*rLR*A1SF + 2*(1.-rLR*
u)*dA1SF)*pow(oom3u,9./2.);
1617 REAL8 df2 = (rLR*
p*A2SF + (1.-rLR*
u)*dA2SF)*pow(oom3u,
p+1);
1620 REAL8 AT2 = - kapA2*u6*( f0 + XA*f1 + XA*XA*f2 ) - kapB2*u6*( f0 + XB*f1 + XB*XB*f2 );
1621 REAL8 AT4 = - kapT4*u10;
1623 REAL8 dAT2 = - kapA2*6.*
u5*( f0 + XA*f1 + XA*XA*f2 ) - kapB2*6.*
u5*( f0 + XB*f1 + XB*XB*f2 ) - kapA2*u6*( df0 + XA*df1 + XA*XA*df2 ) - kapB2*u6*( df0 + XB*df1 + XB*XB*df2 );
1624 REAL8 dAT4 = - kapT4*10.*u9;
1629 const REAL8 C1 = -3.6820949997216643;
1630 const REAL8 C2 = 5.171003322924513;
1631 const REAL8 C3 = -7.639164165720986;
1632 const REAL8 C4 = -8.63278143009751;
1633 const REAL8 C5 = 12.319646912775516;
1634 const REAL8 C6 = 16.36009385150114;
1637 REAL8 A3hat_Sch = (1.0 - 2.0*
u)*( 1.0 + eightthird*
u2*oom3u );
1638 REAL8 dA3hat_Sch = (1.0 - 2.0*
u)*( eightthird*rLR*
u2*oom3u*oom3u + 2.0*eightthird*
u*oom3u ) - 2.0*( 1.0 + eightthird*
u2*oom3u );
1639 REAL8 d2A3hat_Sch = (1.0 - 2.0*
u)*( 2.0*eightthird*rLR*rLR*
u2*oom3u*oom3u*oom3u + 4.0*eightthird*rLR*
u*oom3u*oom3u + 2.0*eightthird*oom3u ) - 4.0*( eightthird*rLR*
u2*oom3u*oom3u + 2.0*eightthird*
u*oom3u );
1645 REAL8 d2A3tilde = 15.*(
C1*(1 + 3*
C4*
u - 3*
C5*pow(
u,2) + 6*
C6*pow(
u,2) -
C4*
C5*pow(
u,3) + 3*
C5*
C6*pow(
u,4) + pow(
C5,2)*
C6*pow(
u,6)) +
C4*(1 - 3*
C5*pow(
u,2) + 10*
C3*pow(
u,3) + 9*
C3*
C5*pow(
u,5) + 3*
C3*pow(
C5,2)*pow(
u,7) +
C2*pow(
u,2)*(6 + 3*
C5*pow(
u,2) + pow(
C5,2)*pow(
u,4))) +
u*(3*(
C6 + 2*
C3*
u + 5*
C3*
C6*pow(
u,3)) +
C5*(-3 -
C6*pow(
u,2) + 3*
C3*pow(
u,3) +17*
C3*
C6*pow(
u,5)) + pow(
C5,2)*(pow(
u,2) +
C3*pow(
u,5) + 6*
C3*
C6*pow(
u,7)) +
C2*(3 + 10*
C6*pow(
u,2) + 3*pow(
C5,2)*
C6*pow(
u,6) +
C5*pow(
u,2)*(-1 + 9*
C6*pow(
u,2)))) )*Denom3*Denom3*Denom3;
1646 REAL8 A3hat1GSFfit = A3tilde*pow(oom3u, 3.5);
1647 REAL8 dA3hat1GSFfit = 3.5*rLR*A3tilde*pow(oom3u, 4.5) + dA3tilde*pow(oom3u, 3.5);
1648 REAL8 d2A3hat1GSFfit = 15.75*rLR*rLR*A3tilde*pow(oom3u, 5.5) + 7.0*rLR*dA3tilde*pow(oom3u, 4.5) + d2A3tilde*pow(oom3u, 3.5);
1651 REAL8 A3hat2GSF = 36.666666666666666667*
u2*pow(oom3u,
p);
1652 REAL8 dA3hat2GSF = 36.666666666666666667*
u*( 2. + (
p - 2.)*rLR*
u ) * pow(oom3u,
p+1);
1653 REAL8 d2A3hat2GSF = 36.666666666666666667*( 2. + 4.*(
p - 1.)*rLR*
u + (2. - 3.*
p + 1.*
p*
p)*rLR*rLR*
u2 ) * pow(oom3u,
p+2);
1656 REAL8 A3hatA = A3hat_Sch + XA*A3hat1GSFfit + XA*XA*A3hat2GSF;
1657 REAL8 dA3hatA = dA3hat_Sch + XA*dA3hat1GSFfit + XA*XA*dA3hat2GSF;
1658 REAL8 A3hatB = A3hat_Sch + XB*A3hat1GSFfit + XB*XB*A3hat2GSF;
1659 REAL8 dA3hatB = dA3hat_Sch + XB*dA3hat1GSFfit + XB*XB*dA3hat2GSF;
1662 REAL8 AT3 = -1.*kapA3*u8*( A3hatA ) - 1.*kapB3*u8*( A3hatB );
1663 REAL8 dAT3 = -1.*kapA3*u7*( 8.*A3hatA + 1.*
u*dA3hatA ) - 1.*kapB3*u7*( 8.*A3hatB + 1.*
u*dA3hatB );
1665 A = AT2 + AT3 + AT4;
1666 dA_u = dAT2 + dAT3 + dAT4;
1669 REAL8 d2f23 = 2*d2*(-1 + 3*d2*
u2 + n1*
u*(-3+d2*
u2))*(Den*Den*Den);
1670 REAL8 d2A1SF = d2Acub*f23 + 2*dAcub*df23 + Acub*d2f23;
1671 REAL8 d2A2SF = 674./28.;
1672 REAL8 d2f0 = 6*(oom3u*oom3u*oom3u);
1673 REAL8 d2f1 = 0.25*(63*(rLR*rLR)*A1SF + 4*(-1+rLR*
u)*(-7*rLR*dA1SF + (-1+rLR*
u)*d2A1SF))*pow(oom3u,11./2.);
1674 REAL8 d2f2 = ( rLR*
p*(1+
p)*rLR*A2SF +(-1+rLR*
u)*( -2.*
p*rLR*dA2SF +(-1.+rLR*
u)*d2A2SF ) )*pow(oom3u,
p+2);
1676 REAL8 d2AT2 = - kapA2*30*
u4*( f0 + XA*f1 + XA*XA*f2 ) - kapB2*30*
u4*( f0 + XB*f1 + XB*XB*f2 ) - 2*kapA2*6*
u5*( df0 + XA*df1 + XA*XA*df2 ) - 2*kapB2*6*
u5*( df0 + XB*df1 + XB*XB*df2 ) - kapA2*u6*( d2f0 + XA*d2f1 + XA*XA*d2f2 ) - kapB2*u6*( d2f0 + XB*d2f1 + XB*XB*d2f2 );
1677 REAL8 d2AT4 = - kapT4*90*u8;
1679 REAL8 d2A3hatA = d2A3hat_Sch + XA*d2A3hat1GSFfit + XA*XA*d2A3hat2GSF;
1680 REAL8 d2A3hatB = d2A3hat_Sch + XB*d2A3hat1GSFfit + XB*XB*d2A3hat2GSF;
1681 REAL8 d2AT3 = -1.*kapA3 * ( 56.*u6*A3hatA + 16.*u7*dA3hatA + 1.*u8*d2A3hatA ) - 1.*kapB3 * ( 56.*u6*A3hatB + 16.*u7*dA3hatB + 1.*u8*d2A3hatB );
1683 d2A_u += d2AT2 + d2AT3 + d2AT4;
1691 A +=-kapT2j*u7*(1. + bar_alph2j_1*
u);
1692 dA_u += -kapT2j*u7*bar_alph2j_1 - 7.*kapT2j*u6*(1. + bar_alph2j_1*
u);
1695 d2A_u += - 14.*kapT2j*
u5*(3. + 4.*bar_alph2j_1*
u);
1701 const REAL8 a1j = 0.728591192;
1702 const REAL8 a2j = 3.100367557;
1703 const REAL8 n1j = -15.04421708;
1704 const REAL8 d2j = 12.55229698;
1706 REAL8 Ahat_Schj = (1.-2.*
u)*oom3u;
1707 REAL8 dAhat_Schj = (rLR-2.)*oom3u*oom3u;
1708 REAL8 d2Ahat_Schj = 2.*rLR*(rLR-2.)*pow(oom3u, 3.);
1711 REAL8 Denomj = 1./(1. + d2j*
u2);
1712 REAL8 Ahat1GSFfitj = elsix*
u*(1. - a1j*
u)*(1. - a2j*
u)*(1. + n1j*
u)*Denomj*pow(oom3u, 3.5);
1713 REAL8 dAhat1GSFfitj = 0.5*elsix * Denomj * Denomj * (2 + 4*n1j*
u + 5*rLR*
u - 2*d2j*pow(
u,2) + 3*n1j*rLR*pow(
u,2) + 9*d2j*rLR*pow(
u,3) + 7*d2j*n1j*rLR*pow(
u,4) - a2j*
u*(4 + rLR*
u*(3 + 7*d2j*pow(
u,2)) + n1j*
u*(6 + 2*d2j*pow(
u,2) + rLR*(
u + 5*d2j*pow(
u,3)))) + a1j*
u*(-4 - 3*rLR*
u - 7*d2j*rLR*pow(
u,3) - n1j*
u*(6 + rLR*
u + 2*d2j*pow(
u,2) + 5*d2j*rLR*pow(
u,3)) + a2j*
u*(6 + rLR*
u + 2*d2j*pow(
u,2) + 5*d2j*rLR*pow(
u,3) + n1j*
u*(8 - rLR*
u + 4*d2j*pow(
u,2) + 3*d2j*rLR*pow(
u,3)))) ) * pow(oom3u, 4.5);
1714 REAL8 d2Ahat1GSFfitj = 0.25*elsix * Denomj * Denomj * Denomj * ( 8*(1 + n1j*
u)*pow(-1 + rLR*
u,2)*pow(1 + d2j*pow(
u,2),2)*(-a2j + a1j*(-1 + 3*a2j*
u)) + 4*(1 - rLR*
u)*(1 + d2j*pow(
u,2))*(1 - 2*a2j*
u + a1j*
u*(-2 + 3*a2j*
u))*(-4*d2j*
u + rLR*(7 + 11*d2j*pow(
u,2)) + n1j*(2 - 2*d2j*pow(
u,2) + rLR*
u*(5 + 9*d2j*pow(
u,2)))) +
u*(-1 + a1j*
u)*(-1 + a2j*
u)*(7*rLR*(9*rLR + n1j*(4 + 5*rLR*
u)) + 2*d2j*(-4 - 20*rLR*
u + 87*pow(rLR,2)*pow(
u,2) +3*n1j*
u*(-4 + 8*rLR*
u + 17*pow(rLR,2)*pow(
u,2))) + pow(d2j,2)*pow(
u,2)*(24 - 104*rLR*
u + 143*pow(rLR,2)*pow(
u,2) + n1j*
u*(8 - 44*rLR*
u + 99*pow(rLR,2)*pow(
u,2)))) ) * pow(oom3u, 5.5);
1717 REAL8 Ahat2GSFj =
u*pow(oom3u,
p);
1718 REAL8 dAhat2GSFj = ( 1.+ (
p-1.)*rLR*
u ) * pow(oom3u,
p+1);
1719 REAL8 d2Ahat2GSFj =
p*rLR * ( 2.+ (
p-1.)*rLR*
u ) * pow(oom3u,
p+2);
1722 REAL8 AhatjA = Ahat_Schj + XA*Ahat1GSFfitj + XA*XA*Ahat2GSFj;
1723 REAL8 dAhatjA = dAhat_Schj + XA*dAhat1GSFfitj + XA*XA*dAhat2GSFj;
1724 REAL8 AhatjB = Ahat_Schj + XB*Ahat1GSFfitj + XB*XB*Ahat2GSFj;
1725 REAL8 dAhatjB = dAhat_Schj + XB*dAhat1GSFfitj + XB*XB*dAhat2GSFj;
1728 REAL8 ATj_2 = -1.*kapA2j*u7*( AhatjA ) - 1.*kapB2j*u7*( AhatjB );
1729 REAL8 dATj_2 = -1.*kapA2j * ( 7.*u6*AhatjA + u7*dAhatjA ) - 1.*kapB2j * ( 7.*u6*AhatjB + u7*dAhatjB );
1735 REAL8 d2AhatjA = d2Ahat_Schj + XA*d2Ahat1GSFfitj + XA*XA*d2Ahat2GSFj;
1736 REAL8 d2AhatjB = d2Ahat_Schj + XB*d2Ahat1GSFfitj + XB*XB*d2Ahat2GSFj;
1737 REAL8 d2ATj_2 = -1.*kapA2j * ( 42.*
u5*AhatjA + 14.*u6*dAhatjA + u7*d2AhatjA ) - 1.*kapB2j * ( 42.*
u5*AhatjB + 14.*u6*dAhatjB + u7*d2AhatjB );
1746 if (d2AT != NULL) *d2AT = d2A_u;
1767 REAL8 Atmp=0., dAtmp_u=0., d2Atmp_u=0.;
1768 REAL8 Btmp=0., dBtmp_r=0.;
1775 REAL8 AT, dAT_u, d2AT_u;
1776 REAL8 UNUSED BT, UNUSED dBT;
1781 #if (USEBTIDALPOTENTIAL)
1784 BT = kT2*(8. - 15.*nu)*u6;
1785 dBT = -kT2*6.*(8. - 15.*nu)*
u4*
u3;
1794 *d2A = 2.*dAtmp_u*
u3 + d2Atmp_u*
u4;
1797 const REAL8 Dp = 1.0 + 6.*nu*
u2 - 2.*(3.0*nu-26.0)*nu*
u3;
1798 const REAL8 D = 1./Dp;
1799 const REAL8 dD = 6.*
u2*(2.*nu*
u-(3.*nu-26.)*nu*
u2)*D*D;
1803 dBtmp_r += (dD*(Atmp) - D*(*dA))/((Atmp)*(Atmp));
1831 REAL8 rc, drc, d2rc;
1832 dyn->
eob_dyn_s_get_rc(
r, nu, a1,
a2, aK2, C_Q1, C_Q2, C_Oct1, C_Oct2, C_Hex1, C_Hex2, usetidal, &rc, &drc, &d2rc);
1835 REAL8 Aorb, dAorb_u, d2Aorb_u;
1840 REAL8 AT, dAT_u, d2AT_u;
1851 REAL8 uc4 = uc2*uc2;
1853 REAL8 dAorb = -dAorb_u*uc2;
1854 REAL8 d2Aorb = 2.*dAorb_u*uc3 + d2Aorb_u*uc4;
1857 REAL8 AKerr_Multipole = (1.+2.*uc)/(1.+2.*
u);
1860 *A = Aorb*AKerr_Multipole*fss;
1861 *dA = dAorb*drc*(1.+2.*uc)/(1.+2.*
u) - 2.*Aorb*drc*uc2/(1.+2.*
u) + 2.*Aorb*(1.+2.*uc)*
u2/((1.+2.*
u)*(1.+2.*
u));
1862 *d2A = d2Aorb*(1.+2.*uc)/(1.+2.*
u) + 4.*dAorb*(
u2*(1.+2.*uc)/((1.+2.*
u)*(1.+2.*
u)) - uc2/(1.+2.*
u)*drc) + Aorb*(-4.*
u3*(1.+2.*uc)/((1.+2.*
u)*(1.+2.*
u)) + 8.*
u4*(1.+2.*uc)/((1.+2.*
u)*(1.+2.*
u)*(1.+2.*
u))+4.*uc3*(1.+2.*
u)*drc*drc - 2.*uc2/(1.+2.*
u)*d2rc);
1865 REAL8 Dp = 1.0 + 6.*nu*uc2 - 2.*(3.0*nu-26.0)*nu*uc3;
1867 REAL8 dD = 6.*uc2*(2.*nu*uc-(3.*nu-26.)*nu*uc2)*D*D;
1870 *
B =
r*
r*uc2*D/(*A);
1871 *dB = (dD*(*A) - D*(*dA))/((*A)*(*A));
1882 const REAL8 nu2 = nu*nu;
1883 const REAL8 nu3 = nu2*nu;
1893 const REAL8 sp2 = 1.-4.*nu;
1894 const REAL8 sp4 = (1-4*nu)*
SQ((1-2*nu));
1895 const REAL8 sp3 = (1.-3.*nu)*(1.-3.*nu);
1896 const REAL8 sp5 = (1.-5.*nu+5.*nu2)*(1.-5.*nu+5.*nu2);
1897 const REAL8 sp6 = (1-4*nu)*(3*nu2-4*nu +1)*(3*nu2-4*nu +1);
1898 const REAL8 sp7 = (1 - 7*nu + 14*nu2 - 7*nu3)*(1 - 7*nu + 14*nu2 - 7*nu3);
1899 const REAL8 sp8 = (1 - 4*nu)*(1 - 6*nu + 10*nu2 - 4*nu3)*(1 - 6*nu + 10*nu2 - 4*nu3);
1903 sp2 * x6, sp3 * x7, sp2 * x6,
1904 sp4 * x8, sp3 * x7, sp4 * x8, sp3 * x7,
1905 sp4 * x8, sp5 * x9, sp4 * x8, sp5 * x9, sp4 * x8,
1906 sp6 * x10, sp5 * x9, sp6 * x10, sp5 * x9, sp6 * x10, sp5 * x9,
1907 sp6 * x10, sp7 * x11, sp6 * x10, sp7 * x11, sp6 * x10, sp7 * x11, sp6 * x10,
1908 sp8 * x12, sp7 * x11, sp8 * x12, sp7 * x11, sp8 * x12, sp7 * x11, sp8 * x12, (7*nu3-14*nu2+7*nu-1)*(7*nu3-14*nu2+7*nu-1) * x11
1912 for (
int k = 0; k <
KMAX; k++) {
1913 Nlm[k] =
CNlm[k] * spx[k];
1923 720., 5040., 40320.,
1924 362880., 3628800., 39916800.,
1925 479001600., 6227020800., 87178291200.};
1931 REAL8 hhatk, x2,
y, prod;
1933 for (k = 0; k <
KMAX; k++ ) {
1935 x2 = 4.*hhatk*hhatk;
1938 prod *= ( j*j + x2 );
1941 y /= ( 1. - exp(-
y) );
1969 const REAL8 FNewt22 = 32./5.*x5;
1972 FlmHLO[1] = 32./5.*(1-4*nu+2*nu2)*x9;
1973 FlmHLO[0] = 32./5.*(1-4*nu+2*nu2)*x10;
1986 c1[1] = (4.-21.*nu + 27.*nu2 - 8.*nu3)/(4.*(1.-4.*nu+2.*nu2));
1991 rhoHlm[1] = 1. +
c1[1]*
x +
c2[1]*x2 + c3[1]*x3 + c4[1]*x4;
1992 rhoHlm[0] = 1. +
c1[0]*
x +
c2[0]*x2 + c3[0]*x3 + c4[0]*x4;
1995 const REAL8 Heff2 = Heff*Heff;
1996 const REAL8 jhat2 = jhat*jhat;
1998 FlmH[1] = FlmHLO[1] * Heff2 * gsl_pow_int(rhoHlm[1],4);
1999 FlmH[0] = FlmHLO[0] * jhat2 * gsl_pow_int(rhoHlm[0],4);
2002 REAL8 hatFH = (FlmH[0]+FlmH[1])/FNewt22;
2023 REAL8 v5 = sqrt(x5);
2029 cv5[0] = -1./4.*chi1*(1.+3.*chi1*chi1)*X1*X1*X1;
2030 cv5[1] = -1./4.*chi2*(1.+3.*chi2*chi2)*X2*X2*X2;
2033 cv8[0] = 0.5*(1.+sqrt(1.-chi1*chi1))*(1.+3.*chi1*chi1)*X1*X1*X1*X1;
2034 cv8[1] = 0.5*(1.+sqrt(1.-chi2*chi2))*(1.+3.*chi2*chi2)*X2*X2*X2*X2;
2036 REAL8 FH22_S = (cv5[0]+cv5[1])*v5;
2037 REAL8 FH22 = (cv8[0]+cv8[1])*x4;
2041 REAL8 hatFH = FH22_S + FH22 + FH21;
2069 const REAL8 X12 = X1-X2;
2070 const REAL8 UNUSED X12sq =
SQ(X12);
2078 jhat, Heff, jhat, Heff,
2079 Heff, jhat, Heff, jhat, Heff,
2080 jhat, Heff, jhat, Heff, jhat, Heff,
2081 Heff, jhat, Heff, jhat, Heff, jhat, Heff,
2082 jhat, Heff, jhat, Heff, jhat, Heff, jhat, Heff};
2084 REAL8 FNewt22, sum_k=0.;
2094 REAL8 x6 = gsl_pow_int(
x, 6);
2095 FNewtlm[0] =
CNlm[0] * x6;
2096 FNewtlm[2] =
CNlm[2] * x6;
2097 FNewtlm[4] =
CNlm[4] * x6;
2099 REAL8 sp4x8 =
SQ((1-2*nu)) * gsl_pow_int(
x, 8);
2100 FNewtlm[5] =
CNlm[5] * sp4x8;
2101 FNewtlm[7] =
CNlm[7] * sp4x8;
2105 REAL8 x6 = gsl_pow_int(
x, 6);
2106 FNewtlm[0] =
CNlm[0] * x6;
2107 FNewtlm[2] =
CNlm[2] * x6;
2108 FNewtlm[4] =
CNlm[4] * x6;
2118 dyn->
eob_wav_flm_s(
x, nu, X1, X2, chi1, chi2, a1,
a2, C_Q1, C_Q2,
dyn->
clm, usetidal, rholm, flm);
2124 FNewt22 = FNewtlm[1];
2127 for (
int k = 0; k <
KMAX; k++) hlmNQC[k] = 1.;
2143 hlmNQC[1] = hNQC.
ampli[1];
2147 for (
int k = 0; k <
KMAX; k++) {
2148 Modhhatlm[k] = prefact[k] * MTlm[k] * flm[k] * hlmNQC[k];
2156 Modhhatlm[0] *= X12;
2157 Modhhatlm[2] *= X12;
2158 Modhhatlm[4] *= X12;
2161 for (
int k = 0; k <
KMAX; k++) {
2162 Modhhatlm[k] += MTlm[k] * hlmTidal[k];
2167 for (
int k =
KMAX; k--;) sum_k +=
SQ(Modhhatlm[k]) * FNewtlm[k];
2170 REAL8 hatf = sum_k/(FNewt22);
2184 return (-32./5. * nu * gsl_pow_int(r_omega,4) * gsl_pow_int(Omega,5) * hatf);
2219 0.1113090643348557, 0.1112593821157397, 0.0424759238428813, 0.0424489015884926, 0.0424717446800903, 0.0215953972500844, 0.0215873812155663, 0.0215776183122621, 0.0216017621863542, 0.0128123696874894, 0.0128097056242375, 0.0128038943888768, 0.0128025242617949, 0.0128202485907368, 0.0083762045692408, 0.0083751913886140, 0.0083724067460769, 0.0083694435961860, 0.0083710364141552, 0.0083834483913443, 0.0058540393221396, 0.0058536069384738, 0.0058522594457692, 0.0058502436535615, 0.0058491157293566, 0.0058514875071582, 0.0058602498033381, 0.0042956812356573, 0.0042954784390887, 0.0042947951664056, 0.0042935886137697, 0.0042923691461384, 0.0042922256848799, 0.0042945927126022, 0.0043009106861259};
2222 0.0004643273300862, 0.0009375605440004, 0.0000597134489198, 0.0002551406918111, 0.0001741036904709, 0.0000124649041611, 0.0000685496215625, 0.0001131160409390, 0.0000419907542591, 0.0000035218982282, 0.0000219211271097, 0.0000473186962874, 0.0000524142634057, 0.0000106823372552, 0.0000012237574387, 0.0000081742188269, 0.0000201940563214, 0.0000295722761753, 0.0000260539631956, 0.0000018994753518, 0.0000004932942990, 0.0000034477210351, 0.0000092294406360, 0.0000155143073237, 0.0000183386499818, 0.0000137922469695, -0.0000007075155453, 0.0000002223410995, 0.0000016045317657, 0.0000045260028113, 0.0000082655700107, 0.0000112393599417, 0.0000115758243113, 0.0000076838709956, -0.0000014020591745};
2224 const REAL8 b3[] = {
2225 -0.0221835462237291, -0.0235386333304348, -0.0042911639711832, -0.0047431560217121, -0.0046577314472149, -0.0013089557502947, -0.0014343968205390, -0.0014978542575474, -0.0014329302934532, -0.0005167994164556, -0.0005573939123058, -0.0005921030407223, -0.0005978284714483, -0.0005673965369076, -0.0002409269302708, -0.0002561516055118, -0.0002723768586352, -0.0002815958312453, -0.0002792078156272, -0.0002646630240693, -0.0001261183503407, -0.0001325622938779, -0.0001403198638518, -0.0001464084186977, -0.0001485971591029, -0.0001459023931717, -0.0001384829633836, -0.0000719062974278, -0.0000749128468013, -0.0000788187384314, -0.0000824202283094, -0.0000846673495936, -0.0000849054394951, -0.0000829269749240, -0.0000788883333858};
2227 const REAL8 b4[] = {
2228 0.0058366730167965, 0.0070452306758401, 0.0006914295465364, 0.0010322294603561, 0.0010057563135650, 0.0001394203795507, 0.0002309706405978, 0.0002596611624417, 0.0002409588083156, 0.0000386949167221, 0.0000679154947896, 0.0000830199015202, 0.0000850120755064, 0.0000780125513602, 0.0000133034384660, 0.0000241813441339, 0.0000311573885555, 0.0000340233089866, 0.0000335167900637, 0.0000307571022927, 0.0000053305073331, 0.0000099143129290, 0.0000132296989826, 0.0000150959309402, 0.0000156304390748, 0.0000151274875147, 0.0000139320508803, 0.0000023959090314, 0.0000045285807761, 0.0000061918979830, 0.0000072894226381, 0.0000078251853305, 0.0000078772667984, 0.0000075606242809, 0.0000069956215270
2235 const REAL8 psi[] = {0.9227843350984671394, 0.9227843350984671394,
2236 1.256117668431800473, 1.256117668431800473, 1.256117668431800473,
2237 1.506117668431800473, 1.506117668431800473, 1.506117668431800473, 1.506117668431800473,
2238 1.706117668431800473, 1.706117668431800473, 1.706117668431800473, 1.706117668431800473, 1.706117668431800473,
2239 1.872784335098467139, 1.872784335098467139, 1.872784335098467139, 1.872784335098467139, 1.872784335098467139, 1.872784335098467139,
2240 2.015641477955609997, 2.015641477955609997, 2.015641477955609997, 2.015641477955609997, 2.015641477955609997, 2.015641477955609997, 2.015641477955609997,
2241 2.140641477955609997, 2.140641477955609997, 2.140641477955609997, 2.140641477955609997, 2.140641477955609997, 2.140641477955609997, 2.140641477955609997, 2.140641477955609997};
2253 num_ang = 1. +
b1[
i]*x2 +
b2[
i]*x3 + b3[
i]*x4 + b4[
i]*x5;
2254 tlm_ang = (- 2. * psi[
i] *
x * num_ang) + 2.*
x* log(2. * k * bphys);
2268 gsl_sf_result num_rad;
2269 gsl_sf_result num_phase;
2270 gsl_sf_result denom_rad;
2271 gsl_sf_result denom_phase;
2283 gsl_sf_lngamma_complex_e((
double)
TEOB_LINDEX[
i] + 1., (
double) (-2.*hhatk), &num_rad, &num_phase);
2284 gsl_sf_lngamma_complex_e((
double)
TEOB_LINDEX[
i] + 1., 0., &denom_rad, &denom_phase);
2286 ratio_rad = (
REAL8) (num_rad.val - denom_rad.val);
2287 ratio_ang = (
REAL8) num_phase.val - 0.;
2289 tlm_rad = ratio_rad +
LAL_PI * hhatk;
2290 tlm_phase = ratio_ang + 2.*hhatk*log(2.*k*bphys);
2292 tlm->
ampli[
i] = exp(tlm_rad);
2293 tlm->
phase[
i] = tlm_phase;
2344 const REAL8 X12 = X1-X2;
2345 const REAL8 a0X12 = a0*X12;
2346 const REAL8 a12X12 = a12*X12;
2350 const REAL8 v3 = v*v2;
2351 const REAL8 v4 = v3*v;
2352 const REAL8 v5 = v4*v;
2356 const REAL8 cSO_lo = (-0.5*a0 - a12X12/6.);
2357 const REAL8 cSO_nlo = (-52./63.-19./504.*nu)*a0 - (50./63.+209./504.*nu)*a12X12;
2362 #if (RC_EXCLUDESPINSPINTIDES)
2372 cSS_lo = 0.5*(C_Q1*a1*a1 + 2.*a1*
a2 + C_Q2*
a2*
a2);
2379 rho22S = cSO_lo*v3 + cSS_lo*v4 + cSO_nlo*v5 ;
2382 rho32S = (a0-a12X12)/(3.*(1.-3.*nu))*v;
2383 rho44S = (-19./30.*a0 - (1.-21.*nu)/(30.-90.*nu)*a12X12)*v3;
2384 rho42S = ( -1./30.*a0 - (19.-39.*nu)/(30.-90.*nu)*a12X12)*v3;
2387 f21S = -1.5*a12*v + ((110./21. + 79./84.*nu)*a12 - 13./84.*a0X12)*v3;
2388 f33S = ((-0.25 + 2.5*nu)*a12 - 1.75*a0X12)*v3;
2389 f31S = ((-2.25 + 6.5*nu)*a12 + 0.25*a0X12)*v3;
2390 f43S = (( 5. -10.*nu)*a12 - 5.*a0X12)/(-4.+8.*nu)*v;
2394 flm[0] =
SQ(rholm[0]);
2395 flm[0] = (X12*flm[0] + f21S);
2397 flm[1] =
SQ(rholm[1]+ rho22S);
2399 flm[2] = gsl_pow_int(rholm[2], 3);
2400 flm[2] = (X12*flm[2] + f31S);
2402 flm[3] = gsl_pow_int(rholm[3]+ rho32S, 3);
2404 flm[4] = gsl_pow_int(rholm[4], 3);
2405 flm[4] = (X12*flm[4] + f33S);
2407 flm[5] = gsl_pow_int(rholm[5], 4);
2408 flm[5] = (X12*flm[5] + f41S);
2410 flm[6] = gsl_pow_int(rholm[6] + rho42S, 4);
2412 flm[7] = gsl_pow_int(rholm[7], 4);
2413 flm[7] = (X12*flm[7] + f43S);
2415 flm[8] = gsl_pow_int(rholm[8] + rho44S, 4);
2463 const REAL8 X12 = X1-X2;
2464 const REAL8 a0X12 = a0*X12;
2465 const REAL8 a12X12 = a12*X12;
2469 const REAL8 v3 = v2*v;
2470 const REAL8 v4 = v3*v;
2471 const REAL8 v5 = v4*v;
2472 const REAL8 v6 = v5*v;
2473 const REAL8 UNUSED v7 = v6*v;
2477 const REAL8 cSO_lo = (-0.5*a0 - a12X12/6.);
2478 const REAL8 cSO_nlo = (-52./63.-19./504.*nu)*a0 - (50./63.+209./504.*nu)*a12X12;
2479 const REAL8 UNUSED cSO_nnlo = (32873./21168 + 477563./42336.*nu + 147421./84672.*nu*nu)*a0 - (23687./63504 - 171791./127008.*nu + 50803./254016.*nu*nu)*a12X12;
2485 #if (RC_EXCLUDESPINSPINTIDES)
2496 cSS_lo = 0.5*(C_Q1*a1*a1 + 2.*a1*
a2 + C_Q2*
a2*
a2);
2497 cSS_nlo = (-85./63. + 383./252.*nu)*a1*
a2 + (-2./3. - 5./18.*nu)*(a1*a1 +
a2*
a2) + (1./7. + 27./56.*nu)*(C_Q1*a1*a1 + C_Q2*
a2*
a2) + 2./9.*X12*(a1*a1 -
a2*
a2) + 55./84.*X12*(C_Q1*a1*a1 - C_Q2*
a2*
a2);
2501 cSS_nlo = 1./504.*(2.*(19. - 70.*nu)*a12*a12 + (-302. + 243.*nu)*a0*a0 + 442.*X12*a0*a12);
2505 rho22S = cSO_lo*v3 + cSS_lo*v4 + cSO_nlo*v5;
2508 rho22S += cSS_nlo*v6;
2511 rho32S = (a0-a12X12)/(3.*(1.-3.*nu))*v;
2512 rho44S = (-19./30.*a0 - (1.-21.*nu)/(30.-90.*nu)*a12X12)*v3;
2513 rho42S = ( -1./30.*a0 - (19.-39.*nu)/(30.-90.*nu)*a12X12)*v3;
2517 f21S = -1.5*a12*v + ((110./21. + 79./84.*nu)*a12 - 13./84.*a0X12)*v3;
2518 f33S = ((-0.25 + 2.5*nu)*a12 - 1.75*a0X12)*v3;
2519 f31S = ((-2.25 + 6.5*nu)*a12 + 0.25*a0X12)*v3;
2520 f43S = (( 5. -10.*nu)*a12 - 5.*a0X12)/(-4.+8.*nu)*v;
2528 #if (RC_EXCLUDESPINSPINTIDES)
2537 c21SS_lo = -19./8.*(a1*a1 -
a2*
a2) - (C_Q1*a1*a1 - C_Q2*
a2*
a2) + 1./8.*(-9.*a1*a1 + 10*a1*
a2 -9.*
a2*
a2 + 12.*(C_Q1*a1*a1 + C_Q2*
a2*
a2))*X12;
2538 c33SS_lo = 3.*(a1*
a2 + 0.5*(C_Q1*a1*a1 + C_Q2*
a2*
a2))*X12;
2539 c31SS_lo = -4.*(C_Q1*a1*a1 - C_Q2*
a2*
a2) + 3.*(a1*
a2 + 0.5*(C_Q1*a1*a1 + C_Q2*
a2*
a2))*X12;
2542 c21SS_lo = 1./8.*(-27.*(a1*a1 -
a2*
a2) + (3.*a1*a1 + 10.*a1*
a2 + 3.*
a2*
a2)*X12);
2543 c33SS_lo = 3./2.*a0*a0*X12;
2544 c31SS_lo = -4.*(a1*a1 -
a2*
a2) + 3./2.*a0*a0*X12;
2548 f21S += c21SS_lo*v4;
2549 f33S += c33SS_lo*v4;
2550 f31S += c31SS_lo*v4;
2553 flm[0] =
SQ(rholm[0]);
2554 flm[0] = (X12*flm[0] + f21S);
2556 flm[1] =
SQ(rholm[1]+ rho22S);
2558 flm[2] = gsl_pow_int(rholm[2], 3);
2559 flm[2] = (X12*flm[2] + f31S);
2561 flm[3] = gsl_pow_int(rholm[3]+ rho32S, 3);
2563 flm[4] = gsl_pow_int(rholm[4], 3);
2564 flm[4] = (X12*flm[4] + f33S);
2566 flm[5] = gsl_pow_int(rholm[5], 4);
2567 flm[5] = (X12*flm[5] + f41S);
2569 flm[6] = gsl_pow_int(rholm[6] + rho42S, 4);
2571 flm[7] = gsl_pow_int(rholm[7], 4);
2572 flm[7] = (X12*flm[7] + f43S);
2574 flm[8] = gsl_pow_int(rholm[8] + rho44S, 4);
2587 const REAL8 nu2 = nu*nu;
2588 const REAL8 nu3 = nu*nu2;
2589 const REAL8 nu4 = nu*nu3;
2591 for (
int k=0; k<
KMAX; k++) clm[k][0] = 1.;
2592 for (
int k=0; k<
KMAX; k++)
for (
int n=1; n<6; n++) clm[k][n] = 0.;
2595 clm[0][1] = (-1.0535714285714286 + 0.27380952380952384 *nu);
2596 clm[0][2] = (-0.8327841553287982 - 0.7789824263038548 *nu + 0.13116496598639457*nu2);
2602 clm[1][1] = (-1.0238095238095237 + 0.6547619047619048*nu);
2603 clm[1][2] = (-1.94208238851096 - 1.5601379440665155*nu + 0.4625614134542706*nu2);
2609 clm[2][1] = (-0.7222222222222222 - 0.2222222222222222*nu);
2610 clm[2][2] = (0.014169472502805836 - 0.9455667789001122*nu - 0.46520763187429853*nu2);
2616 clm[3][1] = (0.003703703703703704*(328. - 1115.*nu + 320.*nu2))/(-1. + 3.*nu);
2617 clm[3][2] = (6.235191420376606e-7*(-1.444528e6 + 8.050045e6*nu - 4.725605e6*nu2 - 2.033896e7*nu3 + 3.08564e6*nu4))/((-1. + 3.*nu)*(-1. + 3.*nu));
2622 clm[4][1] = (-1.1666666666666667 + 0.6666666666666666*nu);
2623 clm[4][2] = (-1.6967171717171716 - 1.8797979797979798*nu + 0.45151515151515154*nu2);
2629 clm[5][1] = (0.001893939393939394*(602. - 1385.*nu + 288.*nu2))/(-1. + 2.*nu);
2630 clm[5][2] = (- 0.36778992787515513);
2635 clm[6][1] = (0.0007575757575757576*(1146. - 3530.*nu + 285.*nu2))/(-1. + 3.*nu);
2636 clm[6][2] = - (3.1534122443213353e-9*(1.14859044e8 - 2.95834536e8*nu - 1.204388696e9*nu2 + 3.04798116e9*nu3 + 3.79526805e8*nu4))/((-1. + 3.*nu)*(-1. + 3.*nu));
2641 clm[7][1] = (0.005681818181818182*(222. - 547.*nu + 160.*nu2))/(-1. + 2.*nu);
2642 clm[7][2] = (- 0.9783218202252293);
2647 clm[8][1] = (0.0007575757575757576*(1614. - 5870.*nu + 2625.*nu2))/(-1. + 3.*nu);
2648 clm[8][2] = (3.1534122443213353e-9*(-5.11573572e8 + 2.338945704e9*nu - 3.13857376e8*nu2 - 6.733146e9*nu3 + 1.252563795e9*nu4))/((-1. + 3.*nu)*(-1. + 3.*nu));
2653 clm[9][1] = (0.002564102564102564*(319. - 626.*nu + 8.*nu2))/(-1. + 2.*nu);
2654 clm[9][2] = (- 0.1047896120973044);
2659 clm[10][1] = (0.00007326007326007326*(-15828. + 84679.*nu - 104930.*nu2 + 21980.*nu3))/(1. - 5.*nu + 5.*nu2);
2664 clm[11][1] = (0.002564102564102564*(375. - 850.*nu + 176.*nu2))/(-1. + 2.*nu);
2665 clm[11][2] = (- 0.5788010707241477);
2670 clm[12][1] = (0.00007326007326007326*(-17448. + 96019.*nu - 127610.*nu2 + 33320.*nu3))/(1. - 5.*nu + 5.*nu2);
2675 clm[13][1] = (0.002564102564102564*(487. - 1298.*nu + 512.*nu2))/(-1. + 2.*nu);
2676 clm[13][2] = (- 1.5749727622804546);
2681 clm[14][1] = (0.006944444444444444*(-161. + 694.*nu - 670.*nu2 + 124.*nu3))/(1. - 4.*nu + 3.*nu2);
2686 clm[15][1] = (0.011904761904761904*(-74. + 378.*nu - 413.*nu2 + 49.*nu3))/(1. - 5.*nu + 5.*nu2);
2687 clm[15][2] = ( - 0.24797525070634313)*
PMTERMS_eps;
2691 clm[16][1] = (0.006944444444444444*(-169. + 742.*nu - 750.*nu2 + 156.*nu3))/(1. - 4.*nu + 3.*nu2);
2696 clm[17][1] = (0.011904761904761904*(-86. + 462.*nu - 581.*nu2 + 133.*nu3))/(1. - 5.*nu + 5.*nu2);
2701 clm[18][1] = (0.006944444444444444*(-185. + 838.*nu - 910.*nu2 + 220.*nu3))/(1. - 4.*nu + 3.*nu2);
2706 clm[19][1] = (0.011904761904761904*(-106. + 602.*nu - 861.*nu2 + 273.*nu3))/(1. - 5.*nu + 5.*nu2);
2711 clm[20][1] = (0.0014005602240896359*(-618. + 2518.*nu - 2083.*nu2 + 228.*nu3))/(1. - 4.*nu + 3.*nu2);
2716 clm[21][1] = (0.00006669334400426837*(16832. - 123489.*nu + 273924.*nu2 - 190239.*nu3 + 32760.*nu4))/(-1. + 7.*nu - 14.*nu2 + 7.*nu3);
2720 clm[22][1] = (0.0014005602240896359*(-666. + 2806.*nu - 2563.*nu2 + 420.*nu3))/(1. - 4.*nu + 3.*nu2);
2725 clm[23][1] = (0.00006669334400426837*(17756. - 131805.*nu + 298872.*nu2 - 217959.*nu3 + 41076.*nu4))/(-1. + 7.*nu - 14.*nu2 + 7.*nu3);
2729 clm[24][1] = (0.0014005602240896359*(-762. + 3382.*nu - 3523.*nu2 + 804.*nu3))/(1. - 4.*nu + 3.*nu2);
2734 clm[25][1] = (0.0006002400960384153*(2144. - 16185.*nu + 37828.*nu2 - 29351.*nu3 + 6104.*nu4))/(-1. + 7.*nu - 14.*nu2 + 7.*nu3);
2738 clm[26][1] = (0.0014005602240896359*(-906. + 4246.*nu - 4963.*nu2 + 1380.*nu3))/(1. - 4.*nu + 3.*nu2);
2743 clm[27][1] = (0.00005482456140350877*(20022. - 126451.*nu + 236922.*nu2 - 138430.*nu3 + 21640.*nu4))/(-1. + 6.*nu - 10.*nu2 + 4.*nu3);
2747 clm[28][1] = (0.0003654970760233918*(2462. - 17598.*nu + 37119.*nu2 - 22845.*nu3 + 3063.*nu4))/(-1. + 7.*nu - 14.*nu2 + 7.*nu3);
2751 clm[29][1] = (0.00005482456140350877*(20598. - 131059.*nu + 249018.*nu2 - 149950.*nu3 + 24520.*nu4))/(-1. + 6.*nu - 10.*nu2 + 4.*nu3);
2755 clm[30][1] = (0.0003654970760233918*(2666. - 19434.*nu + 42627.*nu2 - 28965.*nu3 + 4899.*nu4))/(-1. + 7.*nu - 14.*nu2 + 7.*nu3);
2759 clm[31][1] = (0.00027412280701754384*(4350. - 28055.*nu + 54642.*nu2 - 34598.*nu3 + 6056.*nu4))/(-1. + 6.*nu - 10.*nu2 + 4.*nu3);
2763 clm[32][1] = (0.0010964912280701754*(1002. - 7498.*nu + 17269.*nu2 - 13055.*nu3 + 2653.*nu4))/(-1. + 7.*nu - 14.*nu2 + 7.*nu3);
2767 clm[33][1] = (0.00005482456140350877*(23478. - 154099.*nu + 309498.*nu2 - 207550.*nu3 + 38920.*nu4))/(-1. + 6.*nu - 10.*nu2 + 4.*nu3);
2771 clm[34][1] = (0.0003654970760233918*(3482. - 26778.*nu + 64659.*nu2 - 53445.*nu3 + 12243.*nu4))/(-1. + 7.*nu - 14.*nu2 + 7.*nu3);
2788 const REAL8 nu2 = nu*nu;
2789 const REAL8 nu3 = nu*nu2;
2803 clm[0][3] = (2.9192806270460925 - 1.019047619047619 *el1);
2804 clm[0][4] = (-1.28235780892213 + 1.073639455782313 *el1);
2805 clm[0][5] = (-3.8466571723355227 + 0.8486467106683944 *el1)*
PMTERMS_eps;
2807 clm[1][3] = (12.736034731834051 - 2.902228713904598 *nu - 1.9301558466099282*nu2 + 0.2715020968103451*nu3 - 4.076190476190476*el2);
2808 clm[1][4] = (-2.4172313935587004 + 4.173242630385488 *el2);
2809 clm[1][5] = (-30.14143102836864 + 7.916297736025627 *el2);
2811 clm[2][3] = (1.9098284139598072 - 0.4126984126984127*el1+ (-4.646868015386534 + (0.21354166666666666)*Pi2)*nu + 2.3020866307903347*nu2 - 0.5813492634480288*nu3);
2812 clm[2][4] = (0.5368150316615179 + 0.2980599647266314*el1);
2813 clm[2][5] = (1.4497991763035063 - 0.0058477188106817735*el1)*
PMTERMS_eps;
2815 clm[3][3] = (6.220997955214429 - 1.6507936507936507*el2);
2816 clm[3][4] = (-3.4527288879001268 + 2.005408583186361*el2)*
PMTERMS_eps;
2818 clm[4][3] = (14.10891386831863 - 3.7142857142857144*el3 + (-5.031429681429682 + (0.21354166666666666)*Pi2)*nu - 1.7781727531727531*nu2 + 0.25923767590434255*nu3);
2819 clm[4][4] = (-6.723375314944128 + 4.333333333333333*el3);
2820 clm[4][5] = (-29.568699895427518 + 6.302092352092352*el3)*
PMTERMS_eps;
2822 clm[5][3] = (0.6981550175535535 - 0.2266955266955267*el1);
2823 clm[5][4] = (-0.7931524512893319 + 0.2584672482399755*el1)*
PMTERMS_eps;
2825 clm[6][3] = 4.550378418934105e-12*(8.48238724511e11 - 1.9927619712e11*el2);
2826 clm[6][4] = (-0.6621921297263365 + 0.787251738160829*el2)*
PMTERMS_eps;
2828 clm[7][3] = (8.519456157072423 - 2.0402597402597404*el3)*
PMTERMS_eps;
2829 clm[7][4] = (-5.353216984886716 + 2.5735094451003544*el3)*
PMTERMS_eps;
2831 clm[8][3] = (15.108111214795123 - 3.627128427128427*el4);
2832 clm[8][4] = (-8.857121657199649 + 4.434988849534304*el4)*
PMTERMS_eps;
2834 clm[9][3] = (0.642701885362399 - 0.14414918414918415*el1)*
PMTERMS_eps;
2835 clm[9][4] = (-0.07651588046467575 + 0.11790664036817883*el1)*
PMTERMS_eps;
2837 clm[10][3] = (2.354458371550237 - 0.5765967365967366*el2)*
PMTERMS_eps;
2839 clm[11][3] = (5.733973288504755 - 1.2973426573426574*el3)*
PMTERMS_eps;
2840 clm[11][4] = (-1.9573287625526001 + 1.2474448628294783*el3)*
PMTERMS_eps;
2842 clm[12][3] = (10.252052781721588 - 2.3063869463869464*el4)*
PMTERMS_eps;
2844 clm[13][3] = (15.939827047208668 - 3.6037296037296036*el5)*
PMTERMS_eps;
2845 clm[13][4] = (-10.272578060123237 + 4.500041838503377*el5)*
PMTERMS_eps;
2847 clm[14][3] = (0.21653486654395454 - 0.10001110001110002*el1)*
PMTERMS_eps;
2849 clm[15][3] = (1.7942694138754138 - 0.40004440004440006*el2)*
PMTERMS_eps;
2851 clm[16][3] = (4.002558222882566 - 0.9000999000999002*el3)*
PMTERMS_eps;
2853 clm[17][3] = (7.359388663371044 - 1.6001776001776002*el4)*
PMTERMS_eps;
2855 clm[18][3] = (11.623366217471297 - 2.5002775002775004*el5)*
PMTERMS_eps;
2857 clm[19][3] = (16.645950799433503 - 3.6003996003996006*el6)*
PMTERMS_eps;
2859 clm[20][3] = (0.2581280702019663 - 0.07355557607658449*el1)*
PMTERMS_eps;
2861 clm[22][3] = (3.0835293524055283 - 0.6620001846892604*el3)*
PMTERMS_eps;
2863 clm[24][3] = (8.750589067052443 - 1.838889401914612*el5)*
PMTERMS_eps;
2865 clm[26][3] = (17.255875091408523 - 3.6042232277526396*el7)*
PMTERMS_eps;
2872 const REAL8 xn[] = {1.,
x,x2,x3,x4,x5};
2874 for (
int k=0; k<
KMAX; k++) {
2877 rholm[k] = clm[k][0];
2878 for (
int n=1; n<6; n++) {
2879 rholm[k] += clm[k][n]*xn[n];
2882 rholm[k] = x5*clm[k][5];
2883 for (
int n=5; n-- >1; ) {
2884 rholm[k] += clm[k][n]*xn[n];
2886 rholm[k] += clm[k][0];
2893 for (
int k = 0; k <
KMAX; k++) {
2920 REAL8 delta22LO = 7./3. * y32;
2921 REAL8 delta21LO = 2./3. * y32;
2922 REAL8 delta33LO = 13./10. * y32;
2923 REAL8 delta31LO = 13./30. * y32;
2926 for (
int k = 0; k <
KMAX; k++) {
2936 num = 69020.*nu + 5992.*
LAL_PI*sqrt_y;
2937 den = 5992.*
LAL_PI*sqrt_y + 2456.*nu*(28.+493.*nu*
y);
2938 dlm[0] = delta21LO*num/den;
2940 num = (808920.*nu*
LAL_PI*sqrt(
y) + 137388.*Pi2*
y + 35.*nu2*(136080. + (154975. - 1359276.*nu)*
y));
2941 den = (808920.*nu*
LAL_PI*sqrt(
y) + 137388.*Pi2*
y + 35.*nu2*(136080. + (154975. + 40404.*nu)*
y));
2942 dlm[1] = delta22LO*num/den;
2946 num = 4641.*nu + 1690.*
LAL_PI*sqrt_y;
2947 den = num + 18207.*nu2*
y;
2948 dlm[2] = delta31LO*num/den;
2950 num = 1. + 94770.*
LAL_PI/(566279.*nu)*sqrt_y;
2951 den = num + 80897.* nu/3159.*
y;
2952 dlm[3] = (10.+33.*nu)/(15.*(1.-3.*nu)) * y32 + 52./21.*
LAL_PI*y3;
2954 dlm[4] = delta33LO*num/den;
2957 dlm[5] = (2.+507.*nu)/(10.*(1.-2.*nu))*y32 + 1571./3465.*
LAL_PI*y3;
2958 dlm[6] = 7.*(1.+6.*nu)/(15.*(1.-3.*nu))*y32 + 6284./3465.*
LAL_PI*y3;
2959 dlm[7] = (486.+4961.*nu)/(810.*(1.-2.*nu))*y32 + 1571./385.*
LAL_PI*y3;
2960 dlm[8] = (112.+219.*nu)/(120.*(1.-3.*nu))*y32 + 25136./3465.*
LAL_PI*y3;
2963 dlm[9] = (96875. + 857528.*nu)/(131250.*(1.-2.*nu))*y32;
2982 REAL8 vphi2 = vphi*vphi;
2983 REAL8 vphi3 = vphi*vphi2;
2984 REAL8 vphi4 = vphi*vphi3;
2985 REAL8 vphi5 = vphi*vphi4;
2986 REAL8 vphi6 = vphi*vphi5;
2987 REAL8 vphi7 = vphi*vphi6;
2988 REAL8 vphi8 = vphi*vphi7;
2989 REAL8 vphi9 = vphi*vphi8;
2992 const REAL8 p1 = 1.;
2993 const REAL8 p2 = sqrt(1.-4.*nu);
2994 const REAL8 p3 = (3.*nu-1.);
2995 const REAL8 p4 = (2.*nu-1.)*sqrt(1.-4.*nu);
2996 const REAL8 p5 = 1.-5.*nu+5.*nu2;
2997 const REAL8 p6 = (1.-4.*nu+3.*nu2)*sqrt(1.-4.*nu);
2998 const REAL8 p7 = 7.*nu3 - 14.*nu2 + 7.*nu -1.;
2999 const REAL8 p8 = (4.*nu3 - 10.*nu2 + 6.*nu -1.)*sqrt(1.-4.*nu);
3001 const REAL8 phix2 = 2. * phi;
3002 const REAL8 phix3 = 3. * phi;
3003 const REAL8 phix4 = 4. * phi;
3004 const REAL8 phix5 = 5. * phi;
3005 const REAL8 phix6 = 6. * phi;
3006 const REAL8 phix7 = 7. * phi;
3008 const REAL8 pv23 = p2 * vphi3;
3009 const REAL8 pv34 = p3 * vphi4;
3010 const REAL8 pv45 = p4 * vphi5;
3011 const REAL8 pv56 = p5 * vphi6;
3012 const REAL8 pv67 = p6 * vphi7;
3013 const REAL8 pv78 = p7 * vphi8;
3014 const REAL8 pv89 = p8 * vphi9;
3019 phi,phix2,phix3,phix4,
3020 phi,phix2,phix3,phix4,phix5,
3021 phi,phix2,phix3,phix4,phix5,phix6,
3022 phi,phix2,phix3,phix4,phix5,phix6,phix7,
3023 phi,phix2,phix3,phix4,phix5,phix6,phix7,8.*phi
3029 pv45, pv34, pv45, pv34,
3030 pv45, pv56, pv45, pv56, pv45,
3031 pv67, pv56, pv67, pv56, pv67, pv56,
3032 pv67, pv78, pv67, pv78, pv67, pv78, pv67,
3033 pv89, pv78, pv89, pv78, pv89, pv78, pv89, pv78
3037 for (
int k = 0; k <
KMAX; k++) {
3058 const REAL8 x5 = (
REAL8) gsl_pow_int((
double)
x, 5);
3059 const REAL8 x6 = (
REAL8) gsl_pow_int((
double)
x, 6);
3063 memset(hTidallm, 0.,
KMAX*
sizeof(
REAL8));
3070 hA[1] = 2 * khatA_2 *(XA/XB+3);
3071 hB[1] = 2 * khatB_2 *(XB/XA+3);
3073 betaA1[1] = (-202. + 560*XA - 340*XA*XA + 45*XA*XA*XA)/(42*(3-2*XA));
3074 betaB1[1] = (-202. + 560*XB - 340*XB*XB + 45*XB*XB*XB)/(42*(3-2*XB));
3076 hA[0] = 3 * khatA_2 * (3-4*XA);
3077 hB[0] = 3 * khatB_2 * (3-4*XB);
3080 hA[2] = 12 * khatA_2 * XB;
3081 hB[2] = 12 * khatB_2 * XA;
3083 betaA1[2] = (-6. -5.*XA +131.*XA*XA -130.*XA*XA*XA)/(36.*(1.-XA));
3084 betaB1[2] = (-6. -5.*XB +131.*XB*XB -130.*XB*XB*XB)/(36.*(1.-XB));
3089 betaA1[4] = ( (XA-3.)*(10.*XA*XA - 25.*XA+ 14.) )/(12.*(1.-XA));
3090 betaB1[4] = ( (XB-3.)*(10.*XB*XB - 25.*XB+ 14.) )/(12.*(1.-XB));
3094 hTidallm[0] = ( -hA[0] + hB[0] )*x5;
3096 hTidallm[1] = ( hA[1]*(1. + betaA1[1]*
x) + hB[1]*(1. + betaB1[1]*
x) )*x5;
3100 hTidallm[2] = ( -hA[2]*(1. + betaA1[2]*
x) + hB[2]*(1. + betaB1[2]*
x) )*x5;
3102 hTidallm[3] = 8.*( khatA_2*(1. -2.*XB + 3.*XB*XB) +khatB_2*(1. -2.*XA + 3.*XA*XA) )*x5/(1.-3.*nu);
3104 hTidallm[4] = ( -hA[4]*(1. + betaA1[4]*
x) + hB[4]*(1. + betaB1[4]*
x) )*x5;
3109 const double fourtnine= 1.5555555555555555556;
3110 const double fourthird = 1.3333333333333333333;
3111 hTidallm[0] += 0.5*( -1.*kapA2j/XB + kapB2j/XA )*x5;
3112 hTidallm[1] += fourtnine*kapT2j*x6;
3113 hTidallm[2] += 0.5*( kapA2j*(4. - 17.*XB) - kapB2j*(4. - 17.*XA) )*x6;
3114 hTidallm[3] += fourthird*kapT2j*x5/(1.-3.*nu);
3115 hTidallm[4] += 0.5*( kapA2j*(4. - 9.*XB) - kapB2j*(4. - 9.*XA) )*x6;
3130 const REAL8 xnu = 1-4*nu;
3143 for (
int ki = 0; ki <
KMAX; ki++) {
3153 a1[
k21] = 0.0162387198*(7.32653082*xnu2 + 1.19616248*xnu + 0.73496656);
3154 a2[
k21] = -1.80492460*xnu2 + 1.78172686*xnu + 0.30865284;
3157 b1[
k21] = -0.0647955017*(3.59934444*xnu2 - 4.08628784*xnu + 1.37890907);
3158 b2[
k21] = 1.3410693180*(0.38491989*xnu2 + 0.10969453*xnu + 0.97513971);
3162 a1[
k22] = -0.0805236959*( 1 - 2.00332326*xnu2)/( 1 + 3.08595088*xnu2);
3163 a2[
k22] = 1.5299534255*( 1 + 1.16438929*xnu2)/( 1 + 1.92033923*xnu2);
3166 b1[
k22] = 0.146768094955*( 0.07417121*xnu + 1.01691256);
3167 b2[
k22] = 0.896911234248*(-0.61072011*xnu + 0.94295129);
3171 a1[
k33] = -0.0377680000*(1 - 14.61548907*xnu2)/( 1 + 2.44559263*xnu2);
3172 a2[
k33] = 1.9898000000*(1 + 2.09750346 *xnu2)/( 1 + 2.57489466*xnu2);
3175 b1[
k33] = 0.1418400000*(1.07430512 - 1.23906804*xnu + 4.44910652*xnu2);
3176 b2[
k33] = 0.6191300000*(0.80672432 + 4.07432829*xnu - 7.47270977*xnu2);
3180 for (
int ki = 0; ki <
KMAX; ki++) {
3181 for (
int j = 0; j < 6; j++) {
3187 n[k][0] = (prstar/(
r*Omega))*(prstar/(
r*Omega));
3188 n[k][1] = ddotr/(
r*Omega*Omega);
3189 n[k][2] = n[k][0]*prstar*prstar;
3190 n[k][3] = prstar/(
r*Omega);
3191 n[k][4] = n[k][3]*cbrt(Omega*Omega);
3192 n[k][5] = n[k][4]*prstar*prstar;
3195 n[k][0] = (prstar/(
r*Omega))*(prstar/(
r*Omega));
3196 n[k][1] = ddotr/(
r*Omega*Omega);
3197 n[k][2] = n[k][0]*prstar*prstar;
3198 n[k][3] = prstar/(
r*Omega);
3200 n[k][4] = n[k][3]*(
r*Omega)*(
r*Omega);
3201 n[k][5] = n[k][4]*prstar*prstar;
3204 n[k][0] = (prstar/(
r*Omega))*(prstar/(
r*Omega));
3205 n[k][1] = ddotr/(
r*Omega*Omega);
3206 n[k][2] = n[k][0]*prstar*prstar;
3207 n[k][3] = prstar/(
r*Omega);
3208 n[k][4] = n[k][3]*cbrt(Omega*Omega);
3209 n[k][5] = n[k][4]*prstar*prstar;
3212 for (
int ki = 0; ki <
KMAX; ki++) {
3213 hlmnqc->
ampli[ki] = 1.;
3214 hlmnqc->
phase[ki] = 0.;
3218 hlmnqc->
ampli[k] = 1. + a1[k]*n[k][0] +
a2[k]*n[k][1] + a3[k]*n[k][2];
3219 hlmnqc->
phase[k] =
b1[k]*n[k][3] +
b2[k]*n[k][4] + b3[k]*n[k][5];
3222 hlmnqc->
ampli[k] = 1. + a1[k]*n[k][0] +
a2[k]*n[k][1] + a3[k]*n[k][2];
3223 hlmnqc->
phase[k] =
b1[k]*n[k][3] +
b2[k]*n[k][4] + b3[k]*n[k][5];
3226 hlmnqc->
ampli[k] = 1. + a1[k]*n[k][0] +
a2[k]*n[k][1] + a3[k]*n[k][2];
3227 hlmnqc->
phase[k] =
b1[k]*n[k][3] +
b2[k]*n[k][4] + b3[k]*n[k][5];
3242 const REAL8 n0 = (prstar/(
r*Omega))*(prstar/(
r*Omega));
3243 const REAL8 n1 = ddotr/(
r*Omega*Omega);
3244 const REAL8 n2 = n0*
SQ(prstar);
3245 const REAL8 n3 = prstar/(
r*Omega);
3246 const REAL8 n4 = n3*cbrt(Omega*Omega);
3247 const REAL8 n5 = n4*
SQ(prstar);
3249 const REAL8 n4_k = n3*
SQ((
r*Omega));
3250 const REAL8 n5_k = n4*
SQ(prstar);
3253 for (
int k = 0; k < maxk; k++) {
3266 nqc->
n[
k0][4] = n4_k;
3267 nqc->
n[
k0][5] = n5_k;
3270 for (
int k = 0; k <
KMAX; k++) {
3271 hlmnqc->
ampli[k] = 1.;
3272 hlmnqc->
phase[k] = 0.;
3275 for (
int k = 0; k < maxk; k++) {
3277 hlmnqc->
ampli[k] += nqc->
a1[k]*nqc->
n[k][0] + nqc->
a2[k]*nqc->
n[k][1] + nqc->
a3[k]*nqc->
n[k][2];
3278 hlmnqc->
phase[k] += nqc->
b1[k]*nqc->
n[k][3] + nqc->
b2[k]*nqc->
n[k][4] + nqc->
b3[k]*nqc->
n[k][5];
3289 REAL8 A_tmp, dA_tmp, omg_tmp, domg_tmp;
3299 const REAL8 nu3 = nu2*nu;
3300 const REAL8 X12 = X1 - X2;
3302 const REAL8 aK3 = aK2*aK;
3303 const REAL8 aK4 = aK2*aK2;
3304 const REAL8 a12 = X1*chi1 - X2*chi2;
3305 const REAL8 aeff = aK + 1./3.*a12*X12;
3306 const REAL8 aeff_omg = aK + a12*X12;
3316 REAL8 c_p1, c_p2, c_p3, c_p4;
3317 REAL8 c_pdA1, c_pdA2, c_pdA3, c_pdA4;
3318 REAL8 c_pdomg1, c_pdomg2;
3320 REAL8 a0_omg_tmp, a1_omg_tmp, a2_omg_tmp, b0_omg_tmp, b1_omg_tmp, b2_omg_tmp, a0_domg_tmp, a1_domg_tmp, a2_domg_tmp, b0_domg_tmp, b1_domg_tmp, b2_domg_tmp, a0_A_tmp, a1_A_tmp , a2_A_tmp, b0_A_tmp, b1_A_tmp, b2_A_tmp, a0_dA_tmp, a1_dA_tmp, a2_dA_tmp, b0_dA_tmp, b1_dA_tmp, b2_dA_tmp, omg_tmp_nu, omg_tmp_equal, domg_tmp_nu, domg_tmp_equal, A_tmp_scale_nu, A_tmp_scale_equal, dA_tmp_scale_nu, dA_tmp_scale_equal ;
3322 REAL8 P[2],
M[4], p1[2], p2[2], p3[2], p4[2], pA[5], pdA[5];
3323 REAL8 pomg[5], pdomg[5], pn0[2], pd1[2], ppdomg1[2], ppdomg2[2], pdA1[2],pdA2[2],pdA3[2],pdA4[2];
3329 for (
int i = 0;
i < size;
i++) {
3334 REAL8 *n1,*n2,*n3,*n4,*n5,*n6, *d_n4,*d_n5, *d2_n4,*d2_n5;
3338 for (
int k=0; k<
KMAX; k++) {
3339 omg[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3340 domg[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3341 m11[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3342 m12[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3343 m13[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3344 m21[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3345 m22[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3346 p1tmp[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3347 p2tmp[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3358 d2_n4 = (
REAL8*) calloc (size,
sizeof(
REAL8));
3359 d2_n5 = (
REAL8*) calloc (size,
sizeof(
REAL8));
3365 for (
int k=0; k<
KMAX; k++) {
3366 D0(this->phase->data->data,
dt, size, omg[k]);
3367 D0(omg[k],
dt, size, domg[k]);
3377 pA[3] = -0.00076165;
3379 A_tmp = pA[0]*aK4 + pA[1]*aK3 + pA[2]*aK2 + pA[3]*aK + pA[4];
3381 pdA[0] = 0.00000927;
3382 pdA[1] = -0.00024550;
3383 pdA[2] = 0.00012469;
3384 pdA[3] = 0.00123845;
3385 pdA[4] = -0.00195014;
3386 dA_tmp = pdA[0]*aK4 + pdA[1]*aK3 + pdA[2]*aK2 + pdA[3]*aK + pdA[4];
3388 pomg[0] = 0.00603482;
3389 pomg[1] = 0.01604555;
3390 pomg[2] = 0.02290799;
3391 pomg[3] = 0.07084587;
3392 pomg[4] = 0.38321834;
3393 omg_tmp = pomg[0]*aK4 + pomg[1]*aK3 + pomg[2]*aK2 + pomg[3]*aK + pomg[4];
3395 pdomg[0] = 0.00024066;
3396 pdomg[1] = 0.00038123;
3397 pdomg[2] = -0.00049714;
3398 pdomg[3] = 0.00041219;
3399 pdomg[4] = 0.01190548;
3400 domg_tmp = pdomg[0]*aK4 + pdomg[1]*aK3 + pdomg[2]*aK2 + pdomg[3]*aK + pdomg[4];
3402 }
else if( nu > 0.16) {
3405 p1[1] = -0.00632114;
3407 p2[1] = -0.01180039;
3408 p3[0] = -0.11617413;
3412 c_p1 = p1[0]*nu + p1[1];
3413 c_p2 = p2[0]*nu + p2[1];
3414 c_p3 = p3[0]*nu + p3[1];
3415 c_p4 = p4[0]*nu + p4[1];
3416 A_tmp = c_p1*aK3 + c_p2*aK2 + c_p3*aK + c_p4;
3418 pdA1[0] = -0.00130824;
3419 pdA1[1] = 0.00006202;
3420 pdA2[0] = 0.00199855;
3421 pdA2[1] = -0.00027474;
3422 pdA3[0] = 0.00218838;
3423 pdA3[1] = 0.00071540;
3424 pdA4[0] = -0.00362779;
3425 pdA4[1] = -0.00105397;
3426 c_pdA1 = pdA1[0]*nu + pdA1[1];
3427 c_pdA2 = pdA2[0]*nu + pdA2[1];
3428 c_pdA3 = pdA3[0]*nu + pdA3[1];
3429 c_pdA4 = pdA4[0]*nu + pdA4[1];
3430 dA_tmp = c_pdA1*aK3 + c_pdA2*aK2 + c_pdA3*aK+ c_pdA4;
3432 pn0[0] = 0.46908067;
3433 pn0[1] = 0.27022141;
3434 pd1[0] = 0.64131115;
3435 pd1[1] = -0.37878384;
3436 n0 = pn0[0]*nu + pn0[1];
3437 d1 = pd1[0]*nu + pd1[1];
3438 omg_tmp = n0/(1 + d1*aK);
3440 ppdomg1[0] = 0.00061175;
3441 ppdomg1[1] = 0.00074001;
3442 ppdomg2[0] = 0.02504442;
3443 ppdomg2[1] = 0.00548217;
3444 c_pdomg1 = ppdomg1[0]*nu + ppdomg1[1];
3445 c_pdomg2 = ppdomg2[0]*nu + ppdomg2[1];
3446 domg_tmp = c_pdomg1*aK + c_pdomg2;
3458 a0_omg_tmp = -0.1460961247;
3459 a1_omg_tmp = 0.0998056;
3460 a2_omg_tmp = -0.118098;
3461 b0_omg_tmp = -0.3430184009;
3462 b1_omg_tmp = 0.0921551;
3463 b2_omg_tmp = -0.0740285;
3464 omg_tmp_nu = +0.5427169903*nu2 +0.2512395608*nu +0.2863992248;
3465 omg_tmp_equal =((a2_omg_tmp*X12*X12 + a1_omg_tmp*X12 + a0_omg_tmp)*aeff_omg+1)/((b2_omg_tmp*X12*X12 +b1_omg_tmp*X12 + b0_omg_tmp)*aeff_omg+1);
3466 omg_tmp = omg_tmp_nu*omg_tmp_equal;
3468 a0_domg_tmp = +0.0604556289;
3469 b0_domg_tmp = -0.0299583285;
3470 a1_domg_tmp = 0.0711715;
3471 a2_domg_tmp = -0.0500886;
3472 b1_domg_tmp = 0.0461239;
3473 b2_domg_tmp = -0.0153068;
3475 domg_tmp_nu = ( +0.0045213831*nu +0.0064934920)/( -1.4466409969*nu+1);
3476 domg_tmp_equal = (a2_domg_tmp*X12*X12 +a1_domg_tmp*X12 +b0_domg_tmp)*aeff_omg*aeff_omg +(b2_domg_tmp*X12*X12 +b1_domg_tmp*X12+a0_domg_tmp)*aeff_omg+1;
3477 domg_tmp = domg_tmp_nu*domg_tmp_equal;
3479 a0_A_tmp = -0.2750516062;
3480 b0_A_tmp = -0.4693776065;
3481 a1_A_tmp = 0.143066;
3482 a2_A_tmp = -0.0425947;
3483 b1_A_tmp = 0.176955;
3484 b2_A_tmp = -0.111902;
3486 A_tmp_scale_nu = -0.9862040409*nu3 +0.8167558040*nu2 -0.0427442282*nu+0.2948879452;
3487 A_tmp_scale_equal = ((a2_A_tmp*X12*X12 + a1_A_tmp*X12 +a0_A_tmp)*aeff+1)/((b2_A_tmp*X12*X12 + b1_A_tmp*X12 +b0_A_tmp)*aeff+1);
3488 A_tmp = A_tmp_scale_nu*A_tmp_scale_equal*(1-0.5*omg_tmp*aeff);
3490 a0_dA_tmp = +0.0037461628;
3491 b0_dA_tmp = +0.0636082543;
3492 a1_dA_tmp = 0.00129393;
3493 a2_dA_tmp = -0.00239069;
3494 b1_dA_tmp = -0.0534209;
3495 b2_dA_tmp = -0.186101;
3497 dA_tmp_scale_nu = ( -0.0847947167*nu -0.0042142765)/( +16.1559461812*nu+1);
3498 dA_tmp_scale_equal = ((a2_dA_tmp*X12*X12 + a1_dA_tmp*X12+ a0_dA_tmp)*aeff)/((b2_dA_tmp*X12*X12 + b1_dA_tmp*X12 + b0_dA_tmp)*aeff+1);
3499 dA_tmp = (dA_tmp_scale_nu +dA_tmp_scale_equal)*omg_tmp;
3504 for (
int k=0; k<
KMAX; k++) {
3512 max_omg[1] = omg_tmp;
3513 max_domg[1] = domg_tmp;
3528 for (
int j=0; j<size; j++) {
3529 pr_star2 =
SQ(pr_star[j]);
3532 n1[j] = pr_star2/(
r2*
w2);
3533 n2[j] = ddotr[j]/(
r[j]*
w2);
3534 n4[j] = pr_star[j]/(
r[j]*
w[j]);
3535 n5[j] = n4[j]*
r2*
w2;
3547 D0(n4,
dt,size, d_n4);
3548 D0(n5,
dt,size, d_n5);
3549 D0(d_n4,
dt,size, d2_n4);
3550 D0(d_n5,
dt,size, d2_n5);
3562 int Omgmax_index = 0;
3563 REAL8 Omg_max = Omg_orb[0];
3564 for (
int j=0; j<size; j++) {
3565 if (Omg_orb[j] > Omg_max) {
3566 Omg_max = Omg_orb[j];
3572 REAL8 tOmgOrb_pk = t[Omgmax_index];
3574 REAL8 tNQC = tOmgOrb_pk - DeltaT_nqc;
3584 for (
int j=0; j<size; j++) {
3597 for (
int k=0; k<
KMAX; k++) {
3598 REAL8 nlm = 1./(sqrt( (this->
l+2)*(this->
l+1)*this->
l*(this->
l-1) ) );
3599 for (
int j=0; j<size; j++) {
3600 p1tmp[k][j] = this->ampl->data->data[j] * nlm;
3606 for (
int k=0; k<
KMAX; k++) {
3607 for (
int j=0; j<size; j++) {
3608 m11[k][j] = n1[j] * p1tmp[k][j];
3609 m12[k][j] = n2[j] * p1tmp[k][j];
3614 for (
int k=0; k<
KMAX; k++) {
3615 D0(m11[k],
dt,size, m21[k]);
3616 D0(m12[k],
dt,size, m22[k]);
3617 D0(p1tmp[k],
dt,size, p2tmp[k]);
3630 for (
int k=0; k<
KMAX; k++) {
3632 ai[k][0] = ai[k][1] = 0.;
3633 bi[k][0] = bi[k][1] = 0.;
3636 P[0] = max_A[k] - p1tmp[k][jmax];
3637 P[1] = max_dA[k] - p2tmp[k][jmax];
3639 M[0] = m11[k][jmax];
3640 M[1] = m12[k][jmax];
3641 M[2] = m21[k][jmax];
3642 M[3] = m22[k][jmax];
3648 oodetM = 1.0/(
M[0]*
M[3]-
M[1]*
M[2]);
3649 if (isfinite((
double) oodetM)) {
3650 ai[k][0] = (
M[3]*P[0] -
M[1]*P[1])*oodetM;
3651 ai[k][1] = (
M[0]*P[1] -
M[2]*P[0])*oodetM;
3655 P[0] = omg[k][jmax] - max_omg[k];
3656 P[1] = domg[k][jmax] - max_domg[k];
3667 oodetM = 1.0/(
M[0]*
M[3]-
M[1]*
M[2]);
3668 if (isfinite((
double) oodetM)) {
3669 bi[k][0] = (
M[3]*P[0] -
M[1]*P[1])*oodetM;
3670 bi[k][1] = (
M[0]*P[1] -
M[2]*P[0])*oodetM;
3686 for (
int k=0; k<
KMAX; k++) {
3687 for (
int j=0; j<size; j++) {
3688 this->ampl->data->data[j] = 1. + ai[k][0]*n1[j] + ai[k][1]*n2[j];
3689 this->phase->data->data[j] = bi[k][0]*n4[j] + bi[k][1]*n5[j];
3698 while (
this && that) {
3699 for (
int j=0; j<size; j++) {
3700 this->ampl->data->data[j] *= that->
ampl->
data->
data[j];
3708 for (
int k=0; k<
KMAX; k++) {
3739 REAL8 A_tmp, dA_tmp, omg_tmp, domg_tmp;
3749 const REAL8 nu3 = nu2*nu;
3750 const REAL8 X12 = X1 - X2;
3752 const REAL8 aK3 = aK2*aK;
3753 const REAL8 aK4 = aK2*aK2;
3754 const REAL8 a12 = X1*chi1 - X2*chi2;
3755 const REAL8 aeff = aK + 1./3.*a12*X12;
3756 const REAL8 aeff_omg = aK + a12*X12;
3766 REAL8 c_p1, c_p2, c_p3, c_p4;
3767 REAL8 c_pdA1, c_pdA2, c_pdA3, c_pdA4;
3768 REAL8 c_pdomg1, c_pdomg2;
3770 REAL8 a0_omg_tmp, a1_omg_tmp, a2_omg_tmp, b0_omg_tmp, b1_omg_tmp, b2_omg_tmp, a0_domg_tmp, a1_domg_tmp, a2_domg_tmp, b0_domg_tmp, b1_domg_tmp, b2_domg_tmp, a0_A_tmp, a1_A_tmp , a2_A_tmp, b0_A_tmp, b1_A_tmp, b2_A_tmp, a0_dA_tmp, a1_dA_tmp, a2_dA_tmp, b0_dA_tmp, b1_dA_tmp, b2_dA_tmp, omg_tmp_nu, omg_tmp_equal, domg_tmp_nu, domg_tmp_equal, A_tmp_scale_nu, A_tmp_scale_equal, dA_tmp_scale_nu, dA_tmp_scale_equal ;
3772 REAL8 P[2],
M[4], p1[2], p2[2], p3[2], p4[2], pA[5], pdA[5];
3773 REAL8 pomg[5], pdomg[5], pn0[2], pd1[2], ppdomg1[2], ppdomg2[2], pdA1[2],pdA2[2],pdA3[2],pdA4[2];
3782 REAL8 *n1,*n2,*n3,*n4,*n5,*n6, *d_n4,*d_n5, *d2_n4,*d2_n5;
3786 for (
int k=0; k<
KMAX; k++)
3788 omg[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3789 domg[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3790 m11[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3791 m12[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3792 m13[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3793 m21[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3794 m22[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3795 p1tmp[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3796 p2tmp[k] = (
REAL8*) calloc (size,
sizeof(
REAL8));
3807 d2_n4 = (
REAL8*) calloc (size,
sizeof(
REAL8));
3808 d2_n5 = (
REAL8*) calloc (size,
sizeof(
REAL8));
3814 for (
int k=0; k<
KMAX; k++)
3816 D0(this->phase->data->data,
dt, size, omg[k]);
3817 D0(omg[k],
dt, size, domg[k]);
3828 pA[3] = -0.00076165;
3830 A_tmp = pA[0]*aK4 + pA[1]*aK3 + pA[2]*aK2 + pA[3]*aK + pA[4];
3832 pdA[0] = 0.00000927;
3833 pdA[1] = -0.00024550;
3834 pdA[2] = 0.00012469;
3835 pdA[3] = 0.00123845;
3836 pdA[4] = -0.00195014;
3837 dA_tmp = pdA[0]*aK4 + pdA[1]*aK3 + pdA[2]*aK2 + pdA[3]*aK + pdA[4];
3839 pomg[0] = 0.00603482;
3840 pomg[1] = 0.01604555;
3841 pomg[2] = 0.02290799;
3842 pomg[3] = 0.07084587;
3843 pomg[4] = 0.38321834;
3844 omg_tmp = pomg[0]*aK4 + pomg[1]*aK3 + pomg[2]*aK2 + pomg[3]*aK + pomg[4];
3846 pdomg[0] = 0.00024066;
3847 pdomg[1] = 0.00038123;
3848 pdomg[2] = -0.00049714;
3849 pdomg[3] = 0.00041219;
3850 pdomg[4] = 0.01190548;
3851 domg_tmp = pdomg[0]*aK4 + pdomg[1]*aK3 + pdomg[2]*aK2 + pdomg[3]*aK + pdomg[4];
3853 }
else if( nu > 0.16) {
3856 p1[1] = -0.00632114;
3858 p2[1] = -0.01180039;
3859 p3[0] = -0.11617413;
3863 c_p1 = p1[0]*nu + p1[1];
3864 c_p2 = p2[0]*nu + p2[1];
3865 c_p3 = p3[0]*nu + p3[1];
3866 c_p4 = p4[0]*nu + p4[1];
3867 A_tmp = c_p1*aK3 + c_p2*aK2 + c_p3*aK + c_p4;
3869 pdA1[0] = -0.00130824;
3870 pdA1[1] = 0.00006202;
3871 pdA2[0] = 0.00199855;
3872 pdA2[1] = -0.00027474;
3873 pdA3[0] = 0.00218838;
3874 pdA3[1] = 0.00071540;
3875 pdA4[0] = -0.00362779;
3876 pdA4[1] = -0.00105397;
3877 c_pdA1 = pdA1[0]*nu + pdA1[1];
3878 c_pdA2 = pdA2[0]*nu + pdA2[1];
3879 c_pdA3 = pdA3[0]*nu + pdA3[1];
3880 c_pdA4 = pdA4[0]*nu + pdA4[1];
3881 dA_tmp = c_pdA1*aK3 + c_pdA2*aK2 + c_pdA3*aK+ c_pdA4;
3883 pn0[0] = 0.46908067;
3884 pn0[1] = 0.27022141;
3885 pd1[0] = 0.64131115;
3886 pd1[1] = -0.37878384;
3887 n0 = pn0[0]*nu + pn0[1];
3888 d1 = pd1[0]*nu + pd1[1];
3889 omg_tmp = n0/(1 + d1*aK);
3891 ppdomg1[0] = 0.00061175;
3892 ppdomg1[1] = 0.00074001;
3893 ppdomg2[0] = 0.02504442;
3894 ppdomg2[1] = 0.00548217;
3895 c_pdomg1 = ppdomg1[0]*nu + ppdomg1[1];
3896 c_pdomg2 = ppdomg2[0]*nu + ppdomg2[1];
3897 domg_tmp = c_pdomg1*aK + c_pdomg2;
3909 a0_omg_tmp = -0.1460961247;
3910 a1_omg_tmp = 0.0998056;
3911 a2_omg_tmp = -0.118098;
3912 b0_omg_tmp = -0.3430184009;
3913 b1_omg_tmp = 0.0921551;
3914 b2_omg_tmp = -0.0740285;
3915 omg_tmp_nu = +0.5427169903*nu2 +0.2512395608*nu +0.2863992248;
3916 omg_tmp_equal =((a2_omg_tmp*X12*X12 + a1_omg_tmp*X12 + a0_omg_tmp)*aeff_omg+1)/((b2_omg_tmp*X12*X12 +b1_omg_tmp*X12 + b0_omg_tmp)*aeff_omg+1);
3917 omg_tmp = omg_tmp_nu*omg_tmp_equal;
3919 a0_domg_tmp = +0.0604556289;
3920 b0_domg_tmp = -0.0299583285;
3921 a1_domg_tmp = 0.0711715;
3922 a2_domg_tmp = -0.0500886;
3923 b1_domg_tmp = 0.0461239;
3924 b2_domg_tmp = -0.0153068;
3926 domg_tmp_nu = ( +0.0045213831*nu +0.0064934920)/( -1.4466409969*nu+1);
3927 domg_tmp_equal = (a2_domg_tmp*X12*X12 +a1_domg_tmp*X12 +b0_domg_tmp)*aeff_omg*aeff_omg +(b2_domg_tmp*X12*X12 +b1_domg_tmp*X12+a0_domg_tmp)*aeff_omg+1;
3928 domg_tmp = domg_tmp_nu*domg_tmp_equal;
3930 a0_A_tmp = -0.2750516062;
3931 b0_A_tmp = -0.4693776065;
3932 a1_A_tmp = 0.143066;
3933 a2_A_tmp = -0.0425947;
3934 b1_A_tmp = 0.176955;
3935 b2_A_tmp = -0.111902;
3937 A_tmp_scale_nu = -0.9862040409*nu3 +0.8167558040*nu2 -0.0427442282*nu+0.2948879452;
3938 A_tmp_scale_equal = ((a2_A_tmp*X12*X12 + a1_A_tmp*X12 +a0_A_tmp)*aeff+1)/((b2_A_tmp*X12*X12 + b1_A_tmp*X12 +b0_A_tmp)*aeff+1);
3939 A_tmp = A_tmp_scale_nu*A_tmp_scale_equal*(1-0.5*omg_tmp*aeff);
3941 a0_dA_tmp = +0.0037461628;
3942 b0_dA_tmp = +0.0636082543;
3943 a1_dA_tmp = 0.00129393;
3944 a2_dA_tmp = -0.00239069;
3945 b1_dA_tmp = -0.0534209;
3946 b2_dA_tmp = -0.186101;
3948 dA_tmp_scale_nu = ( -0.0847947167*nu -0.0042142765)/( +16.1559461812*nu+1);
3949 dA_tmp_scale_equal = ((a2_dA_tmp*X12*X12 + a1_dA_tmp*X12+ a0_dA_tmp)*aeff)/((b2_dA_tmp*X12*X12 + b1_dA_tmp*X12 + b0_dA_tmp)*aeff+1);
3950 dA_tmp = (dA_tmp_scale_nu +dA_tmp_scale_equal)*omg_tmp;
3955 for (
int k=0; k<
KMAX; k++) {
3963 max_omg[1] = omg_tmp;
3964 max_domg[1] = domg_tmp;
3979 for (
int j=0; j<size; j++)
3981 pr_star2 =
SQ(pr_star[j]);
3984 n1[j] = pr_star2/(
r2*
w2);
3985 n2[j] = ddotr[j]/(
r[j]*
w2);
3986 n4[j] = pr_star[j]/(
r[j]*
w[j]);
3987 n5[j] = n4[j]*
r2*
w2;
3999 D0(n4,
dt,size, d_n4);
4000 D0(n5,
dt,size, d_n5);
4001 D0(d_n4,
dt,size, d2_n4);
4002 D0(d_n5,
dt,size, d2_n5);
4013 INT4 Omgmax_index = 0;
4014 REAL8 Omg_max = Omg_orb[0];
4015 for (
int j=0; j<size; j++)
4017 if (Omg_orb[j] > Omg_max)
4019 Omg_max = Omg_orb[j];
4037 REAL8 tOmgOrb_pk = t[Omgmax_index];
4039 REAL8 tNQC = tOmgOrb_pk - DeltaT_nqc;
4049 for (
INT4 j=0; j<size; j++)
4064 for (
int k=0; k<
KMAX; k++)
4066 REAL8 nlm = 1./(sqrt( (this->
l+2)*(this->
l+1)*this->
l*(this->
l-1) ) );
4067 for (
int j=0; j<size; j++)
4069 p1tmp[k][j] = this->ampl->data->data[j] * nlm;
4075 for (
int k=0; k<
KMAX; k++)
4077 for (
int j=0; j<size; j++)
4079 m11[k][j] = n1[j] * p1tmp[k][j];
4080 m12[k][j] = n2[j] * p1tmp[k][j];
4085 for (
int k=0; k<
KMAX; k++)
4087 D0(m11[k],
dt,size, m21[k]);
4088 D0(m12[k],
dt,size, m22[k]);
4089 D0(p1tmp[k],
dt,size, p2tmp[k]);
4102 for (
int k=0; k<
KMAX; k++) {
4104 ai[k][0] = ai[k][1] = 0.;
4105 bi[k][0] = bi[k][1] = 0.;
4108 P[0] = max_A[k] - p1tmp[k][jmax];
4109 P[1] = max_dA[k] - p2tmp[k][jmax];
4111 M[0] = m11[k][jmax];
4112 M[1] = m12[k][jmax];
4113 M[2] = m21[k][jmax];
4114 M[3] = m22[k][jmax];
4120 oodetM = 1.0/(
M[0]*
M[3]-
M[1]*
M[2]);
4121 if (isfinite(oodetM)) {
4122 ai[k][0] = (
M[3]*P[0] -
M[1]*P[1])*oodetM;
4123 ai[k][1] = (
M[0]*P[1] -
M[2]*P[0])*oodetM;
4127 P[0] = omg[k][jmax] - max_omg[k];
4128 P[1] = domg[k][jmax] - max_domg[k];
4139 oodetM = 1.0/(
M[0]*
M[3]-
M[1]*
M[2]);
4140 if (isfinite(oodetM)) {
4141 bi[k][0] = (
M[3]*P[0] -
M[1]*P[1])*oodetM;
4142 bi[k][1] = (
M[0]*P[1] -
M[2]*P[0])*oodetM;
4157 for (
int k=0; k<
KMAX; k++) {
4158 for (
int j=0; j<size; j++) {
4159 this->ampl->data->data[j] = 1. + ai[k][0]*n1[j] + ai[k][1]*n2[j];
4160 this->phase->data->data[j] = bi[k][0]*n4[j] + bi[k][1]*n5[j];
4168 while (
this && that) {
4169 for (
int j=0; j<size; j++) {
4170 this->ampl->data->data[j] *= that->
ampl->
data->
data[j];
4197 n1 = (
REAL8*) calloc (fullsize,
sizeof(
REAL8));
4198 n2 = (
REAL8*) calloc (fullsize,
sizeof(
REAL8));
4199 n4 = (
REAL8*) calloc (fullsize,
sizeof(
REAL8));
4200 n5 = (
REAL8*) calloc (fullsize,
sizeof(
REAL8));
4202 for (
int j=0; j<fullsize; j++) {
4203 pr_star2 =
SQ(pr_star[j]);
4206 n1[j] = pr_star2/(
r2*
w2);
4207 n2[j] = ddotr[j]/(
r[j]*
w2);
4208 n4[j] = pr_star[j]/(
r[j]*
w[j]);
4209 n5[j] = n4[j]*
r2*
w2;
4213 for (
int k=0; k<
KMAX; k++) {
4214 for (
int j=0; j<fullsize; j++) {
4215 this->ampl->data->data[j] *= (1. + ai[k][0]*n1[j] + ai[k][1]*n2[j]);
4216 this->phase->data->data[j] -= (bi[k][0]*n4[j] + bi[k][1]*n5[j]);
4222 for (
int k=0; k<
KMAX; k++)
4262 const int usespeedytail = 1;
4263 const REAL8 X12 = X1-X2;
4283 jhat,Heff,jhat,Heff,
4284 Heff,jhat,Heff,jhat,Heff,
4285 jhat,Heff,jhat,Heff,jhat,Heff,
4286 Heff,jhat,Heff,jhat,Heff,jhat,Heff,
4287 jhat,Heff,jhat,Heff,jhat,Heff,jhat,Heff
4299 REAL8 vphi3 = gsl_pow_int(rw*Omega,3);
4303 REAL8 p4_vphi5 = (2.*nu-1) * gsl_pow_int(rw*Omega,5);
4317 REAL8 vphi3 = gsl_pow_int(rw*Omega,3);
4322 REAL8 p4_vphi5 = (2.*nu-1) * gsl_pow_int(rw*Omega,5);
4335 dyn->
eob_wav_flm_s(
x, nu, X1, X2, chi1, chi2, a1,
a2, C_Q1, C_Q2,
dyn->
clm, usetidal, rholm, flm);
4342 #define RTAIL (1.213061319425267e+00)
4345 if (usespeedytail) {
4356 for (
int k = 0; k <
KMAX; k++) {
4371 for (
INT4 k = 0; k < maxk; k++) {
4385 hlm->
ampli[0] *= X12;
4386 hlm->
ampli[2] *= X12;
4387 hlm->
ampli[4] *= X12;
4388 hlm->
ampli[5] *= X12;
4389 hlm->
ampli[7] *= X12;
4390 hlm->
ampli[9] *= X12;
4391 hlm->
ampli[11] *= X12;
4392 hlm->
ampli[13] *= X12;
4395 for (
int k = 0; k <
KMAX; k++) {
4408 const REAL8 nu2 = nu*nu;
4409 const REAL8 nu3 = nu2*nu;
4410 const REAL8 X12 = sqrt(1.-4.*nu);
4415 const REAL8 n1 = -1.174839;
4416 const REAL8 n2 = 0.354064;
4417 const REAL8 d1 = -0.151961;
4419 const REAL8 c3_eq =
c0*(1. + n1*a12 + n2*a12*a12)/(1.+d1*a12);
4422 const REAL8 cnu = 929.579;
4423 const REAL8 cnu2 = -9178.87;
4424 const REAL8 cnu3 = 23632.3;
4425 const REAL8 ca1_a2 = -104.891;
4427 const REAL8 c3_uneq = cnu*a12*nu*X12 + cnu2*a12*nu2*X12 + cnu3*a12*nu3*X12 + ca1_a2*(a1-
a2)*nu2;
4429 return c3_eq + c3_uneq;
4436 REAL8 DeltaT_nqc = 1.0;
4451 if (((chi1 < -0.9) && (nu < 8./81.)) || ((chi1 < -0.8) && (nu < 11./144.))) {
4476 for (
int k = 0; k <
KMAX; k++) {
4477 for (
int j = 0; j < 6; j++) {
4478 nqc->
flx->
n[k][j] = 0.;
4479 nqc->
hlm->
n[k][j] = 0.;
4481 nqc->
flx->
a1[k] = 0.;
4482 nqc->
flx->
a2[k] = 0.;
4483 nqc->
flx->
a3[k] = 0.;
4484 nqc->
flx->
b1[k] = 0.;
4485 nqc->
flx->
b2[k] = 0.;
4486 nqc->
flx->
b3[k] = 0.;
4488 nqc->
hlm->
a1[k] = 0.;
4489 nqc->
hlm->
a2[k] = 0.;
4490 nqc->
hlm->
a3[k] = 0.;
4491 nqc->
hlm->
b1[k] = 0.;
4492 nqc->
hlm->
b2[k] = 0.;
4493 nqc->
hlm->
b3[k] = 0.;
4528 const REAL8 xnu = 1-4*nu;
4540 nqc->
a1[
k21] = 0.0162387198*(7.32653082*xnu2 + 1.19616248*xnu + 0.73496656);
4541 nqc->
a2[
k21] = -1.80492460*xnu2 + 1.78172686*xnu + 0.30865284;
4544 nqc->
b1[
k21] = -0.0647955017*(3.59934444*xnu2 - 4.08628784*xnu + 1.37890907);
4545 nqc->
b2[
k21] = 1.3410693180*(0.38491989*xnu2 + 0.10969453*xnu + 0.97513971);
4549 nqc->
a1[
k22] = -0.0805236959*( 1 - 2.00332326*xnu2)/( 1 + 3.08595088*xnu2);
4550 nqc->
a2[
k22] = 1.5299534255*( 1 + 1.16438929*xnu2)/( 1 + 1.92033923*xnu2);
4553 nqc->
b1[
k22] = 0.146768094955*( 0.07417121*xnu + 1.01691256);
4554 nqc->
b2[
k22] = 0.896911234248*(-0.61072011*xnu + 0.94295129);
4558 nqc->
a1[
k33] = -0.0377680000*(1 - 14.61548907*xnu2)/( 1 + 2.44559263*xnu2);
4559 nqc->
a2[
k33] = 1.9898000000*(1 + 2.09750346 *xnu2)/( 1 + 2.57489466*xnu2);
4562 nqc->
b1[
k33] = 0.1418400000*(1.07430512 - 1.23906804*xnu + 4.44910652*xnu2);
4563 nqc->
b2[
k33] = 0.6191300000*(0.80672432 + 4.07432829*xnu - 7.47270977*xnu2);
4594 const REAL8 a12 = X1*chi1 - X2*chi2;
4595 const REAL8 X12 = X1 - X2;
4596 const REAL8 aeff = aK + 1./3.*a12*X12;
4597 const REAL8 aeff_omg = aK + a12*X12;
4598 const REAL8 af = abh;
4600 const REAL8 nu3 = nu2*nu;
4602 const REAL8 aeff3 = aeff2*aeff;
4604 const REAL8 af3 = af2*af;
4605 const REAL8 aeff_omg2 =
SQ(aeff_omg);
4606 const REAL8 aeff_omg3 = aeff_omg2*aeff_omg;
4607 const REAL8 aeff_omg4 =
SQ(aeff_omg2);
4617 const INT4 UNUSED k44 = 8;
4619 for (
int k=0; k<
KMAX; k++) {
4621 sigmar[k] = sigmai[k] = 0.;
4622 a1[k] =
a2[k] = a3[k] =
a4[k] = 0.;
4623 b1[k] =
b2[k] = b3[k] = b4[k] = 0.;
4637 alpha21[
k22] = -0.3025985041156393 *nu2 + 0.0032794155172817 *nu + 0.1828276903682022;
4638 alpha1[
k22] = -0.1615300454109702 *nu2 + 0.0147030662812516 *nu + 0.0878204175700328;
4639 c3A[
k22] = 0.8118901739129283 *nu - 0.5584875090785957;
4640 c3phi[
k22] = 0.7156419884962878 *nu + 3.8436474282409803;
4641 c4phi[
k22] = 2.2336960710670901 *nu + 1.4736119175780844;
4642 Domg[
k22] = 0.8846304360111242 *nu2 + 0.0872792137250448 *nu + 0.1058414813686749;
4643 Amrg[
k22] = 1.4935750287318139 *nu2 + 0.2157497669089671 *nu + 1.4292027468283439;
4646 alpha21[
k21] = -0.2741607253846813 *nu2 + 0.0079342900879431 *nu + 0.1835522430667348;
4647 alpha1[
k21] = -0.1277546304610336 *nu2 + 0.0093615534859368 *nu + 0.0882855170502398;
4648 c3A[
k21] = -0.9431151070942140 *nu + 0.2569989171628133;
4649 c3phi[
k21] = -3.4479482376671666 *nu + 2.4755856452648359;
4650 c4phi[
k21] = -3.4024504071619841 *nu + 1.0650118588151427;
4651 Domg[
k21] = 0.2660644668923829 *nu2 + 0.2276854484140649 *nu + 0.0884880283627388;
4652 Amrg[
k21] = -5.7236432632743952 *nu2 + 0.0390010969627653 *nu + 0.4291847351869338;
4656 alpha21[
k33] = -0.3620553934265325 *nu2 + 0.0171973908686402 *nu + 0.1865364041200878;
4657 alpha1[
k33] = -0.1821867653548689 *nu2 + 0.0134440240947561 *nu + 0.0916720214797975;
4658 c3A[
k33] = 2.7565431398030675 *nu - 0.5506682334306747;
4659 c3phi[
k33] = -0.2497526471104979 *nu + 2.3737675006958683;
4660 c4phi[
k33] = -2.9538823110315420 *nu + 1.4483501341373066;
4661 Domg[
k33] = 1.3341439550896721 *nu2 - 0.1717105341058959 *nu + 0.1694617455660599;
4662 Amrg[
k33] = -9.3034388918614841 *nu2 + 1.0189351143222705 *nu + 0.4533252110436300;
4676 sigmar[
k21] = -0.208936*nu3 - 0.028103*nu2 - 0.005383*nu + 0.08896;
4677 sigmai[
k21] = 0.733477*nu3 + 0.188359*nu2 + 0.220659*nu + 0.37367;
4679 sigmar[
k22] = -0.364177*nu3 + 0.010951*nu2 - 0.010591*nu + 0.08896;
4680 sigmai[
k22] = 2.392808*nu3 + 0.051309*nu2 + 0.449425*nu + 0.37365;
4682 sigmar[
k33] = -0.319703*nu3 - 0.030076*nu2-0.009034*nu + 0.09270;
4683 sigmai[
k33] = 2.957425*nu3 + 0.178146*nu2 + 0.709560*nu + 0.59944;
4695 REAL8 omega1_c = -0.0598837831 * af3 + 0.8082136788 * af2 - 1.7408467418 * af + 1;
4696 REAL8 omega1_d = -0.2358960279 * af3 + 1.3152369374 * af2 - 2.0764065380 * af + 1;
4697 omega1[
k22] = 0.3736716844 * (omega1_c/omega1_d);
4700 REAL8 alpha1_c = 0.1211263886 * af3 + 0.7015835813 * af2 - 1.8226060896 * af + 1;
4701 REAL8 alpha1_d = 0.0811633377 * af3 + 0.7201166020 * af2 - 1.8002031358 * af + 1;
4702 alpha1[
k22] = 0.0889623157 * (alpha1_c/alpha1_d);
4705 REAL8 alpha21_c = 0.4764196512 * af3 - 0.0593165805 * af2 - 1.4168096833 * af + 1;
4706 REAL8 alpha21_d = 0.4385578151 * af3 - 0.0763529088 * af2 - 1.3595491146 * af + 1;
4707 alpha21[
k22] = 0.1849525596 * (alpha21_c/alpha21_d);
4710 REAL8 a_c3A = 0.0169543;
4711 REAL8 b_c3A = -0.0799343;
4712 REAL8 c_c3A = -0.115928;
4713 REAL8 c3A_nu = 0.8298678603 * nu - 0.5615838975;
4714 REAL8 c3A_eq = (c_c3A * X12 + 0.0907476903) * aeff3 + (b_c3A * X12 + 0.0227344099) * aeff2 + (a_c3A * X12 - 0.1994944332)*aeff;
4715 c3A[
k22] = c3A_nu + c3A_eq;
4718 REAL8 a_c3phi = -0.462321;
4719 REAL8 b_c3phi = -0.904512;
4720 REAL8 c_c3phi = 0.437747;
4721 REAL8 d_c3phi = 1.8275;
4722 REAL8 c3phi_nu = 0.4558467286 * nu + 3.8883812141;
4723 REAL8 c3phi_equal = (d_c3phi*X12-2.0575868122) * aeff_omg4 +(c_c3phi*X12-0.5051534498)*aeff_omg3 +(b_c3phi*X12+2.5742292762)*aeff_omg2 +(a_c3phi*X12+2.5599640181)*aeff_omg;
4724 c3phi[
k22] = c3phi_nu + c3phi_equal;
4727 REAL8 a_c4phi = -0.449976;
4728 REAL8 b_c4phi = -0.980913;
4729 REAL8 c4phi_nu = 2.0822327682 * nu + 1.4996868401;
4730 REAL8 c4phi_equal = (b_c4phi*X12+3.5695199109) * aeff_omg2 + (a_c4phi * X12 + 4.1312404030) * aeff_omg;
4731 c4phi[
k22] = c4phi_nu + c4phi_equal;
4735 REAL8 a2_omgmx = -0.122735;
4736 REAL8 a1_omgmx = 0.0857478;
4737 REAL8 b2_omgmx = -0.0760023;
4738 REAL8 b1_omgmx = 0.0826514;
4739 REAL8 omgmx_eq_c = (a2_omgmx*X12_2 +a1_omgmx*X12 -0.1416002395) * aeff_omg + 1;
4740 REAL8 omgmx_eq_d = (b2_omgmx*X12_2 +b1_omgmx*X12 -0.3484804901) * aeff_omg + 1;
4741 REAL8 omgmx_eq = omgmx_eq_c/omgmx_eq_d;
4742 REAL8 omgmx = (0.481958619443355 * nu2 + 0.223976694441952 * nu + 0.273813064427363) * omgmx_eq;
4746 REAL8 a2_A_scaled = -0.0820894;
4747 REAL8 a1_A_scaled = 0.176126;
4748 REAL8 b2_A_scaled = -0.150239;
4749 REAL8 b1_A_scaled = 0.20491;
4750 REAL8 A_scaled_eq = ((a2_A_scaled*X12*X12 + a1_A_scaled*X12 -0.2935238329)*aeff + 1)/((b2_A_scaled*X12*X12 + b1_A_scaled*X12 -0.4728707630)*aeff + 1);
4751 REAL8 A_scaled = (+1.826573640739664*nu2 +0.100709438291872*nu +1.438424467327531)*A_scaled_eq;
4753 Amrg[
k22] = A_scaled*(1-0.5*omgmx*aeff);
4754 Domg[
k22] = omega1[
k22] - Mbh*omgmx;
4759 sigmar[
k22] = alpha1[
k22];
4760 sigmai[
k22] = omega1[
k22];
4764 for (
int k=0; k<
KMAX; k++) {
4766 c2A[k] = 0.5*alpha21[k];
4768 a1[k] = Amrg[k] * alpha1[k] * cosh_c3A * cosh_c3A / c2A[k];
4771 a4[k] = Amrg[k] - a1[k] * tanh(c3A[k]);
4775 b1[k] = Domg[k] * (1+c3phi[k]+c4phi[k]) / (
b2[k]*(c3phi[k] + 2.*c4phi[k]));
4798 REAL8 phase = -
b1*log((1. + b3*exp(-
b2*
x) + b4*exp(-2.*
b2*
x))/(1.+b3+b4));
4799 psi[0] = amp * exp(-sigmar*
x);
4800 psi[1] = - (phase - sigmai*
x);
4816 const REAL8 xnu = (1.-4.*nu);
4817 const REAL8 ooMbh = 1./Mbh;
4839 UINT4 index_pk = dynsize-1;
4840 REAL8 Omega_pk = Omega[index_pk];
4841 for (
INT4 j = dynsize-2; j-- ; ) {
4842 if (Omega[j] < Omega_pk)
4845 Omega_pk = Omega[j];
4856 REAL8 *Omega_ptr = &Omega[index_pk-3];
4861 if ( (index_pk + (n-1)/2) > (dynsize-1) ) {
4864 REAL8 Omega_pk_grid[7];
4865 const INT4 ni = (index_pk + (n-1)/2) - (dynsize-1) ;
4868 for (
int j = 0; j < (7-ni); j++)
4869 Omega_pk_grid[j] = Omega_ptr[j];
4872 Omega_pk_grid[6] = 2.*Omega_pk_grid[5]-Omega_pk_grid[4];
4875 Omega_pk_grid[5] = 2.*Omega_pk_grid[4]-Omega_pk_grid[3];
4876 Omega_pk_grid[6] = 2.*Omega_pk_grid[5]-Omega_pk_grid[4];
4878 Omega_pk_grid[4] = 2.*Omega_pk_grid[3]-Omega_pk_grid[2];
4879 Omega_pk_grid[5] = 2.*Omega_pk_grid[4]-Omega_pk_grid[3];
4880 Omega_pk_grid[6] = 2.*Omega_pk_grid[5]-Omega_pk_grid[4];
4884 tOmg_pk =
find_max(n,
dt, tmax, Omega_pk_grid, NULL);
4887 tOmg_pk =
find_max(n,
dt, tmax, Omega_ptr, NULL);
4899 REAL8 *Omega_ptr = &Omega[index_pk-3];
4911 REAL8 tmrgA22 = tOmg_pk-(DeltaT_nqc + 2.)/Mbh;
4913 for (
int k=0; k<
KMAX; k++) {
4921 dtmrg[
k21] = 5.70364338 + 1.85804796*xnu + 4.0332262*xnu*xnu;
4922 dtmrg[
k33] = 4.29550934 - 0.85938*xnu;
4927 for (
int k=0; k<
KMAX; k++) {
4928 tmatch[k] = 2.*ooMbh + tmrg[k];
4942 for (
int k=0; k<
KMAX; k++) {
4943 t_lm[k] = malloc ( size *
sizeof(
REAL8) );
4944 for (
int j = 0; j < size; j++ ) {
4945 t_lm[k][j] = t[j] * ooMbh;
4951 for (
int k = 0; k <
KMAX; k++) {
4952 for (
int j = size-1; j-- ; ) {
4953 if (t_lm[k][j] < tmatch[k]) {
4961 REAL8 t0, tm, psi[2];
4964 for (
int k = 0; k <
KMAX; k++) {
4967 t0 = t_lm[k][idx[k]] - tmrg[k];
4968 eob_wav_ringdown_template(t0, a1[k],
a2[k], a3[k],
a4[k],
b1[k],
b2[k], b3[k], b4[k],
sigma[0][k],
sigma[1][k], psi);
4971 for (
int j = idx[k]; j < size ; j++ ) {
4972 tm = t_lm[k][j] - tmrg[k];
4973 eob_wav_ringdown_template(tm, a1[k],
a2[k], a3[k],
a4[k],
b1[k],
b2[k], b3[k], b4[k],
sigma[0][k],
sigma[1][k], psi);
4977 this_hlm = this_hlm->
next;
4982 for (
int k=0; k<
KMAX; k++)
4984 if(t_lm[k]) free(t_lm[k]);
5001 return A + 0.5 *
u * dA_u;
5009 const double epsabs = 0.;
5010 const double epsrel = 1
e-10;
5011 const gsl_root_fsolver_type *T;
5012 double x, x_lo, x_hi;
5018 switch (tidesFlag) {
5043 gsl_root_fsolver *
s;
5048 T = gsl_root_fsolver_brent;
5049 s = gsl_root_fsolver_alloc (T);
5050 gsl_root_fsolver_set (
s, &F, x_lo, x_hi);
5055 status = gsl_root_fsolver_iterate (
s);
5056 x = gsl_root_fsolver_root (
s);
5057 x_lo = gsl_root_fsolver_x_lower (
s);
5058 x_hi = gsl_root_fsolver_x_upper (
s);
5059 status = gsl_root_test_interval (x_lo, x_hi, epsabs, epsrel);
5061 while (
status == GSL_CONTINUE && iter < max_iter);
5062 gsl_root_fsolver_free (
s);
5065 if (isfinite(
x)) *rLR = (
REAL8)
x;
5070 if (
status == GSL_SUCCESS) {
5073 if (iter >= max_iter) {
5076 if (
status != GSL_SUCCESS) {
5179 const REAL8 z3 = 2.0*nu*(4.0-3.0*nu);
5192 for (
int v=0; v <
TEOB_NV; v++)
5193 buffer[v] = malloc(size *
sizeof (
REAL8));
5195 REAL8 *A_vec = buffer[0];
5196 REAL8 *dA_vec = buffer[1];
5197 REAL8 *B_vec = buffer[2];
5198 REAL8 *sqrtAbyB_vec = buffer[3];
5199 REAL8 *rc_vec = buffer[4];
5200 REAL8 *drc_dr_vec = buffer[5];
5201 REAL8 *uc2_vec = buffer[6];
5202 REAL8 *duc_dr_vec = buffer[7];
5203 REAL8 *dAuc2_dr_vec = buffer[8];
5204 REAL8 *dG_dr_vec = buffer[9];
5205 REAL8 *dG_dprstar_vec = buffer[10];
5206 REAL8 *dG_dprstarbyprstar_vec = buffer[11];
5207 REAL8 *G0_vec = buffer[12];
5208 REAL8 *dG_dr0_vec = buffer[13];
5209 REAL8 *E_vec = buffer[14];
5210 REAL8 *Heff_vec = buffer[15];
5211 REAL8 *Heff_orb_vec = buffer[16];
5212 REAL8 *dpphi_dr_vec = buffer[17];
5213 REAL8 *dprstar_dr_vec = buffer[18];
5214 REAL8 *dphi_dr_vec = buffer[19];
5215 REAL8 *dt_dr_vec = buffer[20];
5218 REAL8 a_coeff, b_coeff, c_coeff, Delta, sol_p, sol_m, j02, uc,
u2, prstar2, dHeff_dpphi, dHeff_dprstar, dHeff_dr, dHeff_dprstarbyprstar, d2Heff_dprstar20,
H,
G, pl_hold,
x, jhat, psi, r_omg, v_phi, Fphi, dr_dtbyprstar, prstar4, Heff_orb_f, Heff_f, E_f;
5224 for (
int i = 0;
i < size;
i++)
5234 C_Oct2, C_Hex1, C_Hex2, usetidal, &rc_vec[
i],
5235 &drc_dr_vec[
i], &pl_hold);
5236 eob_dyn_s_GS(
dyn->
r, rc_vec[
i], drc_dr_vec[
i], aK2, 0.0, 0.0, nu, chi1, chi2, X1, X2, c3, ggm);
5238 G = ggm[2] *S+ggm[3] *Sstar;
5239 dG_dr_vec[
i] = ggm[6] *S+ggm[7] *Sstar;
5240 dG_dprstar_vec[
i] = ggm[4] *S+ggm[5] *Sstar;
5241 dG_dprstarbyprstar_vec[
i] = ggm[10]*S+ggm[11]*Sstar;
5252 dG_dprstar_vec[
i] = 0.0;
5253 dG_dprstarbyprstar_vec[
i] = 0.0;
5260 dG_dr0_vec[
i] = dG_dr_vec[
i];
5263 sqrtAbyB_vec[
i] = sqrt(A_vec[
i]/B_vec[
i]);
5266 duc_dr_vec[
i] = -uc2_vec[
i]*drc_dr_vec[
i];
5267 dAuc2_dr_vec[
i] = uc2_vec[
i]*(dA_vec[
i]-2*A_vec[
i]*uc*drc_dr_vec[
i]);
5274 a_coeff =
SQ(dAuc2_dr_vec[
i]) - 4*A_vec[
i]*uc2_vec[
i]*
SQ(dG_dr_vec[
i]);
5275 b_coeff = 2*dA_vec[
i]*dAuc2_dr_vec[
i] - 4*A_vec[
i]*
SQ(dG_dr_vec[
i]);
5276 c_coeff =
SQ(dA_vec[
i]);
5278 Delta =
SQ(b_coeff) - 4*a_coeff*c_coeff ;
5285 sol_p = (-b_coeff + sqrt(Delta))/(2*a_coeff);
5286 sol_m = (-b_coeff - sqrt(Delta))/(2*a_coeff);
5289 if (dG_dr0_vec[
i] > 0) j02 = sol_p;
5296 a_coeff = dAuc2_dr_vec[
i];
5297 b_coeff = dA_vec[
i];
5299 j02 = -b_coeff/a_coeff;
5306 dprstar_dr_vec[
i] = 0.0;
5312 Sstar,chi1,chi2,X1,X2,aK2,c3,A_vec[
i],dA_vec[
i],
5332 d2Heff_dprstar20 = 1/Heff_orb_vec[
i];
5334 Heff_vec[
i] = Heff_orb_vec[
i];
5340 dyn->
Omg = dHeff_dpphi/E_vec[
i];
5346 dyn->
ddotr = -A_vec[
i]/B_vec[
i]*dHeff_dr*d2Heff_dprstar20;
5370 if (n%2==0) parity = 0;
5374 for (
int i = 0;
i < size;
i++) {
5398 Heff_orb_f = sqrt(A_vec[
i]*(1.0 +
SQ(
dyn->
pphi)*uc2_vec[
i]));
5399 Heff_f = G0_vec[
i]*
dyn->
pphi + Heff_orb_f;
5400 E_f = sqrt(1 + 2*nu*(Heff_f - 1));
5401 psi = (duc_dr_vec[
i] + dG_dr0_vec[
i]*rc_vec[
i]*sqrt(A_vec[
i]/(
SQ(
dyn->
pphi)) + A_vec[
i]*uc2_vec[
i])/A_vec[
i])/(-0.5*dA_vec[
i]);
5402 r_omg = 1.0/cbrt(
SQ(((1./sqrt(rc_vec[
i]*rc_vec[
i]*rc_vec[
i]*psi))+G0_vec[
i])/(E_f)));
5405 jhat =
dyn->
pphi/(r_omg*v_phi);
5411 Heff_orb_f = sqrt(A_vec[
i]*(1.0 +
SQ(
dyn->
pphi)*uc2_vec[
i]));
5412 Heff_f = Heff_orb_f;
5413 psi = 2.*(1.0 + 2.0*nu*(Heff_orb_f - 1.0))/(
SQ(
dyn->
r)*dA_vec[
i]);
5414 r_omg =
dyn->
r*cbrt(psi);
5417 jhat =
dyn->
pphi/(r_omg*v_phi);
5424 dHeff_dprstarbyprstar =
dyn->
pphi*dG_dprstarbyprstar_vec[
i] + (1+2*
z3*A_vec[
i]*uc2_vec[
i]*
SQ(
dyn->
prstar))/Heff_orb_vec[
i];
5425 dr_dtbyprstar = sqrtAbyB_vec[
i]/(E_vec[
i])*dHeff_dprstarbyprstar;
5426 dyn->
prstar = Fphi/dpphi_dr_vec[
i]/dr_dtbyprstar;
5433 eob_dyn_s_GS(
dyn->
r, rc_vec[
i], drc_dr_vec[
i], aK2,
dyn->
prstar, 0.0, nu, chi1, chi2, X1, X2, c3, ggm);
5435 dG_dr_vec[
i] = ggm[6] *S+ggm[7] *Sstar;
5436 dG_dprstar_vec[
i] = ggm[4] *S+ggm[5] *Sstar;
5437 dG_dprstarbyprstar_vec[
i] = ggm[10]*S+ggm[11]*Sstar;
5446 a_coeff = dAuc2_dr_vec[
i];
5447 b_coeff = 2*Heff_orb_vec[
i]*(dG_dr_vec[
i] + dG_dprstar_vec[
i]*dprstar_dr_vec[
i]);
5449 Delta =
SQ(b_coeff) - 4*a_coeff*c_coeff;
5452 sol_m = (-b_coeff - sqrt(Delta))/(2*a_coeff);
5465 eob_ham_s(nu,
dyn->
r, rc_vec[
i], drc_dr_vec[
i],
dyn->
pphi,
dyn->
prstar, S, Sstar, chi1, chi2, X1, X2, aK2, c3, A_vec[
i], dA_vec[
i],
5487 d2Heff_dprstar20 = (1. + 2.*A_vec[
i]*
u2*
z3*prstar2)/Heff_orb_vec[
i];
5489 Heff_vec[
i] = Heff_orb_vec[
i];
5495 dyn->
Omg = dHeff_dpphi/E_vec[
i];
5501 dyn->
ddotr = -A_vec[
i]/B_vec[
i]*dHeff_dr*d2Heff_dprstar20;
5504 dt_dr_vec[
i] = E_vec[
i]/(sqrtAbyB_vec[
i]*dHeff_dprstar);
5505 dphi_dr_vec[
i] =
dyn->
Omg*dt_dr_vec[
i];
5538 for (
int v=0; v <
TEOB_NV; v++)
const INT4 TEOB_MINDEX[KMAX]
const INT4 TEOB_LINDEX[KMAX]
GSL routines for ODE integration https://www.gnu.org/software/gsl/doc/html/ode-initval....
void eob_dyn_s_GS(REAL8 r, REAL8 rc, REAL8 drc_dr, REAL8 UNUSED aK2, REAL8 prstar, REAL8 UNUSED pph, REAL8 nu, REAL8 UNUSED chi1, REAL8 UNUSED chi2, REAL8 UNUSED X1, REAL8 UNUSED X2, REAL8 cN3LO, REAL8 *ggm)
void eob_wav_hlm(LALTEOBResumSDynamics *dyn, LALTEOBResumSWaveformModeSingleTime *hlm)
void eob_wav_hlmNQC(REAL8 UNUSED nu, REAL8 r, REAL8 prstar, REAL8 Omega, REAL8 ddotr, NQCcoefs *nqc, LALTEOBResumSWaveformModeSingleTime *hlmnqc)
void eob_ham(REAL8 nu, REAL8 r, REAL8 pphi, REAL8 prstar, REAL8 A, REAL8 dA, REAL8 *H, REAL8 *Heff, REAL8 *dHeff_dr, REAL8 *dHeff_dprstar, REAL8 *dHeff_dpphi)
void eob_nqc_setcoefs(LALTEOBResumSDynamics *dyn)
int eob_dyn_Npostadiabatic(LALTEOBResumSDynamics *dyn, const REAL8 r0)
static const REAL8 ChlmNewt_phase[35]
REAL8 eob_flx_Flux_s(REAL8 x, REAL8 Omega, REAL8 r_omega, REAL8 E, REAL8 Heff, REAL8 jhat, REAL8 r, REAL8 pr_star, REAL8 ddotr, LALTEOBResumSDynamics *dyn)
void eob_metric(REAL8 r, LALTEOBResumSDynamics *dyn, REAL8 *A, REAL8 *B, REAL8 *dA, REAL8 *d2A, REAL8 *dB)
void eob_flx_Tlm(const REAL8 w, REAL8 *MTlm)
REAL8 eob_dyn_bisecOmegaorb0(LALTEOBResumSDynamics *dyn, REAL8 omg_orb0, REAL8 r0_kepl)
void eob_dyn_ic_s(REAL8 r0, LALTEOBResumSDynamics *dyn, REAL8 y_init[])
void eob_wav_flm(REAL8 x, REAL8 nu, REAL8 clm[KMAX][6], REAL8 *rholm, REAL8 *flm)
int eob_dyn_rhs(REAL8 t, const REAL8 y[], REAL8 dy[], void *d)
void eob_dyn_s_get_rc_NNLO(REAL8 r, REAL8 nu, REAL8 at1, REAL8 at2, REAL8 aK2, REAL8 C_Q1, REAL8 C_Q2, REAL8 UNUSED C_Oct1, REAL8 UNUSED C_Oct2, REAL8 UNUSED C_Hex1, REAL8 UNUSED C_Hex2, int usetidal, REAL8 *rc, REAL8 *drc_dr, REAL8 *d2rc_dr2)
static const double CNlm[35]
void eob_wav_flm_s_SSNLO(REAL8 x, REAL8 nu, REAL8 X1, REAL8 X2, REAL8 UNUSED chi1, REAL8 UNUSED chi2, REAL8 a1, REAL8 a2, REAL8 C_Q1, REAL8 C_Q2, REAL8 clm[KMAX][6], int usetidal, REAL8 *rholm, REAL8 *flm)
void eob_wav_deltalm(REAL8 Hreal, REAL8 Omega, REAL8 nu, REAL8 *dlm)
int eob_dyn_adiabLR(LALTEOBResumSDynamics *dyn, REAL8 *rLR, INT4 tidesFlag)
void eob_ham_s(REAL8 nu, REAL8 r, REAL8 rc, REAL8 drc_dr, REAL8 pphi, REAL8 prstar, REAL8 S, REAL8 Sstar, REAL8 chi1, REAL8 chi2, REAL8 X1, REAL8 X2, REAL8 aK2, REAL8 c3, REAL8 A, REAL8 dA, REAL8 *H, REAL8 *Heff, REAL8 *Heff_orb, REAL8 *dHeff_dr, REAL8 *dHeff_dprstar, REAL8 *dHeff_dpphi, REAL8 *d2Heff_dprstar20)
void eob_wav_ringdown_template(REAL8 x, REAL8 a1, REAL8 a2, REAL8 a3, REAL8 a4, REAL8 b1, REAL8 b2, REAL8 b3, REAL8 b4, REAL8 sigmar, REAL8 sigmai, REAL8 *psi)
void eob_wav_speedyTail(REAL8 Omega, REAL8 Hreal, REAL8 bphys, LALTEOBResumSWaveformModeSingleTime *tlm)
void eob_dyn_s_get_rc_LO(REAL8 r, REAL8 nu, REAL8 at1, REAL8 at2, REAL8 aK2, REAL8 C_Q1, REAL8 C_Q2, REAL8 UNUSED C_Oct1, REAL8 UNUSED C_Oct2, REAL8 UNUSED C_Hex1, REAL8 UNUSED C_Hex2, int usetidal, REAL8 *rc, REAL8 *drc_dr, REAL8 *d2rc_dr2)
REAL8 eob_flx_Flux(REAL8 x, REAL8 Omega, REAL8 r_omega, REAL8 E, REAL8 Heff, REAL8 jhat, REAL8 r, REAL8 pr_star, REAL8 ddotr, LALTEOBResumSDynamics *dyn)
void eob_flx_FlmNewt(REAL8 x, REAL8 nu, REAL8 *Nlm)
void eob_metric_Atidal(REAL8 r, LALTEOBResumSDynamics *dyn, REAL8 *AT, REAL8 *dAT, REAL8 *d2AT)
REAL8 eob_flx_HorizonFlux_s(REAL8 x, REAL8 UNUSED Heff, REAL8 UNUSED jhat, REAL8 UNUSED nu, REAL8 X1, REAL8 X2, REAL8 chi1, REAL8 chi2)
void eob_wav_hlmNewt(REAL8 r, REAL8 Omega, REAL8 phi, REAL8 nu, LALTEOBResumSWaveformModeSingleTime *hlmNewt)
void eob_dyn_s_get_rc_NLO(REAL8 r, REAL8 nu, REAL8 at1, REAL8 at2, REAL8 aK2, REAL8 C_Q1, REAL8 C_Q2, REAL8 UNUSED C_Oct1, REAL8 UNUSED C_Oct2, REAL8 UNUSED C_Hex1, REAL8 UNUSED C_Hex2, int usetidal, REAL8 *rc, REAL8 *drc_dr, REAL8 *d2rc_dr2)
void QNMHybridFitCab(REAL8 nu, REAL8 X1, REAL8 X2, REAL8 chi1, REAL8 chi2, REAL8 aK, REAL8 Mbh, REAL8 abh, REAL8 *a1, REAL8 *a2, REAL8 *a3, REAL8 *a4, REAL8 *b1, REAL8 *b2, REAL8 *b3, REAL8 *b4, REAL8 *sigmar, REAL8 *sigmai, int usespins)
void eob_wav_hlmNQC_nospin201602(REAL8 nu, REAL8 r, REAL8 prstar, REAL8 Omega, REAL8 ddotr, LALTEOBResumSWaveformModeSingleTime *hlmnqc)
void eob_wav_flm_coeffs(REAL8 nu, REAL8 clm[KMAX][6])
void eob_wav_hlmNQC_find_a1a2a3_mrg(LALTEOBResumSDynamics *dyn_mrg, SphHarmPolarTimeSeries *hlm_mrg, SphHarmPolarTimeSeries *hnqc, LALTEOBResumSDynamics *dyn, SphHarmPolarTimeSeries *hlm)
REAL8 eob_dyn_bisecHeff0_s(REAL8 nu, REAL8 chi1, REAL8 chi2, REAL8 X1, REAL8 X2, REAL8 c3, REAL8 pph, REAL8 rorb, REAL8 A, REAL8 dA, REAL8 rc, REAL8 drc_dr, REAL8 ak2, REAL8 S, REAL8 Ss)
void eob_metric_s(REAL8 r, LALTEOBResumSDynamics *dyn, REAL8 *A, REAL8 *B, REAL8 *dA, REAL8 *d2A, REAL8 *dB)
REAL8 eob_dyn_r0_eob(REAL8 f0, LALTEOBResumSDynamics *dyn)
void eob_wav_flm_s_SSLO(REAL8 x, REAL8 nu, REAL8 X1, REAL8 X2, REAL8 UNUSED chi1, REAL8 UNUSED chi2, REAL8 a1, REAL8 a2, REAL8 C_Q1, REAL8 C_Q2, REAL8 clm[KMAX][6], int usetidal, REAL8 *rholm, REAL8 *flm)
void eob_dyn_s_get_rc_NOTIDES(REAL8 r, REAL8 nu, REAL8 at1, REAL8 at2, REAL8 aK2, REAL8 UNUSED C_Q1, REAL8 UNUSED C_Q2, REAL8 UNUSED C_Oct1, REAL8 UNUSED C_Oct2, REAL8 UNUSED C_Hex1, REAL8 UNUSED C_Hex2, int usetidal, REAL8 *rc, REAL8 *drc_dr, REAL8 *d2rc_dr2)
void eob_dyn_s_get_rc_NOSPIN(REAL8 r, REAL8 UNUSED nu, REAL8 UNUSED at1, REAL8 UNUSED at2, REAL8 UNUSED aK2, REAL8 UNUSED C_Q1, REAL8 UNUSED C_Q2, REAL8 UNUSED C_Oct1, REAL8 UNUSED C_Oct2, REAL8 UNUSED C_Hex1, REAL8 UNUSED C_Hex2, int UNUSED usetidal, REAL8 UNUSED *rc, REAL8 UNUSED *drc_dr, REAL8 UNUSED *d2rc_dr2)
void eob_nqc_setcoefs_nospin201602(REAL8 nu, NQCcoefs *nqc)
void eob_wav_ringdown(LALTEOBResumSDynamics *dyn, SphHarmPolarTimeSeries *hlm)
REAL8 eob_nqc_timeshift(REAL8 nu, REAL8 chi1)
void eob_wav_hhatlmTail(REAL8 Omega, REAL8 Hreal, REAL8 bphys, LALTEOBResumSWaveformModeSingleTime *tlm)
REAL8 eob_flx_HorizonFlux(REAL8 x, REAL8 Heff, REAL8 jhat, REAL8 nu)
void eob_metric_A5PNlog(REAL8 r, REAL8 nu, REAL8 *A, REAL8 *dA, REAL8 *d2A)
REAL8 eob_c3_fit_global(REAL8 nu, REAL8 UNUSED chi1, REAL8 UNUSED chi2, REAL8 UNUSED X1, REAL8 UNUSED X2, REAL8 a1, REAL8 a2)
void eob_wav_hlmNQC_find_a1a2a3(LALTEOBResumSDynamics *dyn, SphHarmPolarTimeSeries *h, SphHarmPolarTimeSeries *hnqc)
static const REAL8 ChlmNewt_ampli[35]
void eob_dyn_s_get_rc_NNLO_S4(REAL8 r, REAL8 nu, REAL8 at1, REAL8 at2, REAL8 aK2, REAL8 C_Q1, REAL8 C_Q2, REAL8 C_Oct1, REAL8 C_Oct2, REAL8 C_Hex1, REAL8 C_Hex2, int usetidal, REAL8 *rc, REAL8 *drc_dr, REAL8 *d2rc_dr2)
int eob_dyn_rhs_s(REAL8 t, const REAL8 y[], REAL8 dy[], void *d)
REAL8 eob_dyn_Omegaorb0(REAL8 r, void *params)
REAL8 eob_dyn_r0_Kepler(REAL8 f0)
REAL8 eob_dyn_fLR(REAL8 r, void *params)
void eob_wav_hlmTidal(REAL8 x, LALTEOBResumSDynamics *dyn, REAL8 *hTidallm)
REAL8 eob_dyn_DHeff0(REAL8 x, void *params)
void eob_dyn_ic(REAL8 r0, LALTEOBResumSDynamics *dyn, REAL8 y_init[])
REAL8 cumint3(REAL8 *f, REAL8 *x, const INT4 n, REAL8 *sum)
INT4 D0(REAL8 *f, REAL8 dx, INT4 n, REAL8 *df)
XLALWignerdMatrix, XLALWignerDMatrix in <lal/SphericalHarmonics.h>
REAL8 Eulerlog(const REAL8 x, const INT4 m)
INT4 D2(REAL8 *f, REAL8 dx, INT4 n, REAL8 *d2f)
REAL8 find_max(const INT4 n, REAL8 dx, REAL8 x0, REAL8 *f, REAL8 *fmax)
static double sigma(const double a, const double b, const sysq *system)
Internal function that computes the spin-spin couplings.
static vector d(const double L_norm, const double J_norm, const vector roots)
Internal function that returns the coefficients "d_0", "d_2" and "d_4" from 1703.03967 corresponding ...
static REAL8 UNUSED C2(REAL8 e0, REAL8 f0)
static REAL8 UNUSED C1(REAL8 Mtotal)
static REAL8 UNUSED C4(REAL8 e0, REAL8 f0)
static REAL8 UNUSED C5(REAL8 e0, REAL8 f0)
static REAL8 UNUSED C3(REAL8 e0, REAL8 f0)
static REAL8 UNUSED C0(REAL8 eta)
static REAL8 UNUSED z3(REAL8 e0, REAL8 f0, REAL8 inc, REAL8 bet, REAL8 FPlus, REAL8 FCross)
static REAL8 UNUSED C6(REAL8 eta)
This file is part of TEOBResumS.
#define DEQUAL(a, b, eps)
@ LAL_SIM_INSPIRAL_GETIDES_GSF23
@ LAL_SIM_INSPIRAL_GETIDES_NNLO
@ LAL_SIM_INSPIRAL_GETIDES_GSF2
@ LAL_SIM_INSPIRAL_GMTIDES_GSF
@ LAL_SIM_INSPIRAL_GMTIDES_PN
#define XLAL_ERROR_VOID(...)
#define XLAL_ERROR_REAL8(...)
#define XLAL_CHECK_VOID(assertion,...)
EOBDynSGetRCFunc eob_dyn_s_get_rc
EOBWavFlmSFunc eob_wav_flm_s
INT4 use_tidal_gravitomagnetic
REAL8 * data[TEOB_DYNAMICS_NVARS]
Multipolar coefficients for NQC waveform.
NQC data for flux and waveform.
LALTEOBResumSDynamics * dyn
REAL8TimeSeries * ampl
The sequences of mode amplitude.
REAL8TimeSeries * phase
The sequences of mode phase (not modulo 2Pi).
struct tagSphHarmPolarTimeSeries * next
next pointer
REAL8Sequence * tdata
Timestamp values.