47 ( (1.-
Xdust[ipH2])/2.) );
61 double G1[
H2_TOP] = { 0.3 , 0.4 , 0.9 };
62 double G2[
H2_TOP] = { 0.6 , 0.6 , 0.4 };
65 if( (
energy_wn[0][iVib][0]+energy_off) <= Evm )
73 Fv =
sexp(
POW2( (
energy_wn[0][iVib][0]+energy_off - Evm)/(G2[ipH2]*(EH2 - Evm ) ) ) );
83 long int i , iElecHi , iElecLo;
84 long int iVibHi , iVibLo;
85 long int iRotHi , iRotLo;
86 long int iElec, iVib , iRot;
93 realnum sum , sumj , sumv , sumo , sump;
105 fprintf(
ioQQQ,
" H2_Create called in DEBUG mode.\n");
125 for( iVibHi = 0; iVibHi <=
h2.
nVib_hi[iElecHi]; ++iVibHi )
133 for( iVibLo=0; iVibLo<(
h2.
nVib_hi[iElecHi]+1); ++iVibLo )
137 for( iRotLo=0; iRotLo<=
h2.
nRot_hi[iElecHi][iVibLo]; ++iRotLo )
144 CollRateFit.reserve(iVibHi,iRotHi,iVibLo,iRotLo,nColl,3);
167 fprintf(
ioQQQ,
"elec %li highest vib= %li\n", iElec ,
h2.
nVib_hi[iElec] );
181 for( iVib = 0; iVib <=
h2.
nVib_hi[iElec]; ++iVib )
211 for( iVib = 0; iVib <=
h2.
nVib_hi[0]; ++iVib )
230 for( iVib = 0; iVib <=
h2.
nVib_hi[0]; ++iVib )
256 for( iVib = 0; iVib <=
h2.
nVib_hi[0]; ++iVib )
272 for( ipH2=0; ipH2<(int)
H2_TOP; ++ipH2 )
277 for( iVib = 0; iVib <=
h2.
nVib_hi[0]; ++iVib )
300 for( iVib = 0; iVib <=
h2.
nVib_hi[0]; ++iVib )
338 " H2_Create: there are %li electronic levels, in each level there are",
355 for( iVibHi = 0; iVibHi <=
VIB_COLLID; ++iVibHi )
359 for( iVibLo=0; iVibLo<(
VIB_COLLID+1); ++iVibLo )
361 for( iRotLo=0; iRotLo<=
h2.
nRot_hi[iElecHi][iVibLo]; ++iRotLo )
368 for( nColl=0; nColl<N_X_COLLIDER-2; ++nColl )
373 if(
CollRateFit[iVibHi][iRotHi][iVibLo][iRotLo][nColl][0] != 0. )
379 CollRateFit[iVibHi][iRotHi][iVibLo][iRotLo][nColl][0] += r;
384 if(
CollRateFit[iVibHi][iRotHi][iVibLo][iRotLo][N_X_COLLIDER-2][0] != 0. )
390 CollRateFit[iVibHi][iRotHi][iVibLo][iRotLo][N_X_COLLIDER-2][0] *= pow((
realnum)10.f,r);
412 for( i=1; i<nLevels_per_elec[0]; ++i )
423 for( iVibHi = 0; iVibHi <=
h2.
nVib_hi[iElecHi]; ++iVibHi )
434 for( iVibHi=0; iVibHi<=
h2.
nVib_hi[iElecHi]; ++iVibHi )
448 ASSERT( nH2_energies == nEner );
468 if( nEner+1 < nLevels_per_elec[0] )
485 for( nEner=0; nEner<nLevels_per_elec[0]; ++nEner )
509 " The total number of levels used in the matrix solver was set to %li but there are only %li levels in X.\n Sorry.\n",
511 nLevels_per_elec[0]);
522 for( iVibHi=0; iVibHi<=
h2.
nVib_hi[iElecHi]; ++iVibHi )
531 long int lim_elec_lo = 0;
532 H2Lines.reserve(iElecHi,iVibHi,iRotHi,1);
533 for( iElecLo=0; iElecLo<=lim_elec_lo; ++iElecLo )
540 if( iElecLo==iElecHi )
542 H2Lines.reserve(iElecHi,iVibHi,iRotHi,iElecLo,nv+1);
543 for( iVibLo=0; iVibLo<=nv; ++iVibLo )
546 if( iElecLo==iElecHi && iVibHi==iVibLo )
548 nr =
MAX2(1,iRotHi-1);
549 H2Lines.reserve(iElecHi,iVibHi,iRotHi,iElecLo,iVibLo,nr+1);
566 fprintf(
ioQQQ,
" There are a total of %li lines in the entire H2 molecule.\n", nlines );
571 for( iVibHi=0; iVibHi<=
h2.
nVib_hi[iElecHi]; ++iVibHi )
578 long int lim_elec_lo = 0;
579 for( iElecLo=0; iElecLo<=lim_elec_lo; ++iElecLo )
584 if( iElecLo==iElecHi )
586 for( iVibLo=0; iVibLo<=nv; ++iVibLo )
589 if( iElecLo==iElecHi && iVibHi==iVibLo )
591 for( iRotLo=
h2.
Jlowest[iElecLo]; iRotLo<=nr; ++iRotLo )
604 for( iVibHi=0; iVibHi<=
h2.
nVib_hi[iElecHi]; ++iVibHi )
613 long int lim_elec_lo = 0;
614 for( iElecLo=0; iElecLo<=lim_elec_lo; ++iElecLo )
619 if( iElecLo==iElecHi )
621 for( iVibLo=0; iVibLo<=nv; ++iVibLo )
624 if( iElecLo==iElecHi && iVibHi==iVibLo )
626 for( iRotLo=
h2.
Jlowest[iElecLo]; iRotLo<=nr; ++iRotLo )
628 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Hi =
629 H2Lines[iElecHi][iVibHi][iRotHi][0][0][0].Hi;
631 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Lo =
632 H2Lines[iElecLo][iVibLo][iRotLo][0][0][0].Hi;
654 for( iVibHi=0; iVibHi<=
h2.
nVib_hi[iElecHi]; ++iVibHi )
665 H2_stat[iElecHi][iVibHi][iRotHi] = 3.f*(2.f*iRotHi+1.f);
671 H2_stat[iElecHi][iVibHi][iRotHi] = (2.f*iRotHi+1.f);
680 for( iVibHi=0; iVibHi<=
h2.
nVib_hi[iElecHi]; ++iVibHi )
688 long int lim_elec_lo = 0;
689 for( iElecLo=0; iElecLo<=lim_elec_lo; ++iElecLo )
694 if( iElecLo==iElecHi )
696 for( iVibLo=0; iVibLo<=nv; ++iVibLo )
700 if( iElecLo==iElecHi && iVibHi==iVibLo )
702 for( iRotLo=
h2.
Jlowest[iElecLo]; iRotLo<=nr; ++iRotLo )
708 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Hi->nelem =
LIMELM+3;
710 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Hi->IonStg = 1;
713 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Lo->g =
H2_stat[iElecLo][iVibLo][iRotLo];
714 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Hi->g =
H2_stat[iElecHi][iVibHi][iRotHi];
717 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].EnergyWN =
721 if(
H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].EnergyWN >
SMALLFLOAT)
722 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].WLAng =
723 (
realnum)(1.e8f/
H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].EnergyWN /
724 RefIndex(
H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].EnergyWN));
726 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].EnergyK =
730 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].EnergyErg =
739 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Emis->iRedisFun =
ipCRDW;
742 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Emis->TauIn =
opac.
taumin;
743 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Emis->TauCon =
opac.
taumin;
745 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Emis->TauTot = 1e20f;
748 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Emis->dampXvel =
749 (
realnum)(
H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Emis->Aul/
750 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].EnergyWN/
PI4);
752 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Emis->gf =
754 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].EnergyWN,
755 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Hi->g ) );
758 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Emis->opacity = (
realnum)(
759 abscf(
H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Emis->gf,
760 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].EnergyWN,
761 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Lo->g));
778 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Coll.cs = (
realnum)(
779 pow(
H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].WLAng*1e-8,3)*
780 (
H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Hi->g/
781 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Lo->g)*
782 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Emis->Aul*
784 ASSERT(
H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Coll.cs>0.);
793 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Coll.cs = 0.;
800 H2Lines[iElecHi][iVibHi][iRotHi][iElecLo][iVibLo][iRotLo].Coll.cs = 0.;
816 for( ipH2=0; ipH2<(int)
H2_TOP; ++ipH2 )
830 # define T_H2_FORM 50000.
831 for( iVib = 0; iVib <=
h2.
nVib_hi[0]; ++iVib )
838 (1.f+2.f*
H2_lgOrtho[iElec][iVib][iRot]) * (1.f+iVib) *
842 sumv += iVib * H2_X_grain_formation_distribution[ipH2][iVib][iRot];
845 sumo += H2_X_grain_formation_distribution[ipH2][iVib][iRot];
850 sump += H2_X_grain_formation_distribution[ipH2][iVib][iRot];
858 double Xrot[
H2_TOP] = { 0.14 , 0.15 , 0.15 };
859 double Xtrans[
H2_TOP] = { 0.12 , 0.15 , 0.25 };
864 for( iVib = 0; iVib <=
h2.
nVib_hi[0]; ++iVib )
873 for( iVib = 0; iVib <=
h2.
nVib_hi[0]; ++iVib )
885 Erot = (EH2 - Ev) * Xrot[ipH2] / (Xrot[ipH2] + Xtrans[ipH2]);
923 double thermal_dist =
928 double aver = ( gaussian + thermal_dist ) / 2.;
939 (1.f+2.f*
H2_lgOrtho[iElec][iVib][iRot]) * Fv * (2.*iRot+1.) *
aver );
943 sumv += iVib * H2_X_grain_formation_distribution[ipH2][iVib][iRot];
946 sumo += H2_X_grain_formation_distribution[ipH2][iVib][iRot];
950 sump += H2_X_grain_formation_distribution[ipH2][iVib][iRot];
973 # define T_H2_FORM 17329.
974 for( iVib = 0; iVib <=
h2.
nVib_hi[0]; ++iVib )
985 sumv += iVib * H2_X_grain_formation_distribution[ipH2][iVib][iRot];
988 sumo += H2_X_grain_formation_distribution[ipH2][iVib][iRot];
993 sump += H2_X_grain_formation_distribution[ipH2][iVib][iRot];
1002 fprintf(
ioQQQ,
"H2 form grains mean J= %.3f mean v = %.3f ortho/para= %.3f\n",
1003 sumj/sum , sumv/sum , sumo/sump );
1007 for( iVib = 0; iVib <=
h2.
nVib_hi[0]; ++iVib )
1036 for( iVib=0; iVib<=
h2.
nVib_hi[iElec]; ++iVib )
1038 for( iRot=0; iRot<=
h2.
nRot_hi[iElec][iVib]; ++iRot )
1040 fprintf(
ioQQQ,
"%li\t%li\t%li\t%.5e\n",