changa 3.5
Loading...
Searching...
No Matches
cooling_cosmo.h
1#ifndef COOLING_COSMO_HINCLUDED
2#define COOLING_COSMO_HINCLUDED
3
4/*
5 * Cooling code for cosmology simulations.
6 * Originally written by James Wadsley, McMaster University for
7 * GASOLINE.
8 */
9
10/* Global consts */
11
12#include "param.h"
13
14#ifdef __cplusplus
15extern "C" {
16#endif
17
18#include "stiff.h"
19
20/* Constants */
21#define CL_B_gm (6.022e23*(938.7830/931.494))
22#define CL_k_Boltzmann 1.38066e-16
23#define CL_eV_erg 1.60219e-12
24#define CL_eV_per_K (CL_k_Boltzmann/CL_eV_erg)
25/*
26#define CL_RT_FLOAT float
27#define CL_RT_MIN 1e-38
28#define CL_RT_MIN FLT_MIN
29*/
30
31#define CL_RT_FLOAT double
32#define CL_RT_MIN 1e-100
33
34/*
35#define CL_RT_MIN DBL_MIN
36*/
37/*
38 * Work around for Dec ev6 flawed
39 * treatment of sub-normal numbers
40 */
41#define CL_MAX_NEG_EXP_ARG -500.
42
43#define CL_NMAXBYTETABLE 56000
44
46typedef struct CoolingParametersStruct {
47 int bIonNonEqm;
48 int nCoolingTable;
49 int bUV;
50 int bUVTableUsesTime;
51 int bDoIonOutput;
52 int bLowTCool;
53 int bSelfShield;
54 double dMassFracHelium;
55 double dCoolingTmin;
56 double dCoolingTmax;
57 } COOLPARAM;
58
60typedef struct CoolingParticleStruct {
61 double Y_HI,Y_HeI,Y_HeII; /* Abundance of ions */
62 double dLymanWerner; /* Flux of Lyman Werner radiation at the gas particle. Temporary CC */
63 } COOLPARTICLE;
64
66typedef struct {
67 double e,Total;
68 double HI,HII,HeI,HeII,HeIII;
69} PERBARYON;
70
72typedef struct {
73 double zTime;
74
75 double Rate_Phot_HI;
76 double Rate_Phot_HeI;
77 double Rate_Phot_HeII;
78
79 double Heat_Phot_HI;
80 double Heat_Phot_HeI;
81 double Heat_Phot_HeII;
82} UVSPECTRUM;
83
86typedef struct {
87 double Rate_Phot_HI;
88 double Rate_Phot_HeI;
89 double Rate_Phot_HeII;
90
91 double Heat_Phot_HI;
92 double Heat_Phot_HeI;
93 double Heat_Phot_HeII;
94
95 double Cool_Coll_HI;
96 double Cool_Coll_HeI;
97 double Cool_Coll_HeII;
98 double Cool_Diel_HeII;
99 double Cool_Comp;
100 double Tcmb;
101 double Cool_LowTFactor;
102
103} RATES_NO_T;
104
106typedef struct {
107 CL_RT_FLOAT Rate_Coll_HI;
108 CL_RT_FLOAT Rate_Coll_HeI;
109 CL_RT_FLOAT Rate_Coll_HeII;
110 CL_RT_FLOAT Rate_Radr_HII;
111 CL_RT_FLOAT Rate_Radr_HeII;
112 CL_RT_FLOAT Rate_Radr_HeIII;
113 CL_RT_FLOAT Rate_Diel_HeII;
114
115 CL_RT_FLOAT Cool_Brem_1;
116 CL_RT_FLOAT Cool_Brem_2;
117 CL_RT_FLOAT Cool_Radr_HII;
118 CL_RT_FLOAT Cool_Radr_HeII;
119 CL_RT_FLOAT Cool_Radr_HeIII;
120 CL_RT_FLOAT Cool_Line_HI;
121 CL_RT_FLOAT Cool_Line_HeI;
122 CL_RT_FLOAT Cool_Line_HeII;
123 CL_RT_FLOAT Cool_LowT;
124
125} RATES_T;
126
127typedef struct clDerivsDataStruct clDerivsData;
128
130
131typedef struct CoolingPKDStruct {
132 double z; /* Redshift */
133 double dTime;
134 /* Rates independent of Temperature */
135 RATES_NO_T R;
136 /* Table for Temperature dependent rates */
137 int nTable;
138 double TMin;
139 double TMax;
140 double TlnMin;
141 double TlnMax;
142 double rDeltaTln;
143 RATES_T *RT;
144
145 int nTableRead; /* number of Tables read from files */
146
147 int bUV;
148 int nUV;
149 UVSPECTRUM *UV;
150 int bUVTableUsesTime;
151 int bUVTableLinear;
152 int bLowTCool;
153 int bSelfShield;
154
155 double dGmPerCcUnit;
156 double dComovingGmPerCcUnit;
157 double dErgPerGmUnit;
158 double dSecUnit;
159 double dErgPerGmPerSecUnit;
160 double diErgPerGmUnit;
161 double dKpcUnit;
162 double Y_H;
163 double Y_He;
164 double Y_eMAX;
165/* Diagnostic */
166 int its;
167} COOL;
168
170typedef struct {
171 double T, Tln;
172 double Coll_HI;
173 double Coll_HeI;
174 double Coll_HeII;
175 double Radr_HII;
176 double Radr_HeII;
177 double Diel_HeII;
178 double Totr_HeII;
179 double Radr_HeIII;
180
181 double Phot_HI;
182 double Phot_HeI;
183 double Phot_HeII;
184} RATE;
185
187typedef struct {
188 double compton;
189 double bremHII;
190 double bremHeII;
191 double bremHeIII;
192 double radrecHII;
193 double radrecHeII;
194 double radrecHeIII;
195 double collionHI;
196 double collionHeI;
197 double collionHeII;
198 double dielrecHeII;
199 double lineHI;
200 double lineHeI;
201 double lineHeII;
202 double lowT;
203} COOL_ERGPERSPERGM;
204
205
207struct clDerivsDataStruct {
208 STIFF *IntegratorContext;
209 COOL *cl;
210 double rho,ExternalHeating,E,ZMetal;
211 RATE Rate;
212 PERBARYON Y;
213 double Y_Total0, Y_Total1;
214 int its; /* Debug */
215};
216
217
218COOL *CoolInit( );
219void CoolFinalize( COOL *cl );
220clDerivsData *CoolDerivsInit(COOL *cl);
221void CoolDerivsFinalize(clDerivsData *cld ) ;
222
223void clInitConstants( COOL *cl, double dGMPerCcunit, double dComovingGmPerCcUnit,
224 double dErgPerGmUnit, double dSecUnit, double dKpcUnit, COOLPARAM CoolParam);
225void clInitUV(COOL *cl, int nTableColumns, int nTableRows, double *dTableData );
226void clInitRatesTable( COOL *cl, double TMin, double TMax, int nTable );
227void CoolInitRatesTable( COOL *cl, COOLPARAM CoolParam);
228
229void clRatesTableError( COOL *cl );
230void clRatesRedshift( COOL *cl, double z, double dTime );
231double clHeatTotal ( COOL *cl, PERBARYON *Y, RATE *Rate );
232void clRates( COOL *cl, RATE *Rate, double T, double rho );
233double clCoolTotal( COOL *cl, PERBARYON *Y, RATE *Rate, double rho, double ZMetal );
234COOL_ERGPERSPERGM clTestCool ( COOL *cl, PERBARYON *Y, RATE *Rate, double rho );
235void clPrintCool( COOL *cl, PERBARYON *Y, RATE *Rate, double rho );
236void clPrintCoolFile( COOL *cl, PERBARYON *Y, RATE *Rate, double rho, FILE *fp );
237
238void clAbunds( COOL *cl, PERBARYON *Y, RATE *Rate, double rho);
239double clThermalEnergy( double Y_Total, double T );
240double clTemperature( double Y_Total, double E );
241double clRateCollHI( double T );
242double clRateCollHeI( double T );
243double clRateCollHeII( double T );
244double clRateRadrHII( double T );
245double clRateRadrHeII( double T );
246double clRateDielHeII( double T );
247double clRateRadrHeIII( double T );
248double clCoolBrem1( double T );
249double clCoolBrem2( double T );
250double clCoolRadrHII( double T );
251double clCoolRadrHeII( double T );
252double clCoolRadrHeIII( double T );
253double clCoolLineHI( double T );
254double clCoolLineHeI( double T );
255double clCoolLineHeII( double T );
256double clCoolLowT( double T );
257
258double clEdotInstant ( COOL *cl, PERBARYON *Y, RATE *Rate, double rho,
259 double ZMetal, double *EdotHeat, double *EdotCool );
260void clIntegrateEnergy(COOL *cl, clDerivsData *clData, PERBARYON *Y, double *E,
261 double ExternalHeating, double rho, double ZMetal, double dt );
262void clIntegrateEnergyDEBUG(COOL *cl, PERBARYON *Y, double *E,
263 double ExternalHeating, double rho, double ZMetal, double dt );
264
265
266void clDerivs(double x, const double *y, double *dHeat, double *dCool,
267 void *Data) ;
268
269int clJacobn( double x, const double y[], double dfdx[], double *dfdy, void *Data) ;
270
271void CoolAddParams( COOLPARAM *CoolParam, PRM );
272void CoolLogParams( COOLPARAM *CoolParam, FILE *fp );
273void CoolOutputArray( COOLPARAM *CoolParam, int, int *, char * );
274
275#define COOL_ARRAY0_EXT "HI"
276double COOL_ARRAY0(COOL *cl_, COOLPARTICLE *cp,double aa);
277#define COOL_ARRAY0( cl_, cp, aa ) ((cp)->Y_HI)
278double COOL_SET_ARRAY0(COOL *cl_, COOLPARTICLE *cp,double aa, double bb_val);
279#define COOL_SET_ARRAY0( cl_, cp, aa, bb_val ) ((cp)->Y_HI = (bb_val))
280
281#define COOL_ARRAY1_EXT "HeI"
282double COOL_ARRAY1(COOL *cl_, COOLPARTICLE *cp,double aa);
283#define COOL_ARRAY1( cl_, cp, aa ) ((cp)->Y_HeI)
284double COOL_SET_ARRAY1(COOL *cl_, COOLPARTICLE *cp,double aa, double bb_val);
285#define COOL_SET_ARRAY1( cl_, cp, aa, bb_val ) ((cp)->Y_HeI = (bb_val))
286
287#define COOL_ARRAY2_EXT "HeII"
288double COOL_ARRAY2(COOL *cl_, COOLPARTICLE *cp,double aa);
289#define COOL_ARRAY2( cl_, cp, aa ) ((cp)->Y_HeII)
290double COOL_SET_ARRAY2(COOL *cl_, COOLPARTICLE *cp,double aa, double bb_val);
291#define COOL_SET_ARRAY2( cl_, cp, aa, bb_val ) ((cp)->Y_HeII = (bb_val))
292
293#define COOL_ARRAY3_EXT "H2"
294double COOL_ARRAY3(COOL *cl, COOLPARTICLE *cp, double ZMetal);
295#define COOL_ARRAY3(cl_, cp, aa ) (0)
296double COOL_SET_ARRAY3(COOL *cl_, COOLPARTICLE *cp,double aa, double bb_val);
297#define COOL_SET_ARRAY3( cl_, cp, aa, bb_val ) (0)
298
299double COOL_EDOT( COOL *cl_, COOLPARTICLE *cp_, double ECode_, double rhoCode_, double ZMetal_, double *posCode_ );
300#define COOL_EDOT( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_) (CoolCodeWorkToErgPerGmPerSec( cl_, CoolEdotInstantCode( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_ )))
301
302double COOL_COOLING( COOL *cl_, COOLPARTICLE *cp_, double ECode_, double rhoCode_, double ZMetal_, double *posCode_ );
303#define COOL_COOLING( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_) (CoolCodeWorkToErgPerGmPerSec( cl_, CoolCoolingCode( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_ )))
304
305double COOL_HEATING( COOL *cl_, COOLPARTICLE *cp_, double ECode_, double rhoCode_, double ZMetal_, double *posCode_ );
306#define COOL_HEATING( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_) (CoolCodeWorkToErgPerGmPerSec( cl_, CoolHeatingCode( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_ )))
307
308void CoolPARTICLEtoPERBARYON(COOL *cl_, PERBARYON *Y, COOLPARTICLE *cp, double HTotal, double HeTotal);
309
310#define CoolPARTICLEtoPERBARYON(cl_, Y, cp) { \
311 (Y)->HI = (cp)->Y_HI; \
312 (Y)->HII = (cl_)->Y_H - (Y)->HI; \
313 (Y)->HeI = (cp)->Y_HeI; \
314 (Y)->HeII = (cp)->Y_HeII; \
315 (Y)->HeIII = (cl_)->Y_He - (Y)->HeI - (Y)->HeII; \
316 (Y)->e = (Y)->HII + (Y)->HeII + 2*(Y)->HeIII; \
317 (Y)->Total = (Y)->e + (cl_)->Y_H + (cl_)->Y_He; }
318
319void CoolPERBARYONtoPARTICLE(COOL *cl_, PERBARYON *Y, COOLPARTICLE *cp);
320
321#define CoolPERBARYONtoPARTICLE(cl_, Y, cp) { \
322 (cp)->Y_HI = (Y)->HI; \
323 (cp)->Y_HeI = (Y)->HeI; \
324 (cp)->Y_HeII = (Y)->HeII; }
325
326
327double CoolEnergyToTemperature( COOL *Cool, COOLPARTICLE *cp, double E, double fMetal );
328double CoolCodeEnergyToTemperature( COOL *Cool, COOLPARTICLE *cp, double E, double fMetal );
329
330/* Note: nod to cosmology (z parameter) unavoidable unless we want to access cosmo.[ch] from here */
331void CoolSetTime( COOL *Cool, double dTime, double z );
332
333double CoolCodeTimeToSeconds( COOL *Cool, double dCodeTime );
334
335#define CoolCodeTimeToSeconds( Cool, dCodeTime ) ((Cool)->dSecUnit*(dCodeTime))
336
337double CoolSecondsToCodeTime( COOL *Cool, double dTime );
338
339#define CoolSecondsToCodeTime( Cool, dTime ) ((dTime)/(Cool)->dSecUnit)
340
341double CoolCodeEnergyToErgPerGm( COOL *Cool, double dCodeEnergy );
342
343#define CoolCodeEnergyToErgPerGm( Cool, dCodeEnergy ) ((Cool)->dErgPerGmUnit*(dCodeEnergy))
344
345double CoolErgPerGmToCodeEnergy( COOL *Cool, double dEnergy );
346
347#define CoolErgPerGmToCodeEnergy( Cool, dEnergy ) ((Cool)->diErgPerGmUnit*(dEnergy))
348
349double CoolCodeWorkToErgPerGmPerSec( COOL *Cool, double dCodeWork );
350
351#define CoolCodeWorkToErgPerGmPerSec( Cool, dCodeWork ) ((Cool)->dErgPerGmPerSecUnit*(dCodeWork))
352
353double CoolErgPerGmPerSecToCodeWork( COOL *Cool, double dWork );
354
355#define CoolErgPerGmPerSecToCodeWork( Cool, dWork ) ((dWork)/(Cool)->dErgPerGmPerSecUnit)
356
357double CodeDensityToComovingGmPerCc( COOL *Cool, double dCodeDensity );
358
359#define CodeDensityToComovingGmPerCc( Cool, dCodeDensity ) ((Cool)->dComovingGmPerCcUnit*(dCodeDensity))
360
361void CoolIntegrateEnergy(COOL *cl, clDerivsData *cData, COOLPARTICLE *cp, double *E,
362 double ExternalHeating, double rho, double ZMetal, double tStep );
363
364void CoolIntegrateEnergyCode(COOL *cl, clDerivsData *cData, COOLPARTICLE *cp, double *E,
365 double ExternalHeating, double rho, double ZMetal, double *r, double tStep );
366
367void CoolDefaultParticleData( COOLPARTICLE *cp );
368
369void CoolInitEnergyAndParticleData( COOL *cl, COOLPARTICLE *cp, double *E, double dDensity, double dTemp, double fMetal );
370
371/* Deprecated */
372double CoolHeatingRate( COOL *cl, COOLPARTICLE *cp, double E, double dDensity, double ZMetal );
373
374double CoolEdotInstantCode(COOL *cl, COOLPARTICLE *cp, double ECode,
375 double rhoCode, double ZMetal, double *posCode );
376double CoolCoolingCode(COOL *cl, COOLPARTICLE *cp, double ECode,
377 double rhoCode, double ZMetal, double *posCode );
378double CoolHeatingCode(COOL *cl, COOLPARTICLE *cp, double ECode,
379 double rhoCode, double ZMetal, double *posCode );
380
381void CoolCodePressureOnDensitySoundSpeed( COOL *cl, COOLPARTICLE *cp, double uPred, double fDensity, double gamma, double gammam1, double *PoverRho, double *c );
382
383/* Note: gamma should be 5/3 for this to be consistent! */
384#define CoolCodePressureOnDensitySoundSpeed( cl__, cp__, uPred__, fDensity__, gamma__, gammam1__, PoverRho__, c__ ) { \
385 *(PoverRho__) = ((5./3.-1)*(uPred__)); \
386 *(c__) = sqrt((5./3.)*(*(PoverRho__))); }
387
388/*
389double CoolCodePressureOnDensity( COOL *cl, COOLPARTICLE *cp, double uPred, double fDensity, double gammam1 );
390
391#define CoolCodePressureOnDensity( cl, cp, uPred, fDensity, gammam1 ) ((gammam1)*(uPred))
392*/
393
394void CoolTableReadInfo( COOLPARAM *CoolParam, int cntTable, int *nTableColumns, char *suffix );
395
396void CoolTableRead( COOL *Cool, int nData, void *vData);
397
398#ifdef __cplusplus
399}
400#endif
401
402#endif
403
404
COOL_ERGPERSPERGM
return structure for clTestCool()
Definition cooling_cosmo.h:187
CoolingPKDStruct
Heating/Cooling context: parameters and tables.
Definition cooling_boley.h:83
CoolingParametersStruct
Input parameters for cooling.
Definition cooling_boley.h:56
CoolingParticleStruct
per-particle cooling data
Definition cooling_boley.h:68
PERBARYON
abundance of various species in particles/baryon
Definition cooling_boley.h:75
RATES_NO_T
structure to hold Temperature independent cooling and heating rates
Definition cooling_cosmo.h:86
RATES_T
structure to hold Temperature dependent cooling rates
Definition cooling_cosmo.h:106
RATE
Rate information for a given particle.
Definition cooling_cosmo.h:170
UVSPECTRUM
photoionization and heating rates from a uniform UV background
Definition cooling_cosmo.h:72
clDerivsDataStruct
context for calculating cooling derivatives
Definition cooling_boley.h:108
clDerivsDataStruct::cl
COOL * cl
pointer to cooling context
Definition cooling_boley.h:110
Generated by doxygen 1.13.2