00001 #include "../copyright.h" 00002 /*============================================================================*/ 00003 /*! \file cool.c 00004 * \brief Implements various optically thin cooling functions. 00005 * 00006 * These can be 00007 * enrolled by setting CoolingFunc=NAME in the problem generator, where NAME 00008 * is one of the functions in this file. 00009 * 00010 * Each cooling function returns the cooling rate per volume. The total 00011 * (or equivalently the internal) energy then evolves as 00012 * - dE/dt = de/dt = - CoolingFunc 00013 * 00014 * Some of these cooling functions return the cooling rate per volume in 00015 * cgs units [ergs/cm^{3}/s]. Thus, to use these functions, the entire 00016 * calculation must be in cgs, or else the cooling rate has to scaled 00017 * appropriately in the calling function. 00018 * 00019 * To add a new cooling function, implement it below and add the name to 00020 * src/microphysics/prototypes.h. Note the argument list must be (d,P,dt). 00021 * 00022 * CONTAINS PUBLIC FUNCTIONS: 00023 * - KoyInut() - Koyama & Inutsuka cooling function */ 00024 /*============================================================================*/ 00025 00026 #include <math.h> 00027 #include <float.h> 00028 #include "../defs.h" 00029 #include "../athena.h" 00030 #include "../globals.h" 00031 #include "prototypes.h" 00032 #include "../prototypes.h" 00033 00034 /* These constants are in cgs */ 00035 static const Real mbar = (1.37)*(1.6733e-24); 00036 static const Real kb = 1.380658e-16; 00037 static const Real HeatRate = 2.0e-26; 00038 static const Real Tmin = 10; 00039 00040 /*=========================== PUBLIC FUNCTIONS ===============================*/ 00041 /*----------------------------------------------------------------------------*/ 00042 /*! \fn Real KoyInut(const Real dens, const Real Press, const Real dt) 00043 * \brief Analytic fit to cooling in the diffuse ISM given by eq. (4) in 00044 * Koyama & Inutsuka, ApJ 564, L97 (2002); Returns rate in cgs. 00045 */ 00046 00047 #ifndef BAROTROPIC 00048 Real KoyInut(const Real dens, const Real Press, const Real dt) 00049 { 00050 Real n,coolrate=0.0; 00051 Real T,coolratepp,MaxdT,dT; 00052 Real Teq, logn, lognT; 00053 00054 /* Compute number density and Temperature */ 00055 n = dens/mbar; 00056 logn = log10(n); 00057 T = MAX((Press/(n*kb)),Tmin); 00058 00059 /* Compute the minimun Temperature*/ 00060 Teq = Tmin; 00061 00062 /* KI cooling rate per particle */ 00063 coolratepp = HeatRate* 00064 (n*(1.0e7*exp(-1.184e5/(T+1000.)) + 0.014*sqrt(T)*exp(-92.0/T)) - 1.0); 00065 00066 /* Expected dT by KI cooling rate */ 00067 dT = coolratepp*dt*Gamma_1/kb; 00068 00069 if ((T-dT) <= 185.0){ 00070 lognT = 3.9247499 - 1.8479378*logn + 1.5335032*logn*logn 00071 -0.47665872*pow(logn,3) + 0.076789136*pow(logn,4)-0.0049052587*pow(logn,5); 00072 Teq = pow(10.0,lognT) / n; 00073 } 00074 00075 /* Compute maximum change in T allowed to keep T positive, and limit cooling 00076 * rate to this value */ 00077 MaxdT = kb*(T-Teq)/(dt*Gamma_1); 00078 coolrate = MIN(coolratepp,MaxdT); 00079 return n*coolrate; 00080 } 00081 #endif /* BAROTROPIC */
1.7.1