GADGET-4
driftfac.h
Go to the documentation of this file.
1/*******************************************************************************
2 * \copyright This file is part of the GADGET4 N-body/SPH code developed
3 * \copyright by Volker Springel. Copyright (C) 2014-2020 by Volker Springel
4 * \copyright (vspringel@mpa-garching.mpg.de) and all contributing authors.
5 *******************************************************************************/
6
12#ifndef DRIFTFAC_H
13#define DRIFTFAC_H
14
15#include <gsl/gsl_integration.h>
16#include <gsl/gsl_math.h>
17#include <math.h>
18#include <mpi.h>
19#include <stdio.h>
20#include <stdlib.h>
21#include <string.h>
22
23#include "../data/dtypes.h"
24#include "../main/main.h"
25#include "gadgetconfig.h"
26
28{
29 public:
30 void init_drift_table(void);
31 double get_drift_factor(integertime time0, integertime time1);
32 double get_gravkick_factor(integertime time0, integertime time1);
33 double get_hydrokick_factor(integertime time0, integertime time1);
35 double get_comoving_distance_for_scalefactor(double ascale);
36 double get_scalefactor_for_comoving_distance(double dist);
38
39 static double hubble_function(double a)
40 {
41 double hubble_a = All.Omega0 / (a * a * a) + (1 - All.Omega0 - All.OmegaLambda) / (a * a) + All.OmegaLambda;
42
43 hubble_a = All.Hubble * sqrt(hubble_a);
44
45 return hubble_a;
46 }
47
48 private:
49#define DRIFT_TABLE_LENGTH 1000
50
52 double DriftTable[DRIFT_TABLE_LENGTH];
53
55 double GravKickTable[DRIFT_TABLE_LENGTH];
56
58 double HydroKickTable[DRIFT_TABLE_LENGTH];
59
60 double logTimeBegin;
61 double logTimeMax;
62
63 static double drift_integ(double a, void *param)
64 {
65 double h = hubble_function(a);
66
67 return 1 / (h * a * a * a);
68 }
69
70 static double gravkick_integ(double a, void *param)
71 {
72 double h = hubble_function(a);
73
74 return 1 / (h * a * a);
75 }
76
77 static double hydrokick_integ(double a, void *param)
78 {
79 double h = hubble_function(a);
80
81 return 1 / (h * pow(a, 3 * GAMMA_MINUS1) * a);
82 }
83};
84
85extern driftfac Driftfac;
86
87#endif
global_data_all_processes All
Definition: main.cc:40
static double hubble_function(double a)
Definition: driftfac.h:39
double get_drift_factor(integertime time0, integertime time1)
Definition: driftfac.cc:68
double get_comoving_distance_for_scalefactor(double ascale)
Definition: driftfac.cc:199
double get_gravkick_factor(integertime time0, integertime time1)
Definition: driftfac.cc:109
double get_hydrokick_factor(integertime time0, integertime time1)
Definition: driftfac.cc:150
double get_scalefactor_for_comoving_distance(double dist)
Definition: driftfac.cc:209
double get_comoving_distance(integertime time0)
Definition: driftfac.cc:191
integertime get_gravkick_factor_inverse(double fac)
Definition: driftfac.cc:220
void init_drift_table(void)
Definition: driftfac.cc:26
#define GAMMA_MINUS1
Definition: constants.h:110
int integertime
Definition: constants.h:331
driftfac Driftfac
Definition: main.cc:41
#define DRIFT_TABLE_LENGTH
Definition: driftfac.h:49
expr pow(half base, half exp)
Definition: half.hpp:2803
expr sqrt(half arg)
Definition: half.hpp:2777