CLASS MANUAL
transfer.h
Go to the documentation of this file.
1 
3 #ifndef __TRANSFER__
4 #define __TRANSFER__
5 
6 #include "nonlinear.h"
7 #include "hyperspherical.h"
8 #include <sys/shm.h>
9 #include <sys/stat.h>
10 #include "errno.h"
11 
12 /* macro: test if index_tt is in the range between index and index+num, while the flag is true */
13 #define _index_tt_in_range_(index,num,flag) (flag == _TRUE_) && (index_tt >= index) && (index_tt < index+num)
14 
38 struct transfers {
39 
45 
46  double lcmb_rescale;
49  double lcmb_tilt;
52  double lcmb_pivot;
58  short has_nz_file;
60  FileName nz_file_name;
61  int nz_size;
62  double * nz_z;
63  double * nz_nz;
64  double * nz_ddnz;
68  FileName nz_evo_file_name;
70  double * nz_evo_z;
71  double * nz_evo_nz;
72  double * nz_evo_dlog_nz;
73  double * nz_evo_dd_dlog_nz;
76 
80 
81  short has_cls;
84 
88 
89  int md_size;
94  int index_tt_e;
95  int index_tt_b;
110  int * tt_size;
113 
117 
118  int ** l_size_tt;
120  int * l_size;
124  int * l;
126  //int * l_size_bessel; /**< for each wavenumber, maximum value of l at which bessel functions must be evaluated */
127 
131 
135 
136  size_t q_size;
138  double * q;
140  double ** k;
145 
149 
150  double ** transfer;
153 
157 
162  ErrorMsg error_message;
165 };
166 
175 
179 
180  HyperInterpStruct HIS;
184  HyperInterpStruct * pBIS;
186  int l_size;
189 
193 
194  int tau_size;
200  double * sources;
205  double * tau0_minus_tau;
206  double * w_trapz;
207  double * chi;
211  double * cscKgen;
212  double * cotKgen;
215 
219 
220  double K;
221  int sgnK;
224 
227 };
228 
236 typedef enum {SCALAR_TEMPERATURE_0,
237  SCALAR_TEMPERATURE_1,
238  SCALAR_TEMPERATURE_2,
239  SCALAR_POLARISATION_E,
240  VECTOR_TEMPERATURE_1,
241  VECTOR_TEMPERATURE_2,
242  VECTOR_POLARISATION_E,
243  VECTOR_POLARISATION_B,
244  TENSOR_TEMPERATURE_2,
245  TENSOR_POLARISATION_E,
246  TENSOR_POLARISATION_B,
247  NC_RSD} radial_function_type;
248 
249 enum Hermite_Interpolation_Order {HERMITE3, HERMITE4, HERMITE6};
250 
251 /*************************************************************************************************************/
252 /* @cond INCLUDE_WITH_DOXYGEN */
253 /*
254  * Boilerplate for C++
255  */
256 #ifdef __cplusplus
257 extern "C" {
258 #endif
259 
261  struct transfers * ptr,
262  int index_md,
263  int index_ic,
264  int index_type,
265  int index_l,
266  double q,
267  double * ptransfer_local
268  );
269 
270  int transfer_init(
271  struct precision * ppr,
272  struct background * pba,
273  struct thermo * pth,
274  struct perturbs * ppt,
275  struct nonlinear * pnl,
276  struct transfers * ptr
277  );
278 
279  int transfer_free(
280  struct transfers * ptr
281  );
282 
284  struct precision * ppr,
285  struct perturbs * ppt,
286  struct transfers * ptr,
287  double q_period,
288  double K,
289  int sgnK
290  );
291 
292  int transfer_perturbation_copy_sources_and_nl_corrections(
293  struct perturbs * ppt,
294  struct nonlinear * pnl,
295  struct transfers * ptr,
296  double *** sources
297  );
298 
299  int transfer_perturbation_source_spline(
300  struct perturbs * ppt,
301  struct transfers * ptr,
302  double *** sources,
303  double *** sources_spline
304  );
305 
306  int transfer_perturbation_sources_free(
307  struct perturbs * ppt,
308  struct nonlinear * pnl,
309  struct transfers * ptr,
310  double *** sources
311  );
312 
313  int transfer_perturbation_sources_spline_free(
314  struct perturbs * ppt,
315  struct transfers * ptr,
316  double *** sources_spline
317  );
318 
320  struct precision * ppr,
321  struct perturbs * ppt,
322  struct transfers * ptr
323  );
324 
326  struct precision * ppr,
327  struct perturbs * ppt,
328  struct transfers * ptr,
329  double q_period,
330  double K,
331  int sgnK
332  );
333 
334  int transfer_get_q_list_v1(
335  struct precision * ppr,
336  struct perturbs * ppt,
337  struct transfers * ptr,
338  double q_period,
339  double K,
340  int sgnK
341  );
342 
344  struct perturbs * ppt,
345  struct transfers * ptr,
346  double K
347  );
348 
350  struct perturbs * ppt,
351  struct transfers * ptr,
352  int ** tp_of_tt
353  );
354 
355  int transfer_free_source_correspondence(
356  struct transfers * ptr,
357  int ** tp_of_tt
358  );
359 
360  int transfer_source_tau_size_max(
361  struct precision * ppr,
362  struct background * pba,
363  struct perturbs * ppt,
364  struct transfers * ptr,
365  double tau_rec,
366  double tau0,
367  int * tau_size_max
368  );
369 
371  struct precision * ppr,
372  struct background * pba,
373  struct perturbs * ppt,
374  struct transfers * ptr,
375  double tau_rec,
376  double tau0,
377  int index_md,
378  int index_tt,
379  int * tau_size
380  );
381 
383  struct precision * ppr,
384  struct background * pba,
385  struct perturbs * ppt,
386  struct transfers * ptr,
387  int ** tp_of_tt,
388  int index_q,
389  int tau_size_max,
390  double tau_rec,
391  double *** sources,
392  double *** sources_spline,
393  struct transfer_workspace * ptw
394  );
395 
396  int transfer_radial_coordinates(
397  struct transfers * ptr,
398  struct transfer_workspace * ptw,
399  int index_md,
400  int index_q
401  );
402 
404  struct perturbs * ppt,
405  struct transfers * ptr,
406  int index_q,
407  int index_md,
408  int index_ic,
409  int index_type,
410  double * sources,
411  double * source_spline,
412  double * interpolated_sources
413  );
414 
415  int transfer_sources(
416  struct precision * ppr,
417  struct background * pba,
418  struct perturbs * ppt,
419  struct transfers * ptr,
420  double * interpolated_sources,
421  double tau_rec,
422  int index_q,
423  int index_md,
424  int index_tt,
425  double * sources,
426  double * tau0_minus_tau,
427  double * delta_tau,
428  int * tau_size_out
429  );
430 
432  struct precision * ppr,
433  struct perturbs * ppt,
434  struct transfers * ptr,
435  int bin,
436  double z,
437  double * selection);
438 
440  struct transfers * ptr,
441  double z,
442  double * dNdz,
443  double * dln_dNdz_dz);
444 
446  struct precision * ppr,
447  struct background * pba,
448  struct perturbs * ppt,
449  struct transfers * ptr,
450  int bin,
451  double * tau0_minus_tau,
452  int tau_size);
453 
455  struct precision * ppr,
456  struct background * pba,
457  struct perturbs * ppt,
458  struct transfers * ptr,
459  int bin,
460  double tau0,
461  double * tau0_minus_tau,
462  int tau_size);
463 
465  struct precision * ppr,
466  struct background * pba,
467  struct perturbs * ppt,
468  struct transfers * ptr,
469  int bin,
470  double * tau0_minus_tau,
471  int tau_size,
472  int index_md,
473  double tau0,
474  double * interpolated_sources,
475  double * sources);
476 
478  struct precision * ppr,
479  struct background * pba,
480  struct perturbs * ppt,
481  struct transfers * ptr,
482  int bin,
483  double * tau_min,
484  double * tau_mean,
485  double * tau_max);
486 
488  struct precision * ppr,
489  struct background * pba,
490  struct perturbs * ppt,
491  struct transfers * ptr,
492  double * selection,
493  double * tau0_minus_tau,
494  double * delta_tau,
495  int tau_size,
496  double * pvecback,
497  double tau0,
498  int bin);
499 
501  struct transfer_workspace * ptw,
502  struct precision * ppr,
503  struct perturbs * ppt,
504  struct transfers * ptr,
505  int index_q,
506  int index_md,
507  int index_ic,
508  int index_tt,
509  int index_l,
510  double l,
511  double q_max_bessel,
512  radial_function_type radial_type
513  );
514 
515  int transfer_use_limber(
516  struct precision * ppr,
517  struct perturbs * ppt,
518  struct transfers * ptr,
519  double q_max_bessel,
520  int index_md,
521  int index_tt,
522  double q,
523  double l,
524  short * use_limber
525  );
526 
527  int transfer_integrate(
528  struct perturbs * ppt,
529  struct transfers * ptr,
530  struct transfer_workspace *ptw,
531  int index_q,
532  int index_md,
533  int index_tt,
534  double l,
535  int index_l,
536  double q,
537  radial_function_type radial_type,
538  double * trsf
539  );
540 
541  int transfer_limber(
542  struct transfers * ptr,
543  struct transfer_workspace * ptw,
544  int index_md,
545  int index_q,
546  double l,
547  double q,
548  radial_function_type radial_type,
549  double * trsf
550  );
551 
553  struct transfers * ptr,
554  double * tau0_minus_tau,
555  double * sources,
556  int tau_size,
557  double tau0_minus_tau_limber,
558  double * S
559  );
560 
561  int transfer_limber2(
562  int tau_size,
563  struct transfers * ptr,
564  int index_md,
565  int index_q,
566  double l,
567  double q,
568  double * tau0_minus_tau,
569  double * sources,
570  radial_function_type radial_type,
571  double * trsf
572  );
573 
574  int transfer_can_be_neglected(
575  struct precision * ppr,
576  struct perturbs * ppt,
577  struct transfers * ptr,
578  int index_md,
579  int index_ic,
580  int index_tt,
581  double ra_rec,
582  double q,
583  double l,
584  short * neglect
585  );
586 
587  int transfer_late_source_can_be_neglected(
588  struct precision * ppr,
589  struct perturbs * ppt,
590  struct transfers * ptr,
591  int index_md,
592  int index_tt,
593  double l,
594  short * neglect);
595 
596  int transfer_select_radial_function(
597  struct perturbs * ppt,
598  struct transfers * ptr,
599  int index_md,
600  int index_tt,
601  radial_function_type *radial_type
602  );
603 
604  int transfer_radial_function(
605  struct transfer_workspace * ptw,
606  struct perturbs * ppt,
607  struct transfers * ptr,
608  double k,
609  int index_q,
610  int index_l,
611  int x_size,
612  double * radial_function,
613  radial_function_type radial_type
614  );
615 
616  int transfer_init_HIS_from_bessel(
617  struct transfers * ptr,
618  HyperInterpStruct *pHIS
619  );
620 
621  int transfer_global_selection_read(
622  struct transfers * ptr
623  );
624 
625  int transfer_workspace_init(
626  struct transfers * ptr,
627  struct precision * ppr,
628  struct transfer_workspace **ptw,
629  int perturb_tau_size,
630  int tau_size_max,
631  double K,
632  int sgnK,
633  double tau0_minus_tau_cut,
634  HyperInterpStruct * pBIS
635  );
636 
637  int transfer_workspace_free(
638  struct transfers * ptr,
639  struct transfer_workspace *ptw
640  );
641 
642  int transfer_update_HIS(
643  struct precision * ppr,
644  struct transfers * ptr,
645  struct transfer_workspace * ptw,
646  int index_q,
647  double tau0
648  );
649 
650  int transfer_get_lmax(int (*get_xmin_generic)(int sgnK,
651  int l,
652  double nu,
653  double xtol,
654  double phiminabs,
655  double *x_nonzero,
656  int *fevals),
657  int sgnK,
658  double nu,
659  int *lvec,
660  int lsize,
661  double phiminabs,
662  double xmax,
663  double xtol,
664  int *index_l_left,
665  int *index_l_right,
666  ErrorMsg error_message);
667 
668 #ifdef __cplusplus
669 }
670 #endif
671 
672 #endif
673 /* @endcond */
int tau_size
Definition: transfer.h:194
Definition: background.h:25
double selection_magnification_bias[_SELECTION_NUM_MAX_]
Definition: transfer.h:56
int index_tt_e
Definition: transfer.h:94
int transfer_compute_for_each_l(struct transfer_workspace *ptw, struct precision *ppr, struct perturbs *ppt, struct transfers *ptr, int index_q, int index_md, int index_ic, int index_tt, int index_l, double l, double q_max_bessel, radial_function_type radial_type)
Definition: transfer.c:3435
int transfer_lensing_sampling(struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, int bin, double tau0, double *tau0_minus_tau, int tau_size)
Definition: transfer.c:3111
short has_cls
Definition: transfer.h:81
HyperInterpStruct HIS
Definition: transfer.h:180
int transfer_get_l_list(struct precision *ppr, struct perturbs *ppt, struct transfers *ptr)
Definition: transfer.c:786
int transfer_selection_times(struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, int bin, double *tau_min, double *tau_mean, double *tau_max)
Definition: transfer.c:3235
FileName nz_evo_file_name
Definition: transfer.h:68
int tau_size_max
Definition: transfer.h:195
int index_tt_nc_lens
Definition: transfer.h:103
int transfer_free(struct transfers *ptr)
Definition: transfer.c:418
int transfer_functions_at_q(struct transfers *ptr, int index_md, int index_ic, int index_tt, int index_l, double q, double *transfer_function)
Definition: transfer.c:61
FileName nz_file_name
Definition: transfer.h:60
Definition: perturbations.h:95
double * chi
Definition: transfer.h:207
int transfer_sources(struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, double *interpolated_sources, double tau_rec, int index_q, int index_md, int index_tt, double *sources, double *tau0_minus_tau, double *w_trapz, int *tau_size_out)
Definition: transfer.c:2058
int transfer_get_source_correspondence(struct perturbs *ppt, struct transfers *ptr, int **tp_of_tt)
Definition: transfer.c:1293
double * nz_evo_z
Definition: transfer.h:70
Definition: thermodynamics.h:58
int transfer_interpolate_sources(struct perturbs *ppt, struct transfers *ptr, int index_q, int index_md, int index_ic, int index_type, double *pert_source, double *pert_source_spline, double *interpolated_sources)
Definition: transfer.c:1974
int transfer_integrate(struct perturbs *ppt, struct transfers *ptr, struct transfer_workspace *ptw, int index_q, int index_md, int index_tt, double l, int index_l, double k, radial_function_type radial_type, double *trsf)
Definition: transfer.c:3621
double selection_bias[_SELECTION_NUM_MAX_]
Definition: transfer.h:55
double * nz_nz
Definition: transfer.h:63
int index_tt_nc_g5
Definition: transfer.h:108
double ** transfer
Definition: transfer.h:150
int index_tt_nc_g3
Definition: transfer.h:106
int index_tt_density
Definition: transfer.h:97
int nz_size
Definition: transfer.h:61
int * l_size
Definition: transfer.h:120
short has_nz_analytic
Definition: transfer.h:59
double * tau0_minus_tau
Definition: transfer.h:205
int index_q_flat_approximation
Definition: transfer.h:142
int index_tt_lcmb
Definition: transfer.h:96
double * w_trapz
Definition: transfer.h:206
int * l
Definition: transfer.h:124
int l_size
Definition: transfer.h:186
int transfer_dNdz_analytic(struct transfers *ptr, double z, double *dNdz, double *dln_dNdz_dz)
Definition: transfer.c:3000
double tau0_minus_tau_cut
Definition: transfer.h:225
int index_tt_nc_g1
Definition: transfer.h:104
double * nz_evo_dlog_nz
Definition: transfer.h:72
HyperInterpStruct * pBIS
Definition: transfer.h:184
short has_nz_evo_analytic
Definition: transfer.h:67
#define _SELECTION_NUM_MAX_
Definition: perturbations.h:68
ErrorMsg error_message
Definition: transfer.h:162
int transfer_compute_for_each_q(struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, int **tp_of_tt, int index_q, int tau_size_max, double tau_rec, double ***pert_sources, double ***pert_sources_spline, struct transfer_workspace *ptw)
Definition: transfer.c:1659
double * cotKgen
Definition: transfer.h:212
short initialise_HIS_cache
Definition: transfer.h:158
int transfer_get_k_list(struct perturbs *ppt, struct transfers *ptr, double K)
Definition: transfer.c:1216
int index_tt_t2
Definition: transfer.h:93
int transfer_limber(struct transfers *ptr, struct transfer_workspace *ptw, int index_md, int index_q, double l, double q, radial_function_type radial_type, double *trsf)
Definition: transfer.c:3802
int index_tt_t0
Definition: transfer.h:91
int transfer_limber_interpolate(struct transfers *ptr, double *tau0_minus_tau, double *sources, int tau_size, double tau0_minus_tau_limber, double *S)
Definition: transfer.c:3959
double lcmb_pivot
Definition: transfer.h:52
double K
Definition: transfer.h:220
double lcmb_tilt
Definition: transfer.h:49
int index_tt_rsd
Definition: transfer.h:100
int transfer_selection_function(struct precision *ppr, struct perturbs *ppt, struct transfers *ptr, int bin, double z, double *selection)
Definition: transfer.c:2855
int transfer_indices_of_transfers(struct precision *ppr, struct perturbs *ppt, struct transfers *ptr, double q_period, double K, int sgnK)
Definition: transfer.c:475
int transfer_selection_sampling(struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, int bin, double *tau0_minus_tau, int tau_size)
Definition: transfer.c:3043
radial_function_type
Definition: transfer.h:236
Definition: transfer.h:38
double ** k
Definition: transfer.h:140
short has_nz_file
Definition: transfer.h:58
double * interpolated_sources
Definition: transfer.h:196
size_t q_size
Definition: transfer.h:136
int transfer_selection_compute(struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, double *selection, double *tau0_minus_tau, double *w_trapz, int tau_size, double *pvecback, double tau0, int bin)
Definition: transfer.c:3315
int index_tt_lensing
Definition: transfer.h:98
int HIS_allocated
Definition: transfer.h:182
int index_tt_d0
Definition: transfer.h:101
int nz_evo_size
Definition: transfer.h:69
short transfer_verbose
Definition: transfer.h:160
double * cscKgen
Definition: transfer.h:211
double angular_rescaling
Definition: transfer.h:128
double lcmb_rescale
Definition: transfer.h:46
int ** l_size_tt
Definition: transfer.h:118
int transfer_source_resample(struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, int bin, double *tau0_minus_tau, int tau_size, int index_md, double tau0, double *interpolated_sources, double *sources)
Definition: transfer.c:3168
int index_tt_d1
Definition: transfer.h:102
int index_tt_nc_g2
Definition: transfer.h:105
int index_tt_nc_g4
Definition: transfer.h:107
short has_nz_evo_file
Definition: transfer.h:66
Definition: common.h:345
short neglect_late_source
Definition: transfer.h:226
int md_size
Definition: transfer.h:89
int transfer_limber2(int tau_size, struct transfers *ptr, int index_md, int index_k, double l, double k, double *tau0_minus_tau, double *sources, radial_function_type radial_type, double *trsf)
Definition: transfer.c:4046
int sgnK
Definition: transfer.h:221
int transfer_get_q_list(struct precision *ppr, struct perturbs *ppt, struct transfers *ptr, double q_period, double K, int sgnK)
Definition: transfer.c:985
int l_size_max
Definition: transfer.h:122
int transfer_init(struct precision *ppr, struct background *pba, struct thermo *pth, struct perturbs *ppt, struct nonlinear *pnl, struct transfers *ptr)
Definition: transfer.c:115
double * sources
Definition: transfer.h:200
int transfer_source_tau_size(struct precision *ppr, struct background *pba, struct perturbs *ppt, struct transfers *ptr, double tau_rec, double tau0, int index_md, int index_tt, int *tau_size)
Definition: transfer.c:1475
Definition: transfer.h:174
double * nz_evo_dd_dlog_nz
Definition: transfer.h:73
int * tt_size
Definition: transfer.h:110
double * nz_ddnz
Definition: transfer.h:64
double * nz_evo_nz
Definition: transfer.h:71
double * nz_z
Definition: transfer.h:62
int index_tt_t1
Definition: transfer.h:92
int index_tt_b
Definition: transfer.h:95
Definition: nonlinear.h:20
double * q
Definition: transfer.h:138