Loading [MathJax]/extensions/TeX/AMSsymbols.js
LALSimulation 6.2.0.1-5e288d3
All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
LALSimInspiralTaylorT3.c
Go to the documentation of this file.
1/*
2 * Copyright (C) 2011 Drew Keppel, J. Creighton, S. Fairhurst, B. Krishnan, L. Santamaria, Stas Babak, David Churches, B.S. Sathyaprakash, Craig Robinson , Thomas Cokelaer, 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 <lal/LALSimInspiral.h>
23#include <lal/FindRoot.h>
24#include <lal/LALConstants.h>
25#include <lal/LALStdlib.h>
26#include <lal/TimeSeries.h>
27#include <lal/Units.h>
28
31
32#include "check_series_macros.h"
33
34#ifdef __GNUC__
35#define UNUSED __attribute__ ((unused))
36#else
37#define UNUSED
38#endif
39
40typedef struct
41tagexpnCoeffsTaylorT3 {
42 int ieta;
43 /* Taylor expansion coefficents in phi(t)*/
44 REAL8 ptaN, pta2, pta3, pta4, pta5, pta6, pta7, ptl6, pta10, pta12;
45 /* Taylor expansion coefficents in f(t)*/
46 REAL8 ftaN, fta2, fta3, fta4, fta5, fta6, fta7, ftl6, fta10, fta12;
47
48 /* sampling interval*/
50 /* symmetric mass ratio, total mass, component masses*/
51 REAL8 eta, totalmass, m1, m2, chi1, chi2;
52
53 /* initial and final values of frequency, time, velocity; lso
54 values of velocity and frequency; final phase.*/
55 REAL8 f0, fn, t0, theta_lso, v0, vn, vf, vlso, flso, phiC;
56
57 /* last stable orbit and pole defined by various Taylor and P-approximants*/
58 REAL8 vlsoT0, vlsoT2, vlsoT4, vlsoT6;
60
62 REAL8 td,
64
66 REAL8 td,
68
69typedef struct
70tagexpnFuncTaylorT3
71{
75
76typedef struct
77{
78 REAL8 (*func)(REAL8 tC, expnCoeffsTaylorT3 *ak);
80}
82
84{
86 REAL8 freq, f;
87
88 in = (FreqInFromChirptime *) pars;
89 freq = in->func(tC, &(in->ak));
90 if (XLAL_IS_REAL8_FAIL_NAN(freq))
92 f = freq - in->ak.f0;
93
94 /*
95 fprintf(stderr, "Here freq=%e f=%e tc=%e f0=%e\n", freq, *f, tC, in->ak.f0);
96 */
97
98 return f;
99}
100
101static REAL8
103 REAL8 td,
105 )
106{
107 REAL8 theta,theta3;
108 REAL8 frequency;
109
110 if (ak == NULL)
112
113 theta = pow(td,-0.125);
114 theta3 = theta*theta*theta;
115
116 frequency = theta3*ak->ftaN;
117
118 return frequency;
119}
120
121
122static REAL8
124 REAL8 td,
126 )
127{
128 REAL8 theta,theta2,theta3, theta10, theta12;
129 REAL8 frequency;
130
131 if (ak == NULL)
133
134 theta = pow(td,-0.125);
135 theta2 = theta*theta;
136 theta3 = theta2*theta;
137 theta10 = theta2*theta2*theta2*theta2*theta2;
138 theta12 = theta10*theta2;
139
140 frequency = theta3*ak->ftaN * (1.
141 + ak->fta2*theta2
142 + ak->fta10*theta10
143 + ak->fta12*theta12);
144
145 return frequency;
146}
147
148
149static REAL8
151 REAL8 td,
153 )
154{
155 REAL8 theta,theta2,theta3, theta10, theta12;
156 REAL8 frequency;
157
158 if (ak == NULL)
160
161 theta = pow(td,-0.125);
162 theta2 = theta*theta;
163 theta3 = theta2*theta;
164 theta10 = theta3*theta3*theta3*theta;
165 theta12 = theta10*theta2;
166
167 frequency = theta3*ak->ftaN * (1.
168 + ak->fta2*theta2
169 + ak->fta3*theta3
170 + ak->fta10*theta10
171 + ak->fta12*theta12);
172
173 return frequency;
174}
175
176
177static REAL8
179 REAL8 td,
181 )
182{
183 REAL8 theta,theta2,theta3,theta4, theta10, theta12;
184 REAL8 frequency;
185
186 if (ak == NULL)
188
189 theta = pow(td,-0.125);
190 theta2 = theta*theta;
191 theta3 = theta2*theta;
192 theta4 = theta3*theta;
193 theta10 = theta4*theta4*theta2;
194 theta12 = theta10*theta2;
195
196 frequency = theta3*ak->ftaN * (1.
197 + ak->fta2*theta2
198 + ak->fta3*theta3
199 + ak->fta4*theta4
200 + ak->fta10*theta10
201 + ak->fta12*theta12);
202
203 return frequency;
204}
205
206
207static REAL8
209 REAL8 td,
211 )
212{
213 REAL8 theta,theta2,theta3,theta4,theta5, theta10, theta12;
214 REAL8 frequency;
215
216 if (ak == NULL)
218
219 theta = pow(td,-0.125);
220 theta2 = theta*theta;
221 theta3 = theta2*theta;
222 theta4 = theta3*theta;
223 theta5 = theta4*theta;
224 theta10 = theta5*theta5;
225 theta12 = theta10*theta2;
226
227 frequency = theta3*ak->ftaN * (1.
228 + ak->fta2*theta2
229 + ak->fta3*theta3
230 + ak->fta4*theta4
231 + ak->fta5*theta5
232 + ak->fta10*theta10
233 + ak->fta12*theta12);
234
235 return frequency;
236}
237
238
239static REAL8
241 REAL8 td,
243 )
244{
245 REAL8 theta,theta2,theta3,theta4,theta5,theta6, theta10, theta12;
246 REAL8 frequency;
247
248 if (ak == NULL)
250
251 theta = pow(td,-0.125);
252 theta2 = theta*theta;
253 theta3 = theta2*theta;
254 theta4 = theta3*theta;
255 theta5 = theta4*theta;
256 theta6 = theta5*theta;
257 theta10 = theta6*theta4;
258 theta12 = theta10*theta2;
259
260 frequency = theta3*ak->ftaN * (1.
261 + ak->fta2*theta2
262 + ak->fta3*theta3
263 + ak->fta4*theta4
264 + ak->fta5*theta5
265 + (ak->fta6 + ak->ftl6*log(2.*theta))*theta6
266 + ak->fta10*theta10
267 + ak->fta12*theta12);
268
269 return frequency;
270}
271
272
273static REAL8
275 REAL8 td,
277 )
278{
279 REAL8 theta,theta2,theta3,theta4,theta5,theta6,theta7, theta10, theta12;
280 REAL8 frequency;
281
282 if (ak == NULL)
284
285 theta = pow(td,-0.125);
286 theta2 = theta*theta;
287 theta3 = theta2*theta;
288 theta4 = theta3*theta;
289 theta5 = theta4*theta;
290 theta6 = theta5*theta;
291 theta7 = theta6*theta;
292 theta10 = theta7*theta3;
293 theta12 = theta10*theta2;
294
295 frequency = theta3*ak->ftaN * (1.
296 + ak->fta2*theta2
297 + ak->fta3*theta3
298 + ak->fta4*theta4
299 + ak->fta5*theta5
300 + (ak->fta6 + ak->ftl6*log(2.*theta))*theta6
301 + ak->fta7*theta7
302 + ak->fta10*theta10
303 + ak->fta12*theta12);
304
305 return frequency;
306}
307
308
309static REAL8
311 REAL8 td,
313 )
314{
315 REAL8 theta5;
316 REAL8 phase;
317
318 if (ak == NULL)
320
321 theta5 = pow(td,-0.625);
322 phase = (ak->ptaN/theta5);
323
324 return phase;
325}
326
327
328static REAL8
330 REAL8 td,
332 )
333{
334 REAL8 theta,theta2,theta5, theta10, theta12;
335 REAL8 phase;
336
337 if (ak == NULL)
339
340 theta = pow(td,-0.125);
341 theta2 = theta*theta;
342 theta5 = theta2*theta2*theta;
343 theta10 = theta5*theta5;
344 theta12 = theta10*theta2;
345
346 phase = (ak->ptaN/theta5) * (1.
347 + ak->pta2*theta2
348 + ak->pta10*theta10
349 + ak->pta12*theta12);
350
351 return phase;
352}
353
354
355static REAL8
357 REAL8 td,
359 )
360{
361 REAL8 theta,theta2,theta3,theta5, theta10, theta12;
362 REAL8 phase;
363
364 if (ak == NULL)
366
367 theta = pow(td,-0.125);
368 theta2 = theta*theta;
369 theta3 = theta2*theta;
370 theta5 = theta2*theta3;
371 theta10 = theta5*theta5;
372 theta12 = theta10*theta2;
373
374 phase = (ak->ptaN/theta5) * (1.
375 + ak->pta2*theta2
376 + ak->pta3*theta3
377 + ak->pta10*theta10
378 + ak->pta12*theta12);
379
380 return phase;
381}
382
383
384static REAL8
386 REAL8 td,
388 )
389{
390 REAL8 theta,theta2,theta3,theta4,theta5, theta10, theta12;
391 REAL8 phase;
392
393 if (ak == NULL)
395
396 theta = pow(td,-0.125);
397 theta2 = theta*theta;
398 theta3 = theta2*theta;
399 theta4 = theta3*theta;
400 theta5 = theta4*theta;
401 theta10 = theta5*theta5;
402 theta12 = theta10*theta2;
403
404 phase = (ak->ptaN/theta5) * (1.
405 + ak->pta2*theta2
406 + ak->pta3*theta3
407 + ak->pta4*theta4
408 + ak->pta10*theta10
409 + ak->pta12*theta12);
410
411 return phase;
412}
413
414
415static REAL8
417 REAL8 td,
419 )
420{
421 REAL8 theta,theta2,theta3,theta4,theta5, theta10, theta12;
422 REAL8 phase;
423
424 if (ak == NULL)
426
427 theta = pow(td,-0.125);
428 theta2 = theta*theta;
429 theta3 = theta2*theta;
430 theta4 = theta3*theta;
431 theta5 = theta4*theta;
432 theta10 = theta5*theta5;
433 theta12 = theta10*theta2;
434
435 phase = (ak->ptaN/theta5) * (1.
436 + ak->pta2*theta2
437 + ak->pta3*theta3
438 + ak->pta4*theta4
439 + ak->pta5 * log(theta/ak->theta_lso) * theta5
440 + ak->pta10*theta10
441 + ak->pta12*theta12);
442
443 return phase;
444}
445
446
447static REAL8
449 REAL8 td,
451 )
452{
453 REAL8 theta,theta2,theta3,theta4,theta5,theta6, theta10, theta12;
454 REAL8 phase;
455
456 if (ak == NULL)
458
459 theta = pow(td,-0.125);
460 theta2 = theta*theta;
461 theta3 = theta2*theta;
462 theta4 = theta3*theta;
463 theta5 = theta4*theta;
464 theta6 = theta5*theta;
465 theta10 = theta6*theta4;
466 theta12 = theta10*theta2;
467
468 phase = (ak->ptaN/theta5) * (1.
469 + ak->pta2*theta2
470 + ak->pta3*theta3
471 + ak->pta4*theta4
472 + ak->pta5*log(theta/ak->theta_lso)*theta5
473 +(ak->ptl6*log(2.*theta) + ak->pta6)*theta6
474 + ak->pta10*theta10
475 + ak->pta12*theta12);
476
477 return phase;
478}
479
480
481static REAL8
483 REAL8 td,
485 )
486{
487 REAL8 theta,theta2,theta3,theta4,theta5,theta6,theta7, theta10, theta12;
488 REAL8 phase;
489
490 if (ak == NULL)
492
493 theta = pow(td,-0.125);
494 theta2 = theta*theta;
495 theta3 = theta2*theta;
496 theta4 = theta3*theta;
497 theta5 = theta4*theta;
498 theta6 = theta5*theta;
499 theta7 = theta6*theta;
500 theta10 = theta7*theta3;
501 theta12 = theta10*theta2;
502
503 phase = (ak->ptaN/theta5) * (1.
504 + ak->pta2*theta2
505 + ak->pta3*theta3
506 + ak->pta4*theta4
507 + ak->pta5*log(theta/ak->theta_lso)*theta5
508 +(ak->ptl6*log(2.*theta) + ak->pta6)*theta6
509 + ak->pta7*theta7
510 + ak->pta10*theta10
511 + ak->pta12*theta12);
512
513 return phase;
514}
515
516
517/*
518 * Returns the sum of chirp times to a specified order.
519 *
520 * Computes the sum of the chirp times to a specified order. Inputs given in SI
521 * units.
522 */
524 REAL8 m1, /**< mass of companion 1 */
525 REAL8 m2, /**< mass of companion 2 */
526 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
527 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
528 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
529 REAL8 f_min, /**< start frequency */
530 int O /**< twice post-Newtonian order */
531 )
532{
533 REAL8 tN, t2, t3 ,t4, t5, t6, t6l, t7, tC, t10, t12;
534 REAL8 v, v2, v3, v4, v5, v6, v7, v8, v10, v12;
535 REAL8 piFl = LAL_PI * f_min;
536 REAL8 m = m1 + m2;
537 REAL8 chi1 = m1/(m1+m2), chi2 = m2/(m1+m2);
538 REAL8 eta = m1 * m2 / m / m;
539 m *= LAL_G_SI / pow(LAL_C_SI, 3.0); /* convert m from kilograms to seconds */
540
541 /* Should use the coefficients from LALInspiraSetup.c to avoid errors.
542 */
543 v = cbrt(piFl * m);
544 v2 = v*v;
545 v3 = v*v2;
546 v4 = v*v3;
547 v5 = v*v4;
548 v6 = v*v5;
549 v7 = v*v6;
550 v8 = v*v7;
551 v10 = v2*v8;
552 v12 = v2*v10;
553
562
563 /* Tidal co-efficients for t(v). */
564 t10 = 0.;
565 t12 = 0.;
566 switch( tideO )
567 {
572#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
573 __attribute__ ((fallthrough));
574#endif
578#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
579 __attribute__ ((fallthrough));
580#endif
582 break;
583 default:
584 XLALPrintError("XLAL Error - %s: Invalid tidal PN order %d\nSee LALSimInspiralTidalOrder enum in LALSimInspiralWaveformFlags.h for valid tidal orders.\n",
585 __func__, tideO );
587 break;
588 }
589
590 switch (O) {
591 case 0:
592 case 1:
593 t2 = 0.;
594#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
595 __attribute__ ((fallthrough));
596#endif
597 case 2:
598 t3 = 0.;
599#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
600 __attribute__ ((fallthrough));
601#endif
602 case 3:
603 t4 = 0.;
604#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
605 __attribute__ ((fallthrough));
606#endif
607 case 4:
608 t5 = 0.;
609#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
610 __attribute__ ((fallthrough));
611#endif
612 case 5:
613 t6 = 0.;
614 t6l = 0.;
615#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
616 __attribute__ ((fallthrough));
617#endif
618 case 6:
619 t7 = 0.;
620#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
621 __attribute__ ((fallthrough));
622#endif
623 case 7:
624 case -1: // Use the max PN order, move if higher orders implemented
625 break;
626 case 8:
627 XLALPrintError("XLAL Error - %s: Not supported for requested PN order\n", __func__);
629 break;
630 default:
631 XLALPrintError("XLAL Error - %s: Unknown PN order in switch\n", __func__);
633 }
634
635 tC = -tN / v8 * (1.
636 + t2 * v2
637 + t3 * v3
638 + t4 * v4
639 + t5 * v5
640 + (t6 + t6l*log(16*v2)) * v6
641 + t7 * v7
642 + t10 * v10
643 + t12 * v12);
644
645 return tC;
646}
647
648
649/*
650 * Set up the expnCoeffsTaylorT3 and expnFuncTaylorT3 structures for
651 * generating a TaylorT3 waveform.
652 *
653 * Inputs given in SI units.
654 */
656 expnCoeffsTaylorT3 *ak, /**< coefficients for TaylorT3 evolution [modified] */
657 expnFuncTaylorT3 *f, /**< functions for TaylorT3 evolution [modified] */
658 REAL8 deltaT, /**< sampling interval */
659 REAL8 m1, /**< mass of companion 1 */
660 REAL8 m2, /**< mass of companion 2 */
661 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
662 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
663 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
664 REAL8 f_min, /**< start frequency */
665 int O /**< twice post-Newtonian order */
666 )
667{
668 REAL8 eta, tn, chi1, chi2;
669
670 ak->t0 = 0;
671 ak->m1 = m1;
672 ak->m2 = m2;
673 ak->totalmass = ak->m1 + ak->m2;
674 eta = ak->eta = m1 * m2 / (ak->totalmass * ak->totalmass);
675 chi1 = ak->chi1 = m1/ak->totalmass;
676 chi2 = ak->chi2 = m2/ak->totalmass;
677 ak->totalmass *= LAL_G_SI / pow(LAL_C_SI, 3.0); /* convert totalmass from kilograms to seconds */
678
679 ak->f0 = f_min;
681 ak->fn = 1. / (2. * ak->samplinginterval);
682 ak->v0 = cbrt(LAL_PI * ak->totalmass * f_min);
683
692
701
702 /* Tidal co-efficients for phasing and frequency */
703 ak->pta10 = 0.;
704 ak->pta12 = 0.;
705 ak->fta10 = 0.;
706 ak->fta12 = 0.;
707 switch( tideO )
708 {
712 eta, chi1, lambda1)
715 eta, chi1, lambda1)
717#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
718 __attribute__ ((fallthrough));
719#endif
725#if __GNUC__ >= 7 && !defined __INTEL_COMPILER
726 __attribute__ ((fallthrough));
727#endif
729 break;
730 default:
731 XLALPrintError("XLAL Error - %s: Invalid tidal PN order %d\nSee LALSimInspiralTidalOrder enum in LALSimInspiralWaveformFlags.h for valid tidal orders.\n",
732 __func__, tideO );
734 break;
735 }
736
737 switch (O)
738 {
739 case 0:
742 break;
743 case 1:
744 XLALPrintError("XLAL Error - %s: PN approximant not supported for requested PN order\n", __func__);
746 break;
747 case 2:
750 break;
751 case 3:
754 break;
755 case 4:
758 break;
759 case 5:
762 break;
763 case 6:
766 break;
767 case 7:
768 case -1: // Use the highest PN order available, move if higher terms added
771 break;
772 case 8:
773 XLALPrintError("XLAL Error - %s: PN approximant not supported for requested PN order\n", __func__);
775 break;
776 default:
777 XLALPrintError("XLAL Error - %s: Unknown PN order in switch\n", __func__);
779 }
780
781 tn = XLALSimInspiralTaylorLength(deltaT, m1, m2, f_min, O);
782 ak->theta_lso = pow(tn, -0.125);
783
784 return XLAL_SUCCESS;
785}
786
787/**
788 * @addtogroup LALSimInspiralTaylorXX_c
789 * @{
790 * @name Routines for TaylorT3 Waveforms
791 * @sa
792 * Section IIID of Alessandra Buonanno, Bala R Iyer, Evan
793 * Ochsner, Yi Pan, and B S Sathyaprakash, "Comparison of post-Newtonian
794 * templates for compact binary inspiral signals in gravitational-wave
795 * detectors", Phys. Rev. D 80, 084043 (2009), arXiv:0907.0700v1
796 * @{
797 */
798
799/**
800 * Computes a post-Newtonian orbit using the Taylor T3 method.
801 */
803 REAL8TimeSeries **V, /**< post-Newtonian parameter [returned] */
804 REAL8TimeSeries **phi, /**< orbital phase [returned] */
805 REAL8 phiRef, /**< reference orbital phase (rad) */
806 REAL8 deltaT, /**< sampling interval (s) */
807 REAL8 m1, /**< mass of companion 1 (kg) */
808 REAL8 m2, /**< mass of companion 2 (kg) */
809 REAL8 f_min, /**< starting GW frequency (Hz) */
810 REAL8 fRef, /**< reference GW frequency (Hz) */
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 order */
815 )
816{
817
818 const UINT4 blocklen = 1024;
819 const REAL8 visco = sqrt(1.0/6.0);
820 REAL8 m = m1 + m2, VRef = 0.;
821 REAL8 nu = m1 * m2 / m / m;
822 m *= LAL_G_SI / pow(LAL_C_SI, 3.0); /* convert m from kilograms to seconds */
823 REAL8 tmptC, tC, c1, xmin, xmax, xacc, v, phase, fOld, t, td, temp, tempMin = 0, tempMax = 0;
824 REAL8 (*freqfunction)(REAL8, void *);
825 UINT4 j, len, idxRef = 0;
827 REAL8 f;
828 void *pars;
829
830 expnFuncTaylorT3 expnfunc;
832 FreqInFromChirptime timeIn;
833
834 /* allocate memory */
835
836 *V = XLALCreateREAL8TimeSeries("ORBITAL_FREQUENCY_PARAMETER", &t_end, 0.0, deltaT, &lalDimensionlessUnit,
837 blocklen);
838 *phi = XLALCreateREAL8TimeSeries("ORBITAL_PHASE", &t_end, 0.0, deltaT, &lalDimensionlessUnit, blocklen);
839 if (!V || !phi)
841
842
843 /* initialize expnCoeffsTaylorT3 and expnFuncTaylorT3 structures */
844 if (XLALSimInspiralTaylorT3Setup(&ak, &expnfunc, deltaT, m1, m2, lambda1,
845 lambda2, tideO, f_min, O))
847
848 tC = XLALSimInspiralChirpLength(m1, m2, lambda1, lambda2, tideO, f_min, O);
849 c1 = nu/(5.*m);
850
851 /*
852 * In Jan 2003 we realized that the tC determined as a sum of chirp
853 * times is not quite the tC that should enter the definition of Theta
854 * in the expression for the frequency as a function of time (see DIS3,
855 * 2000). This is because chirp times are obtained by inverting t(f).
856 * Rather tC should be obtained by solving the equation f0 - f(tC) = 0.
857 * This is what is implemented below.
858 */
859
860 timeIn.func = expnfunc.frequency3;
861 timeIn.ak = ak;
862 freqfunction = &XLALInspiralFrequency3Wrapper;
863 xmin = c1*tC/2.;
864 xmax = c1*tC*2.;
865 xacc = 1.e-6;
866 pars = (void*) &timeIn;
867 /* tc is the instant of coalescence */
868
869 /* we add 5 so that if tC is small then xmax
870 * is always greater than a given value (here 5)*/
871 xmax = c1*tC*3 + 5.;
872
873 /* for x in [xmin, xmax], we search the value which gives the max
874 * frequency. and keep the corresponding rootIn.xmin. */
875
876 for (tmptC = c1*tC/1000.; tmptC < xmax; tmptC+=c1*tC/1000.){
877 temp = XLALInspiralFrequency3Wrapper(tmptC , pars);
878 if (XLAL_IS_REAL8_FAIL_NAN(temp))
880 if (temp > tempMax) {
881 xmin = tmptC;
882 tempMax = temp;
883 }
884 if (temp < tempMin) {
885 tempMin = temp;
886 }
887 }
888
889 /* if we have found a value positive then everything should be fine in
890 * the BissectionFindRoot function */
891 if (tempMax > 0 && tempMin < 0){
892 tC = XLALDBisectionFindRoot (freqfunction, xmin, xmax, xacc, pars);
895 }
896 else{
897 XLALPrintError("Can't find good bracket for BisectionFindRoot");
899 }
900
901 tC /= c1;
902
903 /* start waveform generation */
904
905 t = 0.;
906 td = c1 * (tC - t);
907 phase = expnfunc.phasing3(td, &ak);
908 if (XLAL_IS_REAL8_FAIL_NAN(phase))
910 f = expnfunc.frequency3(td, &ak);
913
914 v = cbrt(f * LAL_PI * m);
915 (*V)->data->data[0] = v;
916 (*phi)->data->data[0] = phase;
917
918 j = 0;
919 while (1) {
920
921 /* make one step */
922
923 j++;
924 fOld = f;
925 t = j * deltaT;
926 td = c1 * (tC - t);
927 phase = expnfunc.phasing3(td, &ak);
928 if (XLAL_IS_REAL8_FAIL_NAN(phase))
930 f = expnfunc.frequency3(td, &ak);
933 v = cbrt(f * LAL_PI * m);
934
935 /* check termination conditions */
936
937 if (v >= visco) {
938 XLALPrintInfo("XLAL Info - %s: PN inspiral terminated at ISCO\n", __func__);
939 break;
940 }
941 if (f <= fOld) {
942 XLALPrintInfo("XLAL Info - %s: PN inspiral terminated when frequency stalled\n", __func__);
943 break;
944 }
945
946 /* save current values in vectors but first make sure we don't write past end of vectors */
947
948 if ( j >= (*V)->data->length ) {
949 if ( ! XLALResizeREAL8TimeSeries(*V, 0, (*V)->data->length + blocklen) )
951 if ( ! XLALResizeREAL8TimeSeries(*phi, 0, (*phi)->data->length + blocklen) )
953 }
954 (*V)->data->data[j] = v;
955 (*phi)->data->data[j] = phase;
956 }
957
958 /* make the correct length */
959
960 if ( ! XLALResizeREAL8TimeSeries(*V, 0, j) )
962 if ( ! XLALResizeREAL8TimeSeries(*phi, 0, j) )
964
965 /* adjust to correct time */
966
967 XLALGPSAdd(&(*phi)->epoch, -1.0*j*deltaT);
968 XLALGPSAdd(&(*V)->epoch, -1.0*j*deltaT);
969
970 /* Do a constant phase shift to get desired value of phiRef */
971 len = (*phi)->data->length;
972 /* For fRef==0, phiRef is phase of last sample */
973 if( fRef == 0. )
974 phiRef -= (*phi)->data->data[len-1];
975 /* For fRef==fmin, phiRef is phase of first sample */
976 else if( fRef == f_min )
977 phiRef -= (*phi)->data->data[0];
978 /* phiRef is phase when f==fRef */
979 else
980 {
981 VRef = pow(LAL_PI * LAL_G_SI*(m1+m2) * fRef, 1./3.) / LAL_C_SI;
982 j = 0;
983 do {
984 idxRef = j;
985 j++;
986 } while ((*V)->data->data[j] <= VRef);
987 phiRef -= (*phi)->data->data[idxRef];
988 }
989 for (j = 0; j < len; ++j)
990 (*phi)->data->data[j] += phiRef;
991
992 return (int)(*V)->data->length;
993}
994
995
996/**
997 * Driver routine to compute the post-Newtonian inspiral waveform.
998 *
999 * This routine allows the user to specify different pN orders
1000 * for phasing calcuation vs. amplitude calculations.
1001 */
1003 REAL8TimeSeries **hplus, /**< +-polarization waveform */
1004 REAL8TimeSeries **hcross, /**< x-polarization waveform */
1005 REAL8 phiRef, /**< reference orbital phase (rad) */
1006 REAL8 v0, /**< tail-term gauge choice (default = 1) */
1007 REAL8 deltaT, /**< sampling interval (s) */
1008 REAL8 m1, /**< mass of companion 1 (kg) */
1009 REAL8 m2, /**< mass of companion 2 (kg) */
1010 REAL8 f_min, /**< starting GW frequency (Hz) */
1011 REAL8 fRef, /**< reference GW frequency (Hz) */
1012 REAL8 r, /**< distance of source (m) */
1013 REAL8 i, /**< inclination of source (rad) */
1014 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
1015 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
1016 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
1017 int amplitudeO, /**< twice post-Newtonian amplitude order */
1018 int phaseO /**< twice post-Newtonian phase order */
1019 )
1020{
1022 REAL8TimeSeries *phi;
1023 int status;
1024 int n;
1025 n = XLALSimInspiralTaylorT3PNEvolveOrbit(&V, &phi, phiRef, deltaT,
1026 m1, m2, f_min, fRef, lambda1, lambda2, tideO, phaseO);
1027 if ( n < 0 )
1029 status = XLALSimInspiralPNPolarizationWaveforms(hplus, hcross, V, phi,
1030 v0, m1, m2, r, i, amplitudeO);
1033 if ( status < 0 )
1035 return n;
1036}
1037
1038/**
1039 * Driver routine to compute the -2 spin-weighted spherical harmonic modes
1040 * using TaylorT3 phasing.
1041 */
1043 UNUSED REAL8 v0, /**< tail-term gauge choice (default = 1) */
1044 REAL8 deltaT, /**< sampling interval (s) */
1045 REAL8 m1, /**< mass of companion 1 (kg) */
1046 REAL8 m2, /**< mass of companion 2 (kg) */
1047 REAL8 f_min, /**< starting GW frequency (Hz) */
1048 REAL8 fRef, /**< reference GW frequency (Hz) */
1049 REAL8 r, /**< distance of source (m) */
1050 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
1051 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
1052 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
1053 int amplitudeO, /**< twice post-Newtonian amplitude order */
1054 int phaseO, /**< twice post-Newtonian phase order */
1055 int lmax /**< generate all modes with l <= lmax */
1056 )
1057{
1058 SphHarmTimeSeries *hlm = NULL;
1059 /* The Schwarzschild ISCO frequency - for sanity checking fRef */
1060 REAL8 fISCO = pow(LAL_C_SI,3) / (pow(6.,3./2.)*LAL_PI*(m1+m2)*LAL_G_SI);
1061
1062 /* Sanity check fRef value */
1063 if( fRef < 0. )
1064 {
1065 XLALPrintError("XLAL Error - %s: fRef = %f must be >= 0\n",
1066 __func__, fRef);
1068 }
1069 if( fRef != 0. && fRef < f_min )
1070 {
1071 XLALPrintError("XLAL Error - %s: fRef = %f must be > fStart = %f\n",
1072 __func__, fRef, f_min);
1074 }
1075 if( fRef >= fISCO )
1076 {
1077 XLALPrintError("XLAL Error - %s: fRef = %f must be < Schwar. ISCO=%f\n",
1078 __func__, fRef, fISCO);
1080 }
1081
1083 REAL8TimeSeries *phi;
1084 int n;
1086 m1, m2, f_min, fRef, lambda1, lambda2, tideO, phaseO);
1087 if ( n < 0 )
1089 int m, l;
1091 for(l=2; l<=lmax; l++){
1092 for(m=-l; m<=l; m++){
1094 m1, m2, r, amplitudeO, l, m);
1095 if ( !hxx ){
1097 }
1098 hlm = XLALSphHarmTimeSeriesAddMode(hlm, hxx, l, m);
1100 }
1101 }
1104 return hlm;
1105}
1106
1107/**
1108 * Driver routine to compute the -2 spin-weighted spherical harmonic mode
1109 * using TaylorT3 phasing.
1110 */
1112 UNUSED REAL8 v0, /**< tail-term gauge choice (default = 1) */
1113 REAL8 deltaT, /**< sampling interval (s) */
1114 REAL8 m1, /**< mass of companion 1 (kg) */
1115 REAL8 m2, /**< mass of companion 2 (kg) */
1116 REAL8 f_min, /**< starting GW frequency (Hz) */
1117 REAL8 fRef, /**< reference GW frequency (Hz) */
1118 REAL8 r, /**< distance of source (m) */
1119 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
1120 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
1121 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
1122 int amplitudeO, /**< twice post-Newtonian amplitude order */
1123 int phaseO, /**< twice post-Newtonian phase order */
1124 int l, /**< l index of mode */
1125 int m /**< m index of mode */
1126 )
1127{
1129 /* The Schwarzschild ISCO frequency - for sanity checking fRef */
1130 REAL8 fISCO = pow(LAL_C_SI,3) / (pow(6.,3./2.)*LAL_PI*(m1+m2)*LAL_G_SI);
1131
1132 /* Sanity check fRef value */
1133 if( fRef < 0. )
1134 {
1135 XLALPrintError("XLAL Error - %s: fRef = %f must be >= 0\n",
1136 __func__, fRef);
1138 }
1139 if( fRef != 0. && fRef < f_min )
1140 {
1141 XLALPrintError("XLAL Error - %s: fRef = %f must be > fStart = %f\n",
1142 __func__, fRef, f_min);
1144 }
1145 if( fRef >= fISCO )
1146 {
1147 XLALPrintError("XLAL Error - %s: fRef = %f must be < Schwar. ISCO=%f\n",
1148 __func__, fRef, fISCO);
1150 }
1151
1153 REAL8TimeSeries *phi;
1154 int n;
1156 m1, m2, f_min, fRef, lambda1, lambda2, tideO, phaseO);
1157 if ( n < 0 )
1160 m1, m2, r, amplitudeO, l, m);
1163 return hlm;
1164}
1165
1166/**
1167 * Driver routine to compute the post-Newtonian inspiral waveform.
1168 *
1169 * This routine uses the same pN order for phasing and amplitude
1170 * (unless the order is -1 in which case the highest available
1171 * order is used for both of these -- which might not be the same).
1172 *
1173 * Constant log term in amplitude set to 1. This is a gauge choice.
1174 */
1176 REAL8TimeSeries **hplus, /**< +-polarization waveform */
1177 REAL8TimeSeries **hcross, /**< x-polarization waveform */
1178 REAL8 phiRef, /**< reference orbital phase (rad) */
1179 REAL8 deltaT, /**< sampling interval (s) */
1180 REAL8 m1, /**< mass of companion 1 (kg) */
1181 REAL8 m2, /**< mass of companion 2 (kg) */
1182 REAL8 f_min, /**< starting GW frequency (Hz) */
1183 REAL8 fRef, /**< reference GW frequency (Hz) */
1184 REAL8 r, /**< distance of source (m) */
1185 REAL8 i, /**< inclination of source (rad) */
1186 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
1187 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
1188 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
1189 int O /**< twice post-Newtonian order */
1190 )
1191{
1192 /* set v0 to default value 1 */
1193 return XLALSimInspiralTaylorT3PNGenerator(hplus, hcross, phiRef, 1.0,
1194 deltaT, m1, m2, f_min, fRef, r, i, lambda1, lambda2,
1195 tideO, O, O);
1196}
1197
1198
1199/**
1200 * Driver routine to compute the restricted post-Newtonian inspiral waveform.
1201 *
1202 * This routine computes the phasing to the specified order, but
1203 * only computes the amplitudes to the Newtonian (quadrupole) order.
1204 *
1205 * Constant log term in amplitude set to 1. This is a gauge choice.
1206 */
1208 REAL8TimeSeries **hplus, /**< +-polarization waveform */
1209 REAL8TimeSeries **hcross, /**< x-polarization waveform */
1210 REAL8 phiRef, /**< reference orbital phase (rad) */
1211 REAL8 deltaT, /**< sampling interval (s) */
1212 REAL8 m1, /**< mass of companion 1 (kg) */
1213 REAL8 m2, /**< mass of companion 2 (kg) */
1214 REAL8 f_min, /**< starting GW frequency (Hz) */
1215 REAL8 fRef, /**< reference GW frequency (Hz) */
1216 REAL8 r, /**< distance of source (m)*/
1217 REAL8 i, /**< inclination of source (rad) */
1218 REAL8 lambda1, /**< (tidal deformability of body 1)/(mass of body 1)^5 */
1219 REAL8 lambda2, /**< (tidal deformability of body 2)/(mass of body 2)^5 */
1220 LALSimInspiralTidalOrder tideO, /**< flag to control spin and tidal effects */
1221 int O /**< twice post-Newtonian phase order */
1222 )
1223{
1224 /* use Newtonian order for amplitude */
1225 /* set v0 to default value 1 */
1226 return XLALSimInspiralTaylorT3PNGenerator(hplus, hcross, phiRef, 1.0,
1227 deltaT, m1, m2, f_min, fRef, r, i, lambda1, lambda2,
1228 tideO, 0, O);
1229}
1230
1231/** @} */
1232/** @} */
const double c1
REAL8 XLALSimInspiralTaylorLength(REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 f_min, int O)
static REAL8 UNUSED XLALSimInspiralTaylorT2Timing_12PNTidalCoeff(REAL8 eta, REAL8 chi, REAL8 lambda)
static REAL8 UNUSED XLALSimInspiralTaylorT3Phasing_12PNTidalCoeff(REAL8 eta, REAL8 chi, REAL8 lambda)
static REAL8 UNUSED XLALSimInspiralTaylorT2Timing_6PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Frequency_4PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Phasing_4PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Frequency_12PNTidalCoeff(REAL8 eta, REAL8 chi, REAL8 lambda)
static REAL8 UNUSED XLALSimInspiralTaylorT3Frequency_3PNCoeff(REAL8 UNUSED eta)
static REAL8 UNUSED XLALSimInspiralTaylorT2Timing_2PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Phasing_3PNCoeff(REAL8 UNUSED eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Phasing_6PNLogCoeff(REAL8 UNUSED eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Frequency_7PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT2Timing_5PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Phasing_6PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Phasing_7PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Phasing_0PNCoeff(REAL8 eta)
Computes the PN Coefficients for using in the TaylorT3 phasing equation.
static REAL8 UNUSED XLALSimInspiralTaylorT2Timing_4PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Frequency_2PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT2Timing_3PNCoeff(REAL8 UNUSED eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Phasing_2PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Frequency_5PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Phasing_5PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Frequency_0PNCoeff(REAL8 totalmass)
Computes the PN Coefficients for using in the TaylorT3 frequency equation.
static REAL8 UNUSED XLALSimInspiralTaylorT2Timing_7PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Phasing_10PNTidalCoeff(REAL8 chi, REAL8 lambda)
static REAL8 UNUSED XLALSimInspiralTaylorT3Frequency_6PNLogCoeff(REAL8 UNUSED eta)
static REAL8 UNUSED XLALSimInspiralTaylorT2Timing_10PNTidalCoeff(REAL8 chi, REAL8 lambda)
static REAL8 UNUSED XLALSimInspiralTaylorT2Timing_0PNCoeff(REAL8 totalmass, REAL8 eta)
Computes the PN Coefficients for using in the TaylorT2 timing equation.
static REAL8 UNUSED XLALSimInspiralTaylorT3Frequency_6PNCoeff(REAL8 eta)
static REAL8 UNUSED XLALSimInspiralTaylorT3Frequency_10PNTidalCoeff(REAL8 chi, REAL8 lambda)
static REAL8 UNUSED XLALSimInspiralTaylorT2Timing_6PNLogCoeff(REAL8 UNUSED eta)
static REAL8 XLALSimInspiralFrequency3_0PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralFrequency3_2PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static int XLALSimInspiralTaylorT3Setup(expnCoeffsTaylorT3 *ak, expnFuncTaylorT3 *f, REAL8 deltaT, REAL8 m1, REAL8 m2, REAL8 lambda1, REAL8 lambda2, LALSimInspiralTidalOrder tideO, REAL8 f_min, int O)
static REAL8 XLALInspiralFrequency3Wrapper(REAL8 tC, void *pars)
static REAL8 XLALSimInspiralPhasing3_5PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralPhasing3_4PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralPhasing3_3PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralFrequency3_5PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralPhasing3_2PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralFrequency3_3PN(REAL8 td, expnCoeffsTaylorT3 *ak)
REAL8() SimInspiralPhasing3(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralFrequency3_4PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralPhasing3_6PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralPhasing3_0PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralChirpLength(REAL8 m1, REAL8 m2, REAL8 lambda1, REAL8 lambda2, LALSimInspiralTidalOrder tideO, REAL8 f_min, int O)
static REAL8 XLALSimInspiralPhasing3_7PN(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralFrequency3_7PN(REAL8 td, expnCoeffsTaylorT3 *ak)
REAL8() SimInspiralFrequency3(REAL8 td, expnCoeffsTaylorT3 *ak)
static REAL8 XLALSimInspiralFrequency3_6PN(REAL8 td, expnCoeffsTaylorT3 *ak)
Module containing the energy and flux functions for waveform generation.
int l
Definition: bh_qnmode.c:135
double i
Definition: bh_ringdown.c:118
double theta
Definition: bh_sphwf.c:118
#define __attribute__(x)
REAL8 XLALDBisectionFindRoot(REAL8(*y)(REAL8, void *), REAL8 xmin, REAL8 xmax, REAL8 xacc, void *params)
#define LAL_C_SI
#define LAL_PI
#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 * XLALSimInspiralTaylorT3PNModes(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 TaylorT3 phasing.
COMPLEX16TimeSeries * XLALSimInspiralTaylorT3PNMode(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 TaylorT3 phasing.
int XLALSimInspiralTaylorT3PN(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.
int XLALSimInspiralTaylorT3PNGenerator(REAL8TimeSeries **hplus, REAL8TimeSeries **hcross, REAL8 phiRef, 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 XLALSimInspiralTaylorT3PNRestricted(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.
int XLALSimInspiralTaylorT3PNEvolveOrbit(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)
Computes a post-Newtonian orbit using the Taylor T3 method.
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.
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_REAL8(...)
#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
#define XLAL_IS_REAL8_FAIL_NAN(val)
XLAL_SUCCESS
XLAL_EFAULT
XLAL_EFUNC
XLAL_EMAXITER
XLAL_EINVAL
LIGOTimeGPS * XLALGPSAdd(LIGOTimeGPS *epoch, REAL8 dt)
string status
REAL8(* func)(REAL8 tC, expnCoeffsTaylorT3 *ak)
REAL8Sequence * data
Structure to carry a collection of spherical harmonic modes in COMPLEX16 time series.
SimInspiralFrequency3 * frequency3
SimInspiralPhasing3 * phasing3
double V
Definition: unicorn.c:25
double f_min
Definition: unicorn.c:22
double deltaT
Definition: unicorn.c:24