Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALSimulation 6.2.0.1-6c6b863
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
LALSimInspiralTaylorT4.c
Go to the documentation of this file.
1/*
2 * Copyright (C) 2008 J. Creighton, S. Fairhurst, B. Krishnan, L. Santamaria, D. Keppel, Evan Ochsner, Les Wade
3 *
4 * This program is free software; you can redistribute it and/or modify
5 * it under the terms of the GNU General Public License as published by
6 * the Free Software Foundation; either version 2 of the License, or
7 * (at your option) any later version.
8 *
9 * This program is distributed in the hope that it will be useful,
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 * GNU General Public License for more details.
13 *
14 * You should have received a copy of the GNU General Public License
15 * along with with program; see the file COPYING. If not, write to the
16 * Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
17 * MA 02110-1301 USA
18 */
19
20#include <math.h>
21
22#include <gsl/gsl_const.h>
23#include <gsl/gsl_errno.h>
24#include <gsl/gsl_math.h>
25#include <gsl/gsl_odeiv.h>
26
27#include <lal/LALSimInspiral.h>
28#include <lal/LALConstants.h>
29#include <lal/LALStdlib.h>
30#include <lal/TimeSeries.h>
31#include <lal/Units.h>
32
35#include "check_series_macros.h"
36
37#ifdef __GNUC__
38#define UNUSED __attribute__ ((unused))
39#else
40#define UNUSED
41#endif
42
43/*
44 * This structure contains the intrinsic parameters and post-newtonian
45 * co-efficients for the energy and angular acceleration expansions.
46 * These are computed by XLALSimInspiralTaylorT4Setup routine.
47 */
48
49typedef struct
50tagexpnCoeffsTaylorT4 {
51
52 /*Angular velocity coefficient*/
54
55 /* Taylor expansion coefficents in domega/dt*/
56 REAL8 aatN,aat2,aat3,aat4,aat5,aat6,aat6l,aat7,aat10,aat12;
57
58 /* symmetric mass ratio, total mass, component masses*/
59 REAL8 nu,m,m1,m2,mu,chi1,chi2;
61
63 REAL8 v, /**< post-Newtonian parameter */
64 expnCoeffsdEnergyFlux *ak /**< Taylor expansion coefficients */
65);
66
68 REAL8 v, /**< post-Newtonian parameter */
69 expnCoeffsTaylorT4 *ak /**< Taylor expansion coefficients */
70);
71
72/*
73 * This strucuture contains pointers to the functions for calculating
74 * the post-newtonian terms at the desired order. They can be set by
75 * XLALSimInspiralTaylorT4Setup by passing an appropriate PN order.
76 */
77
78typedef struct
79tagexpnFuncTaylorT4
80{
84
85
86/*
87 * Computes the rate of increase of the orbital frequency for a post-Newtonian
88 * inspiral.
89 *
90 * Implements Equation 3.6 of: Alessandra Buonanno, Bala R Iyer, Evan
91 * Ochsner, Yi Pan, and B S Sathyaprakash, "Comparison of post-Newtonian
92 * templates for compact binary inspiral signals in gravitational-wave
93 * detectors", Phys. Rev. D 80, 084043 (2009), arXiv:0907.0700v1
94 */
95
96static REAL8
98 REAL8 v, /**< post-Newtonian parameter */
99 expnCoeffsTaylorT4 *ak /**< PN co-efficients and intrinsic parameters */
100 )
101{
102 REAL8 ans;
103 REAL8 v9;
104
105 v9 = pow(v, 9.0);
106
107 ans = ak->aatN * v9;
108
109 return ans;
110}
111
112static REAL8
114 REAL8 v, /**< post-Newtonian parameter */
115 expnCoeffsTaylorT4 *ak /**< PN co-efficients and intrinsic parameters */
116 )
117{
118 REAL8 ans;
119 REAL8 v2,v9,v10,v12;
120
121 v2 = v*v;
122 v9 = pow(v, 9.0);
123 v10 = v9*v;
124 v12 = v10*v2;
125
126 ans = ak->aatN * (1.
127 + ak->aat2*v2
128 + ak->aat10*v10
129 + ak->aat12*v12);
130
131 ans *= v9;
132
133 return ans;
134}
135
136static REAL8
138 REAL8 v, /**< post-Newtonian parameter */
139 expnCoeffsTaylorT4 *ak /**< PN co-efficients and intrinsic parameters */
140 )
141{
142 REAL8 ans;
143 REAL8 v2,v3,v9,v10,v12;
144
145 v2 = v*v;
146 v3 = v2*v;
147 v9 = v3*v3*v3;
148 v10 = v9*v;
149 v12 = v10*v2;
150
151 ans = ak->aatN * (1.
152 + ak->aat2*v2
153 + ak->aat3*v3
154 + ak->aat10*v10
155 + ak->aat12*v12);
156
157 ans *= v9;
158
159 return ans;
160}
161
162static REAL8
164 REAL8 v, /**< post-Newtonian parameter */
165 expnCoeffsTaylorT4 *ak /**< PN co-efficients and intrinsic parameters */
166 )
167{
168 REAL8 ans;
169 REAL8 v2,v3,v4,v9,v10,v12;
170
171 v2 = v*v;
172 v3 = v2*v;
173 v4 = v3*v;
174 v9 = v4*v4*v;
175 v10 = v9*v;
176 v12 = v10*v2;
177
178 ans = ak->aatN * (1.
179 + ak->aat2*v2
180 + ak->aat3*v3
181 + ak->aat4*v4
182 + ak->aat10*v10
183 + ak->aat12*v12);
184
185 ans *= v9;
186
187 return ans;
188}
189
190static REAL8
192 REAL8 v, /**< post-Newtonian parameter */
193 expnCoeffsTaylorT4 *ak /**< PN co-efficients and intrinsic parameters */
194 )
195{
196 REAL8 ans;
197 REAL8 v2,v3,v4,v5,v9,v10,v12;
198
199 v2 = v*v;
200 v3 = v2*v;
201 v4 = v3*v;
202 v5 = v4*v;
203 v9 = v5*v4;
204 v10 = v9*v;
205 v12 = v10*v2;
206
207 ans = ak->aatN * (1.
208 + ak->aat2*v2
209 + ak->aat3*v3
210 + ak->aat4*v4
211 + ak->aat5*v5
212 + ak->aat10*v10
213 + ak->aat12*v12);
214
215 ans *= v9;
216
217 return ans;
218}
219
220static REAL8
222 REAL8 v, /**< post-Newtonian parameter */
223 expnCoeffsTaylorT4 *ak /**< PN co-efficients and intrinsic parameters */
224 )
225{
226 REAL8 ans;
227 REAL8 v2,v3,v4,v5,v6,v9,v10,v12;
228
229 v2 = v*v;
230 v3 = v2*v;
231 v4 = v3*v;
232 v5 = v4*v;
233 v6 = v5*v;
234 v9 = v6*v3;
235 v10 = v9*v;
236 v12 = v10*v2;
237
238 ans = ak->aatN * (1.
239 + ak->aat2*v2
240 + ak->aat3*v3
241 + ak->aat4*v4
242 + ak->aat5*v5
243 + (ak->aat6 + ak->aat6l*log(v))*v6
244 + ak->aat10*v10
245 + ak->aat12*v12);
246
247 ans *= v9;
248
249 return ans;
250}
251
252static REAL8
254 REAL8 v, /**< post-Newtonian parameter */
255 expnCoeffsTaylorT4 *ak /**< PN co-efficients and intrinsic parameters */
256 )
257{
258 REAL8 ans;
259 REAL8 v2,v3,v4,v5,v6,v7,v9,v10,v12;
260
261 v2 = v*v;
262 v3 = v2*v;
263 v4 = v3*v;
264 v5 = v4*v;
265 v6 = v5*v;
266 v7 = v6*v;
267 v9 = v7*v2;
268 v10 = v9*v;
269 v12 = v10*v2;
270
271 ans = ak->aatN * (1.
272 + ak->aat2*v2
273 + ak->aat3*v3
274 + ak->aat4*v4
275 + ak->aat5*v5
276 + (ak->aat6 + ak->aat6l*log(v))*v6
277 + ak->aat7*v7
278 + ak->aat10*v10
279 + ak->aat12*v12);
280
281 ans *= v9;
282
283 return ans;
284}
285
286typedef struct
287{
291
292/**
293 * This function is used in the call to the GSL integrator.
294 */
295static int
296XLALSimInspiralTaylorT4PNEvolveOrbitIntegrand(double UNUSED t, const double y[], double ydot[], void *params)
297{
299 ydot[0] = p->func(y[0],&p->ak);
300 ydot[1] = y[0]*y[0]*y[0]*p->ak.av;
301 t = 0.0;
302 return GSL_SUCCESS;
303}
304
305
306/*
307 * Set up the expnCoeffsTaylorT4 and expnFuncTaylorT4 structures for
308 * generating a TaylorT4 waveform and select the post-newtonian
309 * functions corresponding to the desired order.
310 *
311 * Inputs given in SI units.
312 */
313static int
315 expnCoeffsTaylorT4 *ak, /**< coefficients for TaylorT4 evolution [modified] */
316 expnFuncTaylorT4 *f, /**< functions for TaylorT4 evolution [modified] */
317 expnCoeffsdEnergyFlux *akdEF, /**< coefficients for Energy calculation [modified] */
318 REAL8 m1, /**< mass of companion 1 */
319 REAL8 m2, /**< mass of companion 2 */
320 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
321 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
322 LALSimInspiralTidalOrder tideO, /**< twice PN order of tidal effects */
323 int O /**< twice post-Newtonian order */
324)
325{
326 ak->m1 = m1;
327 ak->m2 = m2;
328 ak->m = ak->m1 + ak->m2;
329 ak->mu = m1 * m2 / ak->m;
330 ak->nu = ak->mu/ak->m;
331 ak->chi1 = ak->m1/ak->m;
332 ak->chi2 = ak->m2/ak->m;
333
334 /* Angular velocity co-efficient */
335 ak->av = pow(LAL_C_SI, 3.0)/(LAL_G_SI*ak->m);
336
337 /* PN co-efficients for energy */
342
343 /* PN co-efficients for angular acceleration */
352
353 /* Tidal coefficients for energy and angular acceleration */
354 akdEF->ETa5 = 0.;
355 akdEF->ETa6 = 0.;
356 ak->aat10 = 0.;
357 ak->aat12 = 0.;
358 switch( tideO )
359 {
364#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
365 __attribute__ ((fallthrough));
366#endif
370#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
371 __attribute__ ((fallthrough));
372#endif
374 break;
375 default:
376 XLALPrintError("XLAL Error - %s: Invalid tidal PN order %d\nSee LALSimInspiralTidalOrder enum in LALSimInspiralWaveformFlags.h for valid tidal orders.\n",
377 __func__, tideO );
379 break;
380 }
381
382 switch (O)
383 {
384 case 0:
387 break;
388 case 1:
389 XLALPrintError("XLAL Error - %s: PN approximant not supported for PN order %d\n", __func__,O);
391 break;
392 case 2:
395 break;
396 case 3:
399 break;
400 case 4:
403 break;
404 case 5:
407 break;
408 case 6:
411 break;
412 case 7:
413 case -1:
416 break;
417 case 8:
418 XLALPrintError("XLAL Error - %s: PN approximant not supported for PN order %d\n", __func__,O);
420 break;
421 default:
422 XLALPrintError("XLAL Error - %s: Unknown PN order in switch\n", __func__);
424 }
425
426 return 0;
427}
428
429/**
430 * @addtogroup LALSimInspiralTaylorXX_c
431 * @{
432 * @name Routines for TaylorT4 Waveforms
433 * @sa
434 * Section IIIB of Alessandra Buonanno, Bala R Iyer, Evan
435 * Ochsner, Yi Pan, and B S Sathyaprakash, "Comparison of post-Newtonian
436 * templates for compact binary inspiral signals in gravitational-wave
437 * detectors", Phys. Rev. D 80, 084043 (2009), arXiv:0907.0700v1
438 * @{
439 */
440
441/**
442 * Evolves a post-Newtonian orbit using the Taylor T4 method.
443 *
444 * See:
445 * Michael Boyle, Duncan A. Brown, Lawrence E. Kidder, Abdul H. Mroue,
446 * Harald P. Pfeiffer, Mark A. Scheel, Gregory B. Cook, and Saul A. Teukolsky
447 * "High-accuracy comparison of numerical relativity simulations with
448 * post-Newtonian expansions"
449 * <a href="http://arxiv.org/abs/0710.0158v2">arXiv:0710.0158v2</a>.
450 */
452 REAL8TimeSeries **v, /**< post-Newtonian parameter [returned] */
453 REAL8TimeSeries **phi, /**< orbital phase [returned] */
454 REAL8 phiRef, /**< reference orbital phase (rad) */
455 REAL8 deltaT, /**< sampling interval (s) */
456 REAL8 m1, /**< mass of companion 1 (kg) */
457 REAL8 m2, /**< mass of companion 2 (kg) */
458 REAL8 f_min, /**< start frequency (Hz) */
459 REAL8 fRef, /**< reference frequency (Hz) */
460 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
461 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
462 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
463 int O /**< twice post-Newtonian order */
464 )
465{
466 const UINT4 blocklen = 1024;
467 const REAL8 visco = 1./sqrt(6.);
468 REAL8 VRef = 0.;
470 expnFuncTaylorT4 expnfunc;
473
474 if(XLALSimInspiralTaylorT4Setup(&ak,&expnfunc,&akdEF,m1,m2,lambda1,lambda2,
475 tideO,O))
477
478 params.func=expnfunc.angacc4;
479 params.ak=ak;
480
481 REAL8 E;
482 UINT4 j, len, idxRef = 0;
484 double y[2];
485 double yerr[2];
486 const gsl_odeiv_step_type *T = gsl_odeiv_step_rk4;
487 gsl_odeiv_step *s;
488 gsl_odeiv_system sys;
489
490 /* setup ode system */
492 sys.jacobian = NULL;
493 sys.dimension = 2;
494 sys.params = &params;
495
496 /* allocate memory */
497 *v = XLALCreateREAL8TimeSeries( "ORBITAL_VELOCITY_PARAMETER", &tc, 0., deltaT, &lalDimensionlessUnit, blocklen );
498 *phi = XLALCreateREAL8TimeSeries( "ORBITAL_PHASE", &tc, 0., deltaT, &lalDimensionlessUnit, blocklen );
499 if ( !v || !phi )
501
502 y[0] = (*v)->data->data[0] = cbrt(LAL_PI*LAL_G_SI*ak.m*f_min)/LAL_C_SI;
503 y[1] = (*phi)->data->data[0] = 0.;
504 E = expnfunc.energy4(y[0],&akdEF);
505 if (XLALIsREAL8FailNaN(E))
507 j = 0;
508
509 s = gsl_odeiv_step_alloc(T, 2);
510 while (1) {
511 ++j;
512 gsl_odeiv_step_apply(s, j*deltaT, deltaT, y, yerr, NULL, NULL, &sys);
513 /* ISCO termination condition for quadrupole, 1pN, 2.5pN */
514 if ( y[0] > visco ) {
515 XLALPrintInfo("XLAL Info - %s: PN inspiral terminated at ISCO\n", __func__);
516 break;
517 }
518 if ( j >= (*v)->data->length ) {
519 if ( ! XLALResizeREAL8TimeSeries(*v, 0, (*v)->data->length + blocklen) )
521 if ( ! XLALResizeREAL8TimeSeries(*phi, 0, (*phi)->data->length + blocklen) )
523 }
524 (*v)->data->data[j] = y[0];
525 (*phi)->data->data[j] = y[1];
526 }
527 gsl_odeiv_step_free(s);
528
529 /* make the correct length */
530 if ( ! XLALResizeREAL8TimeSeries(*v, 0, j) )
532 if ( ! XLALResizeREAL8TimeSeries(*phi, 0, j) )
534
535 /* adjust to correct time */
536 XLALGPSAdd(&(*v)->epoch, -1.0*j*deltaT);
537 XLALGPSAdd(&(*phi)->epoch, -1.0*j*deltaT);
538
539 /* Do a constant phase shift to get desired value of phiRef */
540 len = (*phi)->data->length;
541 /* For fRef==0, phiRef is phase of last sample */
542 if( fRef == 0. )
543 phiRef -= (*phi)->data->data[len-1];
544 /* For fRef==fmin, phiRef is phase of first sample */
545 else if( fRef == f_min )
546 phiRef -= (*phi)->data->data[0];
547 /* phiRef is phase when f==fRef */
548 else
549 {
550 VRef = pow(LAL_PI * LAL_G_SI*(m1+m2) * fRef, 1./3.) / LAL_C_SI;
551 j = 0;
552 do {
553 idxRef = j;
554 j++;
555 } while ((*v)->data->data[j] <= VRef);
556 phiRef -= (*phi)->data->data[idxRef];
557 }
558 for (j = 0; j < len; ++j)
559 (*phi)->data->data[j] += phiRef;
560
561 return (int)(*v)->data->length;
562}
563
564
565/**
566 * Driver routine to compute the post-Newtonian inspiral waveform.
567 *
568 * This routine allows the user to specify different pN orders
569 * for phasing calcuation vs. amplitude calculations.
570 */
572 REAL8TimeSeries **hplus, /**< +-polarization waveform */
573 REAL8TimeSeries **hcross, /**< x-polarization waveform */
574 REAL8 phiRef, /**< reference orbital phase (rad) */
575 UNUSED REAL8 v0, /**< tail-term gauge choice (default = 1) */
576 REAL8 deltaT, /**< sampling interval (s) */
577 REAL8 m1, /**< mass of companion 1 (kg) */
578 REAL8 m2, /**< mass of companion 2 (kg) */
579 REAL8 f_min, /**< start frequency (Hz) */
580 REAL8 fRef, /**< reference frequency (Hz) */
581 REAL8 r, /**< distance of source (m) */
582 REAL8 i, /**< inclination of source (rad) */
583 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
584 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
585 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
586 int amplitudeO, /**< twice post-Newtonian amplitude order */
587 int phaseO /**< twice post-Newtonian phase order */
588 )
589{
590
591 /* The Schwarzschild ISCO frequency - for sanity checking fRef */
592 REAL8 fISCO = pow(LAL_C_SI,3) / (pow(6.,3./2.)*LAL_PI*(m1+m2)*LAL_G_SI);
593
594 /* Sanity check fRef value */
595 if( fRef < 0. )
596 {
597 XLALPrintError("XLAL Error - %s: fRef = %f must be >= 0\n",
598 __func__, fRef);
600 }
601 if( fRef != 0. && fRef < f_min )
602 {
603 XLALPrintError("XLAL Error - %s: fRef = %f must be > fStart = %f\n",
604 __func__, fRef, f_min);
606 }
607 if( fRef >= fISCO )
608 {
609 XLALPrintError("XLAL Error - %s: fRef = %f must be < Schwar. ISCO=%f\n",
610 __func__, fRef, fISCO);
612 }
613
614
616 REAL8TimeSeries *phi;
617 int status;
618 int n;
619 n = XLALSimInspiralTaylorT4PNEvolveOrbit(&v, &phi, phiRef, deltaT,
620 m1, m2, f_min, fRef, lambda1, lambda2, tideO, phaseO);
621 if ( n < 0 )
623 status = XLALSimInspiralPNPolarizationWaveforms(hplus, hcross, v, phi,
624 v0, m1, m2, r, i, amplitudeO);
627 if ( status < 0 )
629 return n;
630}
631
632/**
633 * Driver routine to compute the -2 spin-weighted spherical harmonic modes
634 * using TaylorT4 phasing.
635 */
637 UNUSED REAL8 v0, /**< tail-term gauge choice (default = 1) */
638 REAL8 deltaT, /**< sampling interval (s) */
639 REAL8 m1, /**< mass of companion 1 (kg) */
640 REAL8 m2, /**< mass of companion 2 (kg) */
641 REAL8 f_min, /**< starting GW frequency (Hz) */
642 REAL8 fRef, /**< reference GW frequency (Hz) */
643 REAL8 r, /**< distance of source (m) */
644 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
645 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
646 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
647 int amplitudeO, /**< twice post-Newtonian amplitude order */
648 int phaseO, /**< twice post-Newtonian phase order */
649 int lmax /**< generate all modes with l <= lmax */
650 )
651{
652 SphHarmTimeSeries *hlm=NULL;
653 /* The Schwarzschild ISCO frequency - for sanity checking fRef */
654 REAL8 fISCO = pow(LAL_C_SI,3) / (pow(6.,3./2.)*LAL_PI*(m1+m2)*LAL_G_SI);
655
656 /* Sanity check fRef value */
657 if( fRef < 0. )
658 {
659 XLALPrintError("XLAL Error - %s: fRef = %f must be >= 0\n",
660 __func__, fRef);
662 }
663 if( fRef != 0. && fRef < f_min )
664 {
665 XLALPrintError("XLAL Error - %s: fRef = %f must be > fStart = %f\n",
666 __func__, fRef, f_min);
668 }
669 if( fRef >= fISCO )
670 {
671 XLALPrintError("XLAL Error - %s: fRef = %f must be < Schwar. ISCO=%f\n",
672 __func__, fRef, fISCO);
674 }
675
677 REAL8TimeSeries *phi;
678 int n;
680 m1, m2, f_min, fRef, lambda1, lambda2, tideO, phaseO);
681 if ( n < 0 )
683 int m, l;
685 for(l=2; l<=lmax; l++){
686 for(m=-l; m<=l; m++){
688 m1, m2, r, amplitudeO, l, m);
689 if ( !hxx ){
691 }
692 hlm = XLALSphHarmTimeSeriesAddMode(hlm, hxx, l, m);
694 }
695 }
698 return hlm;
699}
700
701/**
702 * Driver routine to compute the -2 spin-weighted spherical harmonic mode
703 * using TaylorT4 phasing.
704 */
706 UNUSED REAL8 v0, /**< tail-term gauge choice (default = 1) */
707 REAL8 deltaT, /**< sampling interval (s) */
708 REAL8 m1, /**< mass of companion 1 (kg) */
709 REAL8 m2, /**< mass of companion 2 (kg) */
710 REAL8 f_min, /**< starting GW frequency (Hz) */
711 REAL8 fRef, /**< reference GW frequency (Hz) */
712 REAL8 r, /**< distance of source (m) */
713 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
714 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
715 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
716 int amplitudeO, /**< twice post-Newtonian amplitude order */
717 int phaseO, /**< twice post-Newtonian phase order */
718 int l, /**< l index of mode */
719 int m /**< m index of mode */
720 )
721{
723 /* The Schwarzschild ISCO frequency - for sanity checking fRef */
724 REAL8 fISCO = pow(LAL_C_SI,3) / (pow(6.,3./2.)*LAL_PI*(m1+m2)*LAL_G_SI);
725
726 /* Sanity check fRef value */
727 if( fRef < 0. )
728 {
729 XLALPrintError("XLAL Error - %s: fRef = %f must be >= 0\n",
730 __func__, fRef);
732 }
733 if( fRef != 0. && fRef < f_min )
734 {
735 XLALPrintError("XLAL Error - %s: fRef = %f must be > fStart = %f\n",
736 __func__, fRef, f_min);
738 }
739 if( fRef >= fISCO )
740 {
741 XLALPrintError("XLAL Error - %s: fRef = %f must be < Schwar. ISCO=%f\n",
742 __func__, fRef, fISCO);
744 }
745
747 REAL8TimeSeries *phi;
748 int n;
750 m1, m2, f_min, fRef, lambda1, lambda2, tideO, phaseO);
751 if ( n < 0 )
754 m1, m2, r, amplitudeO, l, m);
757 return hlm;
758}
759
760/**
761 * Driver routine to compute the post-Newtonian inspiral waveform.
762 *
763 * This routine uses the same pN order for phasing and amplitude
764 * (unless the order is -1 in which case the highest available
765 * order is used for both of these -- which might not be the same).
766 *
767 * Constant log term in amplitude set to 1. This is a gauge choice.
768 */
770 REAL8TimeSeries **hplus, /**< +-polarization waveform */
771 REAL8TimeSeries **hcross, /**< x-polarization waveform */
772 REAL8 phiRef, /**< reference orbital phase (rad) */
773 REAL8 deltaT, /**< sampling interval (Hz) */
774 REAL8 m1, /**< mass of companion 1 (kg) */
775 REAL8 m2, /**< mass of companion 2 (kg) */
776 REAL8 f_min, /**< start frequency (Hz) */
777 REAL8 fRef, /**< reference frequency (Hz) */
778 REAL8 r, /**< distance of source (m) */
779 REAL8 i, /**< inclination of source (rad) */
780 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
781 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
782 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
783 int O /**< twice post-Newtonian order */
784 )
785{
786 /* set v0 to default value 1 */
787 return XLALSimInspiralTaylorT4PNGenerator(hplus, hcross, phiRef, 1.0,
788 deltaT, m1, m2, f_min, fRef, r, i, lambda1, lambda2, tideO, O, O);
789}
790
791
792/**
793 * Driver routine to compute the restricted post-Newtonian inspiral waveform.
794 *
795 * This routine computes the phasing to the specified order, but
796 * only computes the amplitudes to the Newtonian (quadrupole) order.
797 *
798 * Constant log term in amplitude set to 1. This is a gauge choice.
799 */
801 REAL8TimeSeries **hplus, /**< +-polarization waveform */
802 REAL8TimeSeries **hcross, /**< x-polarization waveform */
803 REAL8 phiRef, /**< reference orbital phase (rad) */
804 REAL8 deltaT, /**< sampling interval (s) */
805 REAL8 m1, /**< mass of companion 1 (kg) */
806 REAL8 m2, /**< mass of companion 2 (kg) */
807 REAL8 f_min, /**< start frequency (Hz) */
808 REAL8 fRef, /**< reference frequency (Hz) */
809 REAL8 r, /**< distance of source (m) */
810 REAL8 i, /**< inclination of source (rad) */
811 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
812 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
813 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
814 int O /**< twice post-Newtonian phase order */
815 )
816{
817 /* use Newtonian order for amplitude */
818 /* set v0 to default value 1 */
819 return XLALSimInspiralTaylorT4PNGenerator(hplus, hcross, phiRef, 1.0,
820 deltaT, m1, m2, f_min, fRef, r, i, lambda1, lambda2, tideO, 0, O);
821}
822
823/** @} */
824/** @} */
825
826#if 0
827#include <lal/PrintFTSeries.h>
828#include <lal/PrintFTSeries.h>
829int main(void)
830{
831 LIGOTimeGPS tc = { 888888888, 222222222 };
832 REAL8 phic = 1.0;
833 REAL8 deltaT = 1.0/16384.0;
834 REAL8 m1 = 1.4*LAL_MSUN_SI;
835 REAL8 m2 = 1.4*LAL_MSUN_SI;
836 REAL8 r = 1e6*LAL_PC_SI;
837 REAL8 i = 0.5*LAL_PI;
838 REAL8 f_min = 100.0;
839 REAL8 fRef = 0.;
840 int O = -1;
841 REAL8TimeSeries *hplus;
842 REAL8TimeSeries *hcross;
843 XLALSimInspiralTaylorT4PN(&hplus, &hcross, &tc, phic, deltaT, m1, m2, f_min, fRef, r, i, lambda1, lambda2, tideO, O);
844 LALDPrintTimeSeries(hplus, "hp.dat");
845 LALDPrintTimeSeries(hcross, "hc.dat");
849 return 0;
850}
851#endif
void LALCheckMemoryLeaks(void)
static REAL8 UNUSED XLALSimInspiralTaylorT4wdot_12PNTidalCoeff(REAL8 mByM)
static REAL8 UNUSED XLALSimInspiralPNEnergy_4PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT4wdot_7PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT4wdot_4PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralPNEnergy_10PNTidalCoeff(REAL8 mByM)
static REAL8 UNUSED XLALSimInspiralPNEnergy_2PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT4wdot_0PNCoeff(REAL8 eta)
Computes the PN Coefficients for using in the TaylorT4 frequency equation.
static REAL8 UNUSED XLALSimInspiralTaylorT4wdot_2PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT4wdot_6PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralPNEnergy_0PNCoeff(REAL8 eta)
Computes the PN Coefficients for using in the PN energy equation.
static REAL8 UNUSED XLALSimInspiralTaylorT4wdot_10PNTidalCoeff(REAL8 mByM)
static REAL8 UNUSED XLALSimInspiralTaylorT4wdot_6PNLogCoeff(REAL8 UNUSED eta)
static REAL8 UNUSED XLALSimInspiralPNEnergy_12PNTidalCoeff(REAL8 mByM)
static REAL8 UNUSED XLALSimInspiralTaylorT4wdot_5PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralPNEnergy_6PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT4wdot_3PNCoeff(REAL8 UNUSED eta)
static int XLALSimInspiralTaylorT4Setup(expnCoeffsTaylorT4 *ak, expnFuncTaylorT4 *f, expnCoeffsdEnergyFlux *akdEF, REAL8 m1, REAL8 m2, REAL8 lambda1, REAL8 lambda2, LALSimInspiralTidalOrder tideO, int O)
static REAL8 XLALSimInspiralAngularAcceleration4_2PN(REAL8 v, expnCoeffsTaylorT4 *ak)
static REAL8 XLALSimInspiralAngularAcceleration4_3PN(REAL8 v, expnCoeffsTaylorT4 *ak)
static REAL8 XLALSimInspiralAngularAcceleration4_7PN(REAL8 v, expnCoeffsTaylorT4 *ak)
REAL8() SimInspiralEnergy4(REAL8 v, expnCoeffsdEnergyFlux *ak)
static REAL8 XLALSimInspiralAngularAcceleration4_4PN(REAL8 v, expnCoeffsTaylorT4 *ak)
static REAL8 XLALSimInspiralAngularAcceleration4_0PN(REAL8 v, expnCoeffsTaylorT4 *ak)
REAL8() SimInspiralAngularAcceleration4(REAL8 v, expnCoeffsTaylorT4 *ak)
static int XLALSimInspiralTaylorT4PNEvolveOrbitIntegrand(double UNUSED t, const double y[], double ydot[], void *params)
This function is used in the call to the GSL integrator.
static REAL8 XLALSimInspiralAngularAcceleration4_5PN(REAL8 v, expnCoeffsTaylorT4 *ak)
static REAL8 XLALSimInspiralAngularAcceleration4_6PN(REAL8 v, expnCoeffsTaylorT4 *ak)
Module containing the energy and flux functions for waveform generation.
static REAL8 UNUSED XLALSimInspiralEt6(REAL8 v, expnCoeffsdEnergyFlux *ak)
static REAL8 UNUSED XLALSimInspiralEt0(REAL8 v, expnCoeffsdEnergyFlux *ak)
static REAL8 UNUSED XLALSimInspiralEt4(REAL8 v, expnCoeffsdEnergyFlux *ak)
static REAL8 UNUSED XLALSimInspiralEt2(REAL8 v, expnCoeffsdEnergyFlux *ak)
int main(int argc, char *argv[])
Definition: bh_qnmode.c:226
int s
Definition: bh_qnmode.c:137
int l
Definition: bh_qnmode.c:135
double i
Definition: bh_ringdown.c:118
#define __attribute__(x)
#define LAL_C_SI
#define LAL_MSUN_SI
#define LAL_PI
#define LAL_MTSUN_SI
#define LAL_PC_SI
#define LAL_G_SI
#define LIGOTIMEGPSZERO
double REAL8
uint32_t UINT4
int XLALSimInspiralPNPolarizationWaveforms(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, REAL8TimeSeries *V, REAL8TimeSeries *Phi, REAL8 v0, REAL8 m1, REAL8 m2, REAL8 r, REAL8 i, int ampO)
Given time series for a binary's orbital dynamical variables, construct the waveform polarizations h+...
LALSimInspiralTidalOrder
Enumeration of allowed PN orders of tidal effects.
@ LAL_SIM_INSPIRAL_TIDAL_ORDER_5PN
@ LAL_SIM_INSPIRAL_TIDAL_ORDER_6PN
@ LAL_SIM_INSPIRAL_TIDAL_ORDER_ALL
@ LAL_SIM_INSPIRAL_TIDAL_ORDER_0PN
COMPLEX16TimeSeries * XLALCreateSimInspiralPNModeCOMPLEX16TimeSeriesLALConvention(REAL8TimeSeries *v, REAL8TimeSeries *phi, REAL8 m1, REAL8 m2, REAL8 r, int O, int l, int m)
Computes h(l,m) mode timeseries of spherical harmonic decomposition of the post-Newtonian inspiral wa...
SphHarmTimeSeries * XLALSimInspiralTaylorT4PNModes(UNUSED REAL8 v0, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, REAL8 fRef, REAL8 r, REAL8 lambda1, REAL8 lambda2, LALSimInspiralTidalOrder tideO, int amplitudeO, int phaseO, int lmax)
Driver routine to compute the -2 spin-weighted spherical harmonic modes using TaylorT4 phasing.
int XLALSimInspiralTaylorT4PNEvolveOrbit(REAL8TimeSeries **v, REAL8TimeSeries **phi, REAL8 phiRef, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, REAL8 fRef, REAL8 lambda1, REAL8 lambda2, LALSimInspiralTidalOrder tideO, int O)
Evolves a post-Newtonian orbit using the Taylor T4 method.
int XLALSimInspiralTaylorT4PNGenerator(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, REAL8 phiRef, UNUSED REAL8 v0, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, REAL8 fRef, REAL8 r, REAL8 i, REAL8 lambda1, REAL8 lambda2, LALSimInspiralTidalOrder tideO, int amplitudeO, int phaseO)
Driver routine to compute the post-Newtonian inspiral waveform.
int XLALSimInspiralTaylorT4PN(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, REAL8 phiRef, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, REAL8 fRef, REAL8 r, REAL8 i, REAL8 lambda1, REAL8 lambda2, LALSimInspiralTidalOrder tideO, int O)
Driver routine to compute the post-Newtonian inspiral waveform.
COMPLEX16TimeSeries * XLALSimInspiralTaylorT4PNMode(UNUSED REAL8 v0, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, REAL8 fRef, REAL8 r, REAL8 lambda1, REAL8 lambda2, LALSimInspiralTidalOrder tideO, int amplitudeO, int phaseO, int l, int m)
Driver routine to compute the -2 spin-weighted spherical harmonic mode using TaylorT4 phasing.
int XLALSimInspiralTaylorT4PNRestricted(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, REAL8 phiRef, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, REAL8 fRef, REAL8 r, REAL8 i, REAL8 lambda1, REAL8 lambda2, LALSimInspiralTidalOrder tideO, int O)
Driver routine to compute the restricted post-Newtonian inspiral waveform.
SphHarmTimeSeries * XLALSphHarmTimeSeriesAddMode(SphHarmTimeSeries *appended, const COMPLEX16TimeSeries *inmode, UINT4 l, INT4 m)
Prepend a node to a linked list of SphHarmTimeSeries, or create a new head.
void LALDPrintTimeSeries(REAL8TimeSeries *series, const CHAR *filename)
static const INT4 r
static const INT4 m
void XLALDestroyCOMPLEX16TimeSeries(COMPLEX16TimeSeries *series)
void XLALDestroyREAL8TimeSeries(REAL8TimeSeries *series)
REAL8TimeSeries * XLALResizeREAL8TimeSeries(REAL8TimeSeries *series, int first, size_t length)
REAL8TimeSeries * XLALCreateREAL8TimeSeries(const CHAR *name, const LIGOTimeGPS *epoch, REAL8 f0, REAL8 deltaT, const LALUnit *sampleUnits, size_t length)
const LALUnit lalDimensionlessUnit
#define XLAL_ERROR_NULL(...)
int int int XLALPrintInfo(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
#define XLAL_ERROR(...)
int XLALPrintError(const char *fmt,...) _LAL_GCC_PRINTF_FORMAT_(1
static _LAL_INLINE_ int XLALIsREAL8FailNaN(REAL8 val)
XLAL_EFUNC
XLAL_EINVAL
LIGOTimeGPS * XLALGPSAdd(LIGOTimeGPS *epoch, REAL8 dt)
list y
list mu
string status
p
REAL8Sequence * data
Structure to carry a collection of spherical harmonic modes in COMPLEX16 time series.
SimInspiralEnergy4 * energy4
SimInspiralAngularAcceleration4 * angacc4
Definition: burst.c:245
double V
Definition: unicorn.c:25
double f_min
Definition: unicorn.c:22
double deltaT
Definition: unicorn.c:24