changa 3.5
Loading...
Searching...
No Matches
cooling_metal.h
1
2#ifndef COOLING_METAL_HINCLUDED
3#define COOLING_METAL_HINCLUDED
4
5/*
6 * Cooling code for cosmology simulations.
7 * Originally written by James Wadsley and Sijing Shen, McMaster
8 * University for GASOLINE.
9 */
10
11/* Global consts */
12#include "param.h"
13#include "rpc/xdr.h"
14
15#ifdef __cplusplus
16extern "C" {
17#endif
18
19#include "stiff.h"
20
21/* Constants */
22#define CL_B_gm (6.022e23*(938.7830/931.494))/*Avegadro's Number * Mass_Hydrogen/Energy_AMU */
23#define CL_k_Boltzmann 1.38066e-16
24#define CL_eV_erg 1.60219e-12
25#define CL_eV_per_K (CL_k_Boltzmann/CL_eV_erg)
26/*
27#define CL_RT_FLOAT float
28#define CL_RT_MIN 1e-38
29#define CL_RT_MIN FLT_MIN
30*/
31
32#define CL_RT_FLOAT double
33#define CL_RT_MIN 1e-100
34
35/*
36#define CL_RT_MIN DBL_MIN
37*/
38/*
39 * Work around for Dec ev6 flawed
40 * treatment of sub-normal numbers
41 */
42#define CL_MAX_NEG_EXP_ARG -500.
43
44#define CL_NMAXBYTETABLE 56000
45#define MU_METAL 17.6003
46#define ZSOLAR 0.0130215
47
48typedef struct CoolingParametersStruct {
49 int bIonNonEqm;
50 int nCoolingTable;
51 int bUV;
52 int bMetal;
53 char *CoolInFile;
54 int bUVTableUsesTime;
55 int bDoIonOutput;
56 int bLowTCool;
57 int bSelfShield;
58 double dMassFracHelium;
59 double dCoolingTmin;
60 double dCoolingTmax;
61} COOLPARAM;
62
63typedef struct CoolingParticleStruct {
64 double f_HI,f_HeI,f_HeII;
65} COOLPARTICLE;
66
67typedef struct {
68 double e,Total;
69 double HI,HII,HeI,HeII,HeIII;
70} PERBARYON;
71
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
84typedef struct {
85 double Rate_Phot_HI;
86 double Rate_Phot_HeI;
87 double Rate_Phot_HeII;
88
89 double Heat_Phot_HI;
90 double Heat_Phot_HeI;
91 double Heat_Phot_HeII;
92
93 double Cool_Coll_HI;
94 double Cool_Coll_HeI;
95 double Cool_Coll_HeII;
96 double Cool_Diel_HeII;
97
98 double Cool_Comp;
99 double Tcmb;
100 double Cool_LowTFactor;
101
102} RATES_NO_T;
103
104typedef struct {
105 CL_RT_FLOAT Rate_Coll_HI;
106 CL_RT_FLOAT Rate_Coll_HeI;
107 CL_RT_FLOAT Rate_Coll_HeII;
108 CL_RT_FLOAT Rate_Radr_HII;
109 CL_RT_FLOAT Rate_Radr_HeII;
110 CL_RT_FLOAT Rate_Radr_HeIII;
111 CL_RT_FLOAT Rate_Diel_HeII;
112 CL_RT_FLOAT Rate_Chtr_HeII;
113
114 CL_RT_FLOAT Cool_Brem_1;
115 CL_RT_FLOAT Cool_Brem_2;
116 CL_RT_FLOAT Cool_Radr_HII;
117 CL_RT_FLOAT Cool_Radr_HeII;
118 CL_RT_FLOAT Cool_Radr_HeIII;
119 CL_RT_FLOAT Cool_Line_HI;
120 CL_RT_FLOAT Cool_Line_HeI;
121 CL_RT_FLOAT Cool_Line_HeII;
122 CL_RT_FLOAT Cool_LowT;
123} RATES_T;
124
125typedef struct clDerivsDataStruct clDerivsData;
126
127/* Heating Cooling Context */
128
129typedef struct CoolingPKDStruct {
130 double z; /* Redshift */
131 double dTime;
132 /* Rates independent of Temperature */
133 RATES_NO_T R;
134 /* Table for Temperature dependent rates */
135 int nTable;
136 double TMin;
137 double TMax;
138 double TlnMin;
139 double TlnMax;
140 double rDeltaTln;
141 RATES_T *RT;
142
143 int bMetal;
144 int nzMetalTable;
145 int nnHMetalTable;
146 int nTMetalTable;
147 double MetalTMin;
148 double MetalTMax;
149 double MetalTlogMin;
150 double MetalTlogMax;
151 double rDeltaTlog;
152 double MetalnHMin;
153 double MetalnHMax;
154 double MetalnHlogMin;
155 double MetalnHlogMax;
156 double rDeltanHlog;
157 double MetalzMin;
158 double MetalzMax;
159 double rDeltaz;
160 float ***MetalCoolln;
161 float ***MetalHeatln;
162
163 int nTableRead; /* number of Tables read from files */
164
165 int bUV;
166 int nUV;
167 UVSPECTRUM *UV;
168 int bUVTableUsesTime;
169 int bUVTableLinear;
170 int bLowTCool;
171 int bSelfShield;
172
173 double dGmPerCcUnit;
174 double dComovingGmPerCcUnit;
175 double dErgPerGmUnit;
176 double dSecUnit;
177 double dErgPerGmPerSecUnit;
178 double diErgPerGmUnit;
179 double dKpcUnit;
180 double dMassFracHelium;
181
182/* Diagnostic */
183 int its;
184} COOL;
185
186typedef struct {
187 double T, Tln;
188 double Coll_HI;
189 double Coll_HeI;
190 double Coll_HeII;
191 double Radr_HII;
192 double Radr_HeII;
193 double Diel_HeII;
194 double Chtr_HeII;
195 double Totr_HeII;
196 double Radr_HeIII;
197 double Cool_Metal;
198 double Heat_Metal;
199
200 double Phot_HI;
201 double Phot_HeI;
202 double Phot_HeII;
203} RATE;
204
205typedef struct {
206 double compton;
207 double bremHII;
208 double bremHeII;
209 double bremHeIII;
210 double radrecHII;
211 double radrecHeII;
212 double radrecHeIII;
213 double collionHI;
214 double collionHeI;
215 double collionHeII;
216 double dielrecHeII;
217 double lineHI;
218 double lineHeI;
219 double lineHeII;
220 double lowT;
221 double NetMetalCool;
222} COOL_ERGPERSPERGM;
223
224
225struct clDerivsDataStruct {
226 STIFF *IntegratorContext;
227 COOL *cl;
228 double rho,ExternalHeating,E,ZMetal;
229/* double Y_H, Y_He; */ /* will be needed -- also for temperature , Y_MetalIon, Y_eMetal */
230 RATE Rate;
231 PERBARYON Y;
232 double Y_H, Y_He, Y_eMax;
233 double Y_Total0, Y_Total1;
234 int its; /* Debug */
235 int bCool;
236};
237
238COOL *CoolInit( );
239void CoolFinalize( COOL *cl );
240clDerivsData *CoolDerivsInit(COOL *cl);
241void CoolDerivsFinalize(clDerivsData *cld ) ;
242
243void clInitConstants( COOL *cl, double dGMPerCcunit, double dComovingGmPerCcUnit,
244 double dErgPerGmUnit, double dSecUnit, double dKpcUnit, COOLPARAM CoolParam);
245void clInitUV(COOL *cl, int nTableColumns, int nTableRows, double *dTableData );
246void clInitRatesTable( COOL *cl, double TMin, double TMax, int nTable );
247void clReadMetalTable(COOL *cl, COOLPARAM clParam);
248void clRateMetalTable(COOL *cl, RATE *Rate, double T, double rho, double Y_H, double ZMetal);
249void clHHeTotal(COOL *cl, double ZMetal);
250void CoolInitRatesTable( COOL *cl, COOLPARAM CoolParam);
251
252void clRatesTableError( COOL *cl );
253void clRatesRedshift( COOL *cl, double z, double dTime );
254double clHeatTotal ( COOL *cl, PERBARYON *Y, RATE *Rate );
255void clRates( COOL *cl, RATE *Rate, double T, double rho);
256double clCoolTotal( COOL *cl, PERBARYON *Y, RATE *Rate, double rho, double ZMetal );
257COOL_ERGPERSPERGM clTestCool ( COOL *cl, PERBARYON *Y, RATE *Rate, double rho );
258void clPrintCool( COOL *cl, PERBARYON *Y, RATE *Rate, double rho );
259void clPrintCoolFile( COOL *cl, PERBARYON *Y, RATE *Rate, double rho, FILE *fp );
260
261void clAbunds( COOL *cl, PERBARYON *Y, RATE *Rate, double rho, double ZMetal);
262double clThermalEnergy( double Y_Total, double T );
263double clTemperature( double Y_Total, double E );
264double clRateCollHI( double T );
265double clRateCollHeI( double T );
266double clRateCollHeII( double T );
267double clRateRadrHII( double T );
268double clRateRadrHeII( double T );
269double clRateDielHeII( double T );
270double clRateChtrHeII(double T);
271double clRateRadrHeIII( double T );
272double clCoolBrem1( double T );
273double clCoolBrem2( double T );
274double clCoolRadrHII( double T );
275double clCoolRadrHeII( double T );
276double clCoolRadrHeIII( double T );
277double clCoolLineHI( double T );
278double clCoolLineHeI( double T );
279double clCoolLineHeII( double T );
280double clCoolLowT( double T );
281double clEdotInstant ( COOL *cl, PERBARYON *Y, RATE *Rate, double rho,
282 double ZMetal, double *dEdotHeat, double *EdotCool );
283 void clIntegrateEnergy(COOL *cl, clDerivsData *clData, PERBARYON *Y, double *E,
284 double ExternalHeating, double rho, double ZMetal, double dt );
285void clIntegrateEnergyDEBUG(COOL *cl, PERBARYON *Y, double *E,
286 double ExternalHeating, double rho, double ZMetal, double dt );
287
288
289void clDerivs(double x, const double *y, double *yheat,
290 double *ycool, void *Data) ;
291
292void CoolAddParams( COOLPARAM *CoolParam, PRM );
293void CoolLogParams( COOLPARAM *CoolParam, FILE *fp );
294void CoolOutputArray( COOLPARAM *CoolParam, int, int *, char * );
295
296#define COOL_ARRAY0_EXT "HI"
297double COOL_ARRAY0(COOL *cl, COOLPARTICLE *cp, double ZMetal);
298double COOL_SET_ARRAY0(COOL *cl, COOLPARTICLE *cp,double ZMetal, double val);
299
300#define COOL_ARRAY1_EXT "HeI"
301double COOL_ARRAY1(COOL *cl, COOLPARTICLE *cp, double ZMetal);
302double COOL_SET_ARRAY1(COOL *cl, COOLPARTICLE *cp,double ZMetal, double val);
303
304#define COOL_ARRAY2_EXT "HeII"
305double COOL_ARRAY2(COOL *cl, COOLPARTICLE *cp, double ZMetal);
306double COOL_SET_ARRAY2(COOL *cl, COOLPARTICLE *cp,double ZMetal, double val);
307
308#define COOL_ARRAY3_EXT "H2"
309#define COOL_ARRAY3(cl_, cp, aa ) (0)
310double COOL_SET_ARRAY3(COOL *cl_, COOLPARTICLE *cp,double aa, double bb_val);
311#define COOL_SET_ARRAY3( cl_, cp, aa, bb_val ) (0)
312
313double COOL_EDOT( COOL *cl_, COOLPARTICLE *cp_, double ECode_, double rhoCode_, double ZMetal_, double *posCode_ );
314#define COOL_EDOT( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_) (CoolCodeWorkToErgPerGmPerSec( cl_, CoolEdotInstantCode( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_ )))
315
316double COOL_COOLING( COOL *cl_, COOLPARTICLE *cp_, double ECode_, double rhoCode_, double ZMetal_, double *posCode_ );
317#define COOL_COOLING( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_) (CoolCodeWorkToErgPerGmPerSec( cl_, CoolCoolingCode( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_ )))
318
319double COOL_HEATING( COOL *cl_, COOLPARTICLE *cp_, double ECode_, double rhoCode_, double ZMetal_, double *posCode_ );
320#define COOL_HEATING( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_) (CoolCodeWorkToErgPerGmPerSec( cl_, CoolHeatingCode( cl_, cp_, ECode_, rhoCode_, ZMetal_, posCode_ )))
321
322void clSetAbundanceTotals(COOL *cl, double ZMetal, double *Y_H, double *Y_He, double *Y_eMAX);
323void CoolPARTICLEtoPERBARYON(COOL *cl_, PERBARYON *Y, COOLPARTICLE *cp, double ZMetal);
324void CoolPERBARYONtoPARTICLE(COOL *cl_, PERBARYON *Y, COOLPARTICLE *cp, double ZMetal);
325
326
327double CoolEnergyToTemperature( COOL *Cool, COOLPARTICLE *cp, double E, double ZMetal);
328double CoolCodeEnergyToTemperature( COOL *Cool, COOLPARTICLE *cp, double E, double ZMetal);
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 ZMetal);
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