cloudy trunk
|
00001 /* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and 00002 * others. For conditions of distribution and use see copyright notice in license.txt */ 00003 /*PunchSpecial generate output for the punch special command */ 00004 #include "cddefines.h" 00005 #include "taulines.h" 00006 #include "radius.h" 00007 #include "phycon.h" 00008 //#include "h2.h" 00009 #include "punch.h" 00010 00011 /*PunchSpecial generate output for the punch special command */ 00012 void PunchSpecial(FILE* ioPUN , 00013 const char *chTime) 00014 { 00015 /*long int i;*/ 00016 00017 DEBUG_ENTRY( "PunchSpecial()" ); 00018 00019 if( strncmp(chTime,"LAST",4) == 0 ) 00020 { 00021 /* code to execute only after last zone */ 00022 # if 0 00023 long ipISO , nelem , limit , i; 00024 double EdenAbund , fach; 00025 # include "physconst.h" 00026 # include "hydrogenic.h" 00027 PunFeII( ioPUN );*/ 00028 ipISO = ipHYDROGEN; 00029 nelem = ipHYDROGEN; 00030 00031 /* in all following the factor of two is because a single 00032 * decay produces two photons */ 00033 EdenAbund = StatesElem[ipH_LIKE][nelem][ipH2s].Pop*dense.xIonDense[nelem][nelem+1]*8.226*pow(1.+nelem,6); 00034 fprintf(ioPUN," 2s = %.3e\n", EdenAbund); 00035 00036 /* upper limit to H-like 2-phot is energy of La, which is in ipCont-1 cell */ 00037 limit = Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].ipCont-1; 00038 /* remember sum of rates, this will add up to twice the real rate since 00039 * each transition makes two photons */ 00040 for( i=0; i < limit; i++ ) 00041 { 00042 /*>>chng 01 jan 23, previous change had doubled cross section for H two-photon, 00043 * so here we divide by 2 to get old answer */ 00045 fach = iso.As2nu[ipISO][nelem][i]/2.f; 00046 fach *= rfield.anu2[i]/rfield.widflx[i]*EN1RYD; 00047 fprintf(ioPUN,"%.3e\t%.3e\t%.3e\n", 00048 RYDLAM/1e4/rfield.anu[i] , fach , fach*(realnum)EdenAbund ); 00049 } 00050 # endif 00051 00052 } 00053 else 00054 { 00055 fprintf(ioPUN,"%.5e\t%.3e\t%.2e\n", 00056 radius.depth , phycon.te , C12O16Rotate[0].Emis->TauIn ); 00057 # if 0 00058 long int iElecHi=1 , iVibHi=0, iRotHi=0; 00059 /* code to execute after every zone */ 00060 if( h2.lgH2ON ) 00061 { 00062 ASSERT( H2Lines[iElecHi][iVibHi][iRotHi][0][0][1].ipCont > 0 ); 00063 00064 fprintf(ioPUN,"DEBUG Oion\t%li\t%.2f\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\n", 00065 iteration,fnzone, 00066 radius.depth, 00067 H2Lines[iElecHi][iVibHi][iRotHi][0][0][1].Lo->Pop, 00068 H2Lines[iElecHi][iVibHi][iRotHi][0][0][1].Emis->PopOpc, 00069 H2Lines[iElecHi][iVibHi][iRotHi][0][0][1].Emis->opacity, 00070 H2Lines[iElecHi][iVibHi][iRotHi][0][0][1].Emis->TauCon, 00071 H2Lines[iElecHi][iVibHi][iRotHi][0][0][1].Emis->TauIn, 00072 H2Lines[iElecHi][iVibHi][iRotHi][0][0][1].Emis->pump, 00073 radius.drad 00074 ); 00075 } 00076 # endif 00077 /*DumpLine(&Transitions[ipHE_LIKE][ipHELIUM][ipHe2p1P][0] );*/ 00078 00079 } 00080 return; 00081 }