00001 /* This file is part of Cloudy and is copyright (C)1978-2007 by Gary J. Ferland 00002 * For conditions of distribution and use see copyright notice in license.txt */ 00013 /*#include "physconst.h"*/ 00014 00015 /********************************************************************* 00016 * first come math constants * 00017 *********************************************************************/ 00018 00020 const double E = 2.718281828459045235360287; 00021 00023 const double EULER = 0.577215664901532860606512090082; 00024 00026 const double PI = 3.141592653589793238462643; 00027 00029 const double PI2 = 6.283185307179586476925287; 00030 00032 const double PI4 = 12.56637061435917295385057; 00033 00035 const double PI8 = 25.13274122871834590770115; 00036 00038 const double SQRT2 = 1.414213562373095048801689; 00039 00041 const double SQRTPI = 1.772453850905516027298167; 00042 00044 const double SQRTPIBY2 = 1.253314137315500251207883; 00045 00047 const double LN_TWO = 0.6931471805599453094172321; 00048 00050 const double LN_TEN = 2.302585092994045684017991; 00051 00053 const double LOG10_E = 0.4342944819032518276511289; 00054 00057 const double OPTDEP2EXTIN = 1.085736204758129569127822; 00058 00060 const double RADIAN = 57.29577951308232087679815; 00061 00062 /********************************************************************* 00063 * astronomical constants go here * 00064 *********************************************************************/ 00065 00068 const double SOLAR_MASS = 1.9889e33; 00069 00072 const double SOLAR_LUMINOSITY = 3.846e33; 00073 00076 /* >>refer phys const McCarthy, D.D., Petit, G., eds., IERS Technical Note 32, 00077 * >>refercon Frankfurt am Main: Verlag des Bundesamts fuer Kartographie und Geodaesie, 2004, 12 */ 00078 const double AU = 1.49597870691e13; 00079 00080 /********************************************************************* 00081 * fundamental constants go next, eventually rest should be defined * 00082 * in terms of these, these are Codata2002 values. * 00083 *********************************************************************/ 00084 00086 const double ATOMIC_MASS_UNIT = 1.66053886e-24; 00087 00089 const double ELECTRON_MASS = 9.1093826e-28; 00090 00092 const double PROTON_MASS = 1.67262171e-24; 00093 00095 const double BOLTZMANN = 1.3806505e-16; 00096 00098 const double SPEEDLIGHT = 2.99792458e10; 00099 00101 const double HPLANCK = 6.6260693e-27; 00102 00104 const double ELEM_CHARGE = 1.60217653e-19; 00105 00107 const double RYD_INF = 1.0973731568525e5; 00108 00111 const double HIONPOT = 0.999466508219; 00112 00113 /********************************************************************* 00114 * below here should be derived constants * 00115 * * 00116 * NB - explicit values in comments are approximate * 00117 * and are not maintained ! * 00118 *********************************************************************/ 00119 00121 const double PARSEC = AU*RADIAN*3600.; 00122 00124 const double H_BAR = HPLANCK/(2.*PI); 00125 00127 const double ELEM_CHARGE_ESU = ELEM_CHARGE*SPEEDLIGHT/10.; 00128 00130 const double ELECTRIC_CONST = 1.e11/(PI4*pow2(SPEEDLIGHT)); 00131 00137 const double HION_LTE_POP = pow2(HPLANCK)/(PI2*BOLTZMANN*ELECTRON_MASS); 00138 00141 const double SAHA = sqrt(pow3(HION_LTE_POP)); 00142 00144 const double ERG1CM = HPLANCK*SPEEDLIGHT; 00145 00147 const double T1CM = HPLANCK*SPEEDLIGHT/BOLTZMANN; 00148 00150 const double WAVNRYD = 1./RYD_INF; 00151 00153 const double RYDLAM = 1.e8/RYD_INF; 00154 00156 const double EN1RYD = HPLANCK*SPEEDLIGHT*RYD_INF; 00157 00160 const double TE1RYD = HPLANCK*SPEEDLIGHT*RYD_INF/BOLTZMANN; 00161 00163 const double EVDEGK = ELEM_CHARGE*1.e7/BOLTZMANN; 00164 00166 const double EVRYD = HPLANCK*SPEEDLIGHT*RYD_INF/ELEM_CHARGE*1.e-7; 00167 00169 const double EN1EV = EN1RYD/EVRYD; 00170 00172 const double FR1RYD = SPEEDLIGHT*RYD_INF; 00173 00175 const double HNU3C2 = 2.*HPLANCK*SPEEDLIGHT*pow3(RYD_INF); 00176 00178 const double FR1RYDHYD = SPEEDLIGHT*RYD_INF*HIONPOT; 00179 00181 const double HBAReV = H_BAR/EN1EV; 00182 00184 const double RYDLAMHYD = RYDLAM/HIONPOT; 00185 00187 const double STEFAN_BOLTZ = pow2(PI*pow2(BOLTZMANN))/(60.*pow3(H_BAR)*pow2(SPEEDLIGHT)); 00188 00190 const double FREQ_1EV = SPEEDLIGHT*RYD_INF/EVRYD; 00191 00193 const double FINE_STRUCTURE = pow2(ELEM_CHARGE_ESU)/SPEEDLIGHT/H_BAR; 00194 00196 const double FINE_STRUCTURE2 = pow2(FINE_STRUCTURE); 00197 00199 const double BOHR_RADIUS_CM = FINE_STRUCTURE/(PI4*RYD_INF); 00200 00202 const double TWO_PHOT_CONST = 9.*pow3(FINE_STRUCTURE2)*FR1RYD/2048.; 00203 00206 const double COLL_CONST = SAHA*BOLTZMANN/HPLANCK; 00207 00210 const double MILNE_CONST = SPEEDLIGHT*sqrt(pow3(FINE_STRUCTURE2)*pow3(TE1RYD)/PI); 00211 00214 const double TRANS_PROB_CONST = PI4*HPLANCK*FINE_STRUCTURE/ELECTRON_MASS;