51STATIC void newreact(
const char label[],
const char fun[],
double a,
double b,
double c);
83 template<
class T>
void newfunc()
89 ratefunc findfunc(
const char name[]);
92 class mole_reaction_hmrate;
94 class mole_reaction_th85rate;
96 class mole_reaction_crnurate_noalbedo;
100 class mole_reaction_grn_react;
102 class mole_reaction_h2_collh_excit;
104 class mole_reaction_h2_collh2_excit;
106 class mole_reaction_h2_collh_deexc;
108 class mole_reaction_h2_collh2_deexc;
110 class mole_reaction_rh2g_dis_h;
112 class mole_reaction_rh2s_dis_h;
114 class mole_reaction_rh2g_dis_h2;
116 class mole_reaction_rh2s_dis_h2;
118 class mole_reaction_rh2s_dis_h2_nodeexcit;
120 class mole_reaction_bh2g_dis_h;
122 class mole_reaction_bh2s_dis_h;
124 class mole_reaction_bh2g_dis_h2;
126 class mole_reaction_hneut;
128 class mole_reaction_cionhm;
145 typedef mole_reaction_null T;
147 virtual T* Create()
const {
return new T;}
148 virtual const char* name() {
return "null";}
186 for(n=0;n<nreact;n++)
208 te =
phycon.te+noneq_offset(rate);
213 return pow(te/300.,rate->
b)*exp(-rate->
c/te);
218 typedef mole_reaction_hmrate_exo T;
220 virtual T* Create()
const {
return new T;}
222 virtual const char* name() {
return "hmrate_exo";}
233 te =
MIN2( te, -10. * this->c );
235 return pow(te/300.,this->b)*exp(-te/this->c);
242 typedef mole_reaction_hmrate T;
244 virtual T* Create()
const {
return new T;}
245 virtual const char* name() {
return "hmrate";}
254 typedef mole_reaction_constrate T;
256 virtual T* Create()
const {
return new T;}
257 virtual const char* name() {
return "constrate";}
287 rk =
hmi.UV_Cont_rel2_Habing_TH85_depth/1.66;
291 rk =
hmi.UV_Cont_rel2_Habing_TH85_face/1.66*exp(-(rate->
c*
rfield.extin_mag_V_point));
298 typedef mole_reaction_th85rate T;
300 virtual T* Create()
const {
return new T;}
301 virtual const char* name() {
return "th85rate";}
304 return th85rate(
this);
330 typedef mole_reaction_crnurate_noalbedo T;
332 virtual T* Create()
const {
return new T;}
334 virtual const char* name() {
return "crnurate_noalbedo";}
337 return crnurate_noalbedo(
this);
343 typedef mole_reaction_crnurate T;
345 virtual T* Create()
const {
return new T;}
346 virtual const char* name() {
return "crnurate";}
349 return crnurate_noalbedo(
this)/(1.0-
albedo);
358 typedef mole_reaction_cryden_ov_bg T;
360 virtual T* Create()
const {
return new T;}
361 virtual const char* name() {
return "cryden_ov_bg";}
370 typedef mole_reaction_co_lnu_c_o_lnu T;
372 virtual T* Create()
const {
return new T;}
373 virtual const char* name() {
return "co_lnu_c_o_lnu";}
385 for( ns=0; ns<2; ++ns )
431 typedef mole_reaction_vib_evap T;
433 virtual T* Create()
const {
return new T;}
434 virtual const char* name() {
return "vib_evap";}
437 double binding_energy, exponent, vib_freq;
443 binding_energy = this->b;
444 double bin_total=0.0;
445 for(
size_t nd=0; nd <
gv.
bin.size() ; nd++ )
447 double bin_area =
gv.
bin[nd]->IntArea*
gv.
bin[nd]->cnv_H_pCM3;
448 exponent += exp(-binding_energy/
gv.
bin[nd]->tedust)*bin_area;
449 bin_total += bin_area;
451 exponent /= bin_total;
452 const double surface_density_of_sites = 1.5e15;
455 vib_freq = sqrt(2*surface_density_of_sites*(0.3*
BOLTZMANN)*binding_energy/(
PI*
PI*this->reactants[0]->mole_mass));
475 typedef mole_reaction_grn_abs T;
477 virtual T* Create()
const {
return new T;}
478 virtual const char* name() {
return "grn_abs";}
486 if (this->reactants[0]->n_nuclei() != 0)
487 mass = this->reactants[0]->mole_mass;
489 mass = this->reactants[1]->mole_mass;
497 typedef mole_reaction_grn_react T;
499 virtual T* Create()
const {
return new T;}
500 virtual const char* name() {
return "grn_react";}
503 return grn_react(
this);
527 double quant_barrier = 1e-8;
528 double surface_density_of_sites = 1.5e15;
533 double activ_barrier = rate->
c;
542 double dust_density = 0.0;
547 for(
size_t nd=0; nd <
gv.bin.size(); nd++ )
549 double bin_density =
gv.bin[nd]->IntArea*
gv.bin[nd]->cnv_H_pCM3;
550 Exp_i += exp(-E_i/
gv.bin[nd]->tedust)*bin_density;
551 Exp_j += exp(-E_j/
gv.bin[nd]->tedust)*bin_density;
552 dust_density += bin_density/(4*1e-10);
557 double total_sites = 4.*
mole.grain_area*surface_density_of_sites;
560 double diff_i = vib_freq_i*Exp_i/total_sites;
561 double diff_j = vib_freq_j*Exp_j/total_sites;
568 return scale*(diff_i + diff_j)/
SDIV(dust_density);
573 typedef mole_reaction_grn_photo T;
575 virtual T* Create()
const {
return new T;}
576 virtual const char* name() {
return "grn_photo";}
604 typedef mole_reaction_th85rate_co T;
606 virtual T* Create()
const {
return new T;}
608 virtual const char* name() {
return "th85rate_co";}
612 double esc_co, column;
632 if (this->reactants[0]->n_nuclei() != 0)
633 column =
mole.
species[ this->reactants[0]->index ].column;
635 column =
mole.
species[ this->reactants[1]->index ].column;
637 esc_co = 4.4e-15 * column /
642 return esca0k2(esc_co)*th85rate(
this);
648 typedef mole_reaction_oh_c2h2_co_ch3 T;
650 virtual T* Create()
const {
return new T;}
651 virtual const char* name() {
return "oh_c2h2_co_ch3";}
670 typedef mole_reaction_h_hnc_hcn_h T;
672 virtual T* Create()
const {
return new T;}
674 virtual const char* name() {
return "h_hnc_hcn_h";}
693 if(
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
695 return hmi.H2star_forms_hminus /
696 SDIV(
hmi.H2star_forms_hminus+
hmi.H2_forms_hminus);
710 return 1. - 4.938e-6;
715 double frac_H2star_grains(
void)
718 if(
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
720 return hmi.H2star_forms_grains /
721 SDIV(
hmi.H2star_forms_grains+
hmi.H2_forms_grains);
741 typedef mole_reaction_h2ph3p T;
743 virtual T* Create()
const {
return new T;}
744 virtual const char* name() {
return "h2ph3p";}
766 typedef mole_reaction_hopexch T;
768 virtual T* Create()
const {
return new T;}
769 virtual const char* name() {
return "hopexch";}
778 typedef mole_reaction_hpoexch T;
780 virtual T* Create()
const {
return new T;}
781 virtual const char* name() {
return "hpoexch";}
790 typedef mole_reaction_hmattach T;
792 virtual T* Create()
const {
return new T;}
793 virtual const char* name() {
return "hmattach";}
802 typedef mole_reaction_hmihph2p T;
804 virtual T* Create()
const {
return new T;}
805 virtual const char* name() {
return "hmihph2p";}
829 typedef mole_reaction_hmphoto T;
831 virtual T* Create()
const {
return new T;}
832 virtual const char* name() {
return "hmphoto";}
862 typedef mole_reaction_cionhm T;
864 virtual T* Create()
const {
return new T;}
865 virtual const char* name() {
return "cionhm";}
872 double assoc_detach(
void)
887 y=545969508.1323510+x*71239.23653059864;
893 typedef mole_reaction_c3bod T;
895 virtual T* Create()
const {
return new T;}
896 virtual const char* name() {
return "c3bod";}
905 typedef mole_reaction_asdfg T;
907 virtual T* Create()
const {
return new T;}
908 virtual const char* name() {
return "asdfg";}
917 typedef mole_reaction_asdfs T;
919 virtual T* Create()
const {
return new T;}
920 virtual const char* name() {
return "asdfs";}
929 typedef mole_reaction_asdbg T;
931 virtual T* Create()
const {
return new T;}
932 virtual const char* name() {
return "asdbg";}
942 typedef mole_reaction_asdbs T;
944 virtual T* Create()
const {
return new T;}
945 virtual const char* name() {
return "asdbs";}
955 typedef mole_reaction_bhneut T;
957 virtual T* Create()
const {
return new T;}
958 virtual const char* name() {
return "bhneut";}
969 return (hneut(
this)*ratio*
988 return 1.4e-7*pow(
phycon.te/300,-0.487)*exp(
phycon.te/29300);
992 return 3.4738192887404660e-008;
997 typedef mole_reaction_hneut T;
999 virtual T* Create()
const {
return new T;}
1000 virtual const char* name() {
return "hneut";}
1010 typedef mole_reaction_h2_spon_diss T;
1012 virtual T* Create()
const {
return new T;}
1013 virtual const char* name() {
return "h2_spon_diss";}
1034 virtual const char* name() {
return "grnh2tot";}
1037 return grnh2tot(
this);
1043 typedef mole_reaction_grnh2 T;
1045 virtual T* Create()
const {
return new T;}
1046 virtual const char* name() {
return "grnh2";}
1049 return grnh2tot(
this)*(1.-frac_H2star_grains());
1055 typedef mole_reaction_grnh2s T;
1057 virtual T* Create()
const {
return new T;}
1058 virtual const char* name() {
return "grnh2s";}
1061 return grnh2tot(
this)*frac_H2star_grains();
1067 typedef mole_reaction_radasc T;
1069 virtual T* Create()
const {
return new T;}
1070 virtual const char* name() {
return "radasc";}
1078 return hmrate(
this) *
1090 typedef mole_reaction_assoc_ion T;
1092 virtual T* Create()
const {
return new T;}
1093 virtual const char* name() {
return "assoc_ion";}
1100 return hmrate(
this) *
1113 if(
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
1115 return h2.Average_collH_dissoc_g;
1119 double corr =
MIN2(6.,14.44-
phycon.alogte*3.08);
1131 typedef mole_reaction_rh2g_dis_h T;
1133 virtual T* Create()
const {
return new T;}
1134 virtual const char* name() {
return "rh2g_dis_h";}
1137 return rh2g_dis_h(
this);
1143 if(
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
1145 return h2.Average_collH_dissoc_s;
1150 return HMRATE(4.67e-7,-1.,5.5e4);
1155 typedef mole_reaction_rh2s_dis_h T;
1157 virtual T* Create()
const {
return new T;}
1158 virtual const char* name() {
return "rh2s_dis_h";}
1161 return rh2s_dis_h(
this);
1167 if(
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
1169 return h2.Average_collH2_dissoc_g;
1175 return HMRATE(5.5e-29*0.5/(
SAHA*3.634e-5)*sqrt(300.),0.5,5.195e4);
1180 typedef mole_reaction_rh2g_dis_h2 T;
1182 virtual T* Create()
const {
return new T;}
1183 virtual const char* name() {
return "rh2g_dis_h2";}
1186 return rh2g_dis_h2(
this);
1192 if(
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
1194 return h2.Average_collH2_dissoc_s;
1199 return HMRATE(1e-11,0.,0.);
1204 typedef mole_reaction_rh2s_dis_h2 T;
1206 virtual T* Create()
const {
return new T;}
1207 virtual const char* name() {
return "rh2s_dis_h2";}
1210 return rh2s_dis_h2(
this);
1216 if(
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
1218 return h2.Average_collH2_dissoc_s;
1223 return HMRATE(1e-11,0.,2.18e4);
1226 class mole_reaction_rh2s_dis_h2_nodeexcit :
public mole_reaction
1228 typedef mole_reaction_rh2s_dis_h2_nodeexcit T;
1230 virtual T* Create()
const {
return new T;}
1231 virtual const char* name() {
return "rh2s_dis_h2_nodeexcit";}
1234 return rh2s_dis_h2_nodeexcit(
this);
1240 return rh2g_dis_h(rate)*
hmi.rel_pop_LTE_H2g;
1244 typedef mole_reaction_bh2g_dis_h T;
1246 virtual T* Create()
const {
return new T;}
1247 virtual const char* name() {
return "bh2g_dis_h";}
1250 return bh2g_dis_h(
this);
1256 return rh2s_dis_h(rate)*
hmi.rel_pop_LTE_H2s;
1260 typedef mole_reaction_bh2s_dis_h T;
1262 virtual T* Create()
const {
return new T;}
1263 virtual const char* name() {
return "bh2s_dis_h";}
1266 return bh2s_dis_h(
this);
1272 return rh2g_dis_h2(rate)*
hmi.rel_pop_LTE_H2g;
1276 typedef mole_reaction_bh2g_dis_h2 T;
1278 virtual T* Create()
const {
return new T;}
1279 virtual const char* name() {
return "bh2g_dis_h2";}
1282 return bh2g_dis_h2(
this);
1288 typedef mole_reaction_bh2s_dis_h2 T;
1290 virtual T* Create()
const {
return new T;}
1291 virtual const char* name() {
return "bh2s_dis_h2";}
1300 typedef mole_reaction_h2photon T;
1302 virtual T* Create()
const {
return new T;}
1303 virtual const char* name() {
return "h2photon";}
1312 typedef mole_reaction_h2crphot T;
1314 virtual T* Create()
const {
return new T;}
1315 virtual const char* name() {
return "h2crphot";}
1324 typedef mole_reaction_h2crphh T;
1326 virtual T* Create()
const {
return new T;}
1327 virtual const char* name() {
return "h2crphh";}
1343 typedef mole_reaction_h2scrphh T;
1345 virtual T* Create()
const {
return new T;}
1346 virtual const char* name() {
return "h2scrphh";}
1362 typedef mole_reaction_radath T;
1364 virtual T* Create()
const {
return new T;}
1365 virtual const char* name() {
return "radath";}
1374 typedef mole_reaction_gamtwo T;
1376 virtual T* Create()
const {
return new T;}
1377 virtual const char* name() {
return "gamtwo";}
1387 typedef mole_reaction_hlike_phot T;
1389 virtual T* Create()
const {
return new T;}
1390 virtual const char* name() {
return "hlike_phot";}
1406 typedef mole_reaction_h2s_sp_decay T;
1408 virtual T* Create()
const {
return new T;}
1409 virtual const char* name() {
return "h2s_sp_decay";}
1427 if(
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
1429 return h2.Average_collH2_deexcit;
1438 typedef mole_reaction_h2_collh2_deexc T;
1440 virtual T* Create()
const {
return new T;}
1441 virtual const char* name() {
return "h2_collh2_deexc";}
1444 return h2_collh2_deexc(
this);
1450 if (
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
1452 return h2.Average_collH_deexcit;
1461 typedef mole_reaction_h2_collh_deexc T;
1463 virtual T* Create()
const {
return new T;}
1464 virtual const char* name() {
return "h2_collh_deexc";}
1467 return h2_collh_deexc(
this);
1474 if (
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
1476 return h2.Average_collH2_excit;
1480 return h2_collh2_deexc(rate)*
sexp( 30172./
phycon.te);
1485 typedef mole_reaction_h2_collh2_excit T;
1487 virtual T* Create()
const {
return new T;}
1488 virtual const char* name() {
return "h2_collh2_excit";}
1491 return h2_collh2_excit(
this);
1498 if (
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
1500 return h2.Average_collH_excit;
1504 return h2_collh_deexc(rate)*
sexp( 30172./
phycon.te);
1509 typedef mole_reaction_h2_collh_excit T;
1511 virtual T* Create()
const {
return new T;}
1512 virtual const char* name() {
return "h2_collh_excit";}
1515 return h2_collh_excit(
this);
1521 typedef mole_reaction_h2gexcit T;
1523 virtual T* Create()
const {
return new T;}
1524 virtual const char* name() {
return "h2gexcit";}
1533 typedef mole_reaction_h2sdissoc T;
1535 virtual T* Create()
const {
return new T;}
1536 virtual const char* name() {
return "h2sdissoc";};
1551 typedef mole_reaction_h2gdissoc T;
1553 virtual T* Create()
const {
return new T;}
1554 virtual const char* name() {
return "h2gdissoc";}
1569 typedef mole_reaction_hd_photodissoc T;
1571 virtual T* Create()
const {
return new T;}
1572 virtual const char* name() {
return "hd_photodissoc";}
1598 else if( te < 1200. )
1602 else if( te < 3800. )
1607 else if( te <= 1.4e4 )
1630 typedef mole_reaction_gamheh T;
1632 virtual T* Create()
const {
return new T;}
1633 virtual const char* name() {
return "gamheh";}
1645 for( i=
hmi.
iheh1-1; i < limit; i++ )
1658 enum {exclude, base, umisthack, federman, lithium, deuterium, ti, misc};
1662 static bool lgReactInitialized =
false;
1672 if( lgReactInitialized )
1675 lgReactInitialized =
true;
1684 newfunc<mole_reaction_null>();
1685 newfunc<mole_reaction_hmrate>();
1686 newfunc<mole_reaction_hmrate_exo>();
1687 newfunc<mole_reaction_constrate>();
1688 newfunc<mole_reaction_th85rate>();
1689 newfunc<mole_reaction_crnurate>();
1690 newfunc<mole_reaction_crnurate_noalbedo>();
1691 newfunc<mole_reaction_co_lnu_c_o_lnu>();
1692 newfunc<mole_reaction_vib_evap>();
1693 newfunc<mole_reaction_th85rate_co>();
1694 newfunc<mole_reaction_grn_abs>();
1696 newfunc<mole_reaction_grn_react>();
1698 newfunc<mole_reaction_grn_photo>();
1699 newfunc<mole_reaction_oh_c2h2_co_ch3>();
1700 newfunc<mole_reaction_h_hnc_hcn_h>();
1702 newfunc<mole_reaction_gamheh>();
1703 newfunc<mole_reaction_hd_photodissoc>();
1704 newfunc<mole_reaction_h2gdissoc>();
1705 newfunc<mole_reaction_h2sdissoc>();
1706 newfunc<mole_reaction_h2gexcit>();
1707 newfunc<mole_reaction_h2_collh_excit>();
1708 newfunc<mole_reaction_h2_collh2_excit>();
1709 newfunc<mole_reaction_h2_collh_deexc>();
1710 newfunc<mole_reaction_h2_collh2_deexc>();
1711 newfunc<mole_reaction_h2s_sp_decay>();
1712 newfunc<mole_reaction_hlike_phot>();
1713 newfunc<mole_reaction_gamtwo>();
1714 newfunc<mole_reaction_h2ph3p>();
1715 newfunc<mole_reaction_radath>();
1716 newfunc<mole_reaction_cryden_ov_bg>();
1717 newfunc<mole_reaction_h2scrphh>();
1718 newfunc<mole_reaction_h2crphh>();
1719 newfunc<mole_reaction_h2photon>();
1720 newfunc<mole_reaction_h2crphot>();
1721 newfunc<mole_reaction_rh2g_dis_h>();
1722 newfunc<mole_reaction_rh2s_dis_h>();
1723 newfunc<mole_reaction_rh2g_dis_h2>();
1724 newfunc<mole_reaction_rh2s_dis_h2>();
1725 newfunc<mole_reaction_rh2s_dis_h2_nodeexcit>();
1726 newfunc<mole_reaction_bh2g_dis_h>();
1727 newfunc<mole_reaction_bh2s_dis_h>();
1728 newfunc<mole_reaction_bh2g_dis_h2>();
1729 newfunc<mole_reaction_bh2s_dis_h2>();
1730 newfunc<mole_reaction_radasc>();
1731 newfunc<mole_reaction_assoc_ion>();
1732 newfunc<mole_reaction_grnh2>();
1733 newfunc<mole_reaction_grnh2s>();
1734 newfunc<mole_reaction_h2_spon_diss>();
1735 newfunc<mole_reaction_bhneut>();
1736 newfunc<mole_reaction_hneut>();
1737 newfunc<mole_reaction_asdbs>();
1738 newfunc<mole_reaction_asdbg>();
1739 newfunc<mole_reaction_asdfs>();
1740 newfunc<mole_reaction_asdfg>();
1741 newfunc<mole_reaction_c3bod>();
1742 newfunc<mole_reaction_cionhm>();
1743 newfunc<mole_reaction_hmphoto>();
1744 newfunc<mole_reaction_hmihph2p>();
1745 newfunc<mole_reaction_hmattach>();
1780 newreact(
"C+,OH=>CO,H+",
"hmrate",0.,0.,0.);
1783 newreact(
"H,H,grn=>H2,grn",
"grnh2",1.,0.,0.);
1784 newreact(
"H,H,grn=>H2*,grn",
"grnh2s",1.,0.,0.);
1785 newreact(
"H-,PHOTON=>H,e-",
"hmphoto",1.,0.,0.);
1794 long nelem = (*atom)->el->Z-1;
1798 sprintf(react,
"H-,%s+=>H,%s", (*atom)->label().c_str(), (*atom)->label().c_str() );
1799 newreact(react,
"hmrate",4e-6/sqrt(300.),-0.5,0.);
1803 newreact(
"H,e-=>H-,PHOTON",
"hmattach",1.,0.,0.);
1804 newreact(
"H-,H+=>H2+,e-",
"hmihph2p",1.,0.,0.);
1805 newreact(
"H-,e-=>H,e-,e-",
"cionhm",1.,0.,0.);
1806 newreact(
"H,e-,e-=>H-,e-",
"c3bod",1.,0.,0.);
1807 newreact(
"H,H-=>H2,e-",
"asdfg",1.,0.,0.);
1808 newreact(
"H,H-=>H2*,e-",
"asdfs",1.,0.,0.);
1809 newreact(
"H2,e-=>H,H-",
"asdbg",1.,0.,0.);
1810 newreact(
"H2*,e-=>H,H-",
"asdbs",1.,0.,0.);
1811 newreact(
"H-,H+=>H,H",
"hneut",1.,0.,0.);
1812 newreact(
"H,H=>H-,H+",
"bhneut",1.,0.,0.);
1813 newreact(
"H2*=>H,H",
"h2_spon_diss",1.,0.,0.);
1816 newreact(
"H,H=>H2,PHOTON",
"radasc",3e-14,0.,0.);
1817 newreact(
"H,H=>H2+,e-",
"assoc_ion",3.27e-10,-0.35,17829.);
1818 newreact(
"H2,H=>H,H,H",
"rh2g_dis_h",1.,0.,0.);
1819 newreact(
"H2,H2=>H,H,H2",
"rh2g_dis_h2",1.,0.,0.);
1820 newreact(
"H,H,H=>H2,H",
"bh2g_dis_h",1.,0.,0.);
1821 newreact(
"H,H,H2=>H2,H2",
"bh2g_dis_h2",1.,0.,0.);
1823 newreact(
"H2,PHOTON=>H2+,e-",
"h2photon",1.,0.,0.);
1824 newreact(
"H2,CRPHOT=>H2+,e-",
"h2crphot",1.,0.,0.);
1825 newreact(
"H2*,PHOTON=>H2+,e-",
"h2photon",1.,0.,0.);
1826 newreact(
"H2*,CRPHOT=>H2+,e-",
"h2crphot",1.,0.,0.);
1827 newreact(
"H2,CRPHOT=>H,H",
"h2crphh",1.,0.,0.);
1828 newreact(
"H2,CRPHOT=>H+,H-",
"cryden_ov_bg",3.9e-21,0.,0.);
1829 newreact(
"H2*,CRPHOT=>H+,H-",
"cryden_ov_bg",3.9e-21,0.,0.);
1830 newreact(
"H2,CRPHOT=>H+,H,e-",
"cryden_ov_bg",2.2e-19,0.,0.);
1831 newreact(
"H2*,CRPHOT=>H+,H,e-",
"cryden_ov_bg",2.2e-19,0.,0.);
1832 newreact(
"H+,H=>H2+,PHOTON",
"radath",1.,0.,0.);
1833 newreact(
"H2+,PHOTON=>H,H+",
"gamtwo",1.,0.,0.);
1834 newreact(
"H2+,CRPHOT=>H,H+",
"hlike_phot",1.,0.,0.);
1837 newreact(
"H3+,CRPHOT=>H2+,H+,e-",
"hlike_phot",2.,0.,0.);
1838 newreact(
"H,H+,e-=>H2+,e-",
"hmrate",2e-7 *
SAHA * (4./(2.*1.)) * 3.634e-5 * pow(300.,-1.50),-1.50,0.);
1842 newreact(
"H2,CRPHOT=>H2+,e-",
"hmrate",4.4e-17,0.,0.);
1843 newreact(
"H2,CRPHOT=>H,H+,e-",
"hmrate",1e-19,0.,0.);
1844 newreact(
"H2*,CRPHOT=>H,H+,e-",
"hmrate",1e-19,0.,0.);
1845 newreact(
"H2*,CRPHOT=>H2+,e-",
"hmrate",4.4e-17,0.,0.);
1846 newreact(
"H2,CRPHOT=>H,H",
"hmrate",5e-18,0.,0.);
1847 newreact(
"e-,H3+=>H2,H",
"hmrate",2.5e-8,-0.3,0.);
1848 newreact(
"e-,H3+=>H,H,H",
"hmrate",7.5e-8,-0.3,0.);
1849 newreact(
"H+,H=>H2+,PHOTON",
"hmrate",5.3e-19,1.85,0);
1853 newreact(
"H2+,e-=>H,H",
"hmrate",1.6e-8,-0.43,0.);
1854 newreact(
"H2+,PHOTON=>H,H+",
"th85rate",5.7e-10,0.,1.9);
1857 newreact(
"H2*,CRPHOT=>H,H",
"h2scrphh",1.,0.,0.);
1858 newreact(
"H2,H2+=>H,H3+",
"h2ph3p",1.,0.,0.);
1859 newreact(
"H2*=>H2,PHOTON",
"h2s_sp_decay",1.,0.,0.);
1860 newreact(
"H2*,H2=>H2,H2",
"h2_collh2_deexc",1.,0.,0.);
1861 newreact(
"H2*,H=>H2,H",
"h2_collh_deexc",1.,0.,0.);
1862 newreact(
"H2,H2=>H2*,H2",
"h2_collh2_excit",1.,0.,0.);
1863 newreact(
"H2,H=>H2*,H",
"h2_collh_excit",1.,0.,0.);
1866 newreact(
"H2*,H=>H,H,H",
"rh2s_dis_h",1.,0.,0.);
1867 newreact(
"H2*,H2=>H2,H,H",
"rh2s_dis_h2",1.,0.,0.);
1868 newreact(
"H2*,H2*=>H2,H,H",
"rh2s_dis_h2",1.,0.,0.);
1869 newreact(
"H2*,H2*=>H2*,H,H",
"rh2s_dis_h2_nodeexcit",1.,0.,0.);
1873 newreact(
"H2,He=>He,H,H",
"rh2g_dis_h",1.,0.,0.);
1874 newreact(
"H2,H+=>H+,H,H",
"rh2g_dis_h",1.,0.,0.);
1875 newreact(
"H2,H3+=>H3+,H,H",
"rh2g_dis_h",1.,0.,0.);
1876 newreact(
"H2,e-=>e-,H,H",
"rh2g_dis_h",1.,0.,0.);
1877 newreact(
"H2*,He=>He,H,H",
"rh2s_dis_h",1.,0.,0.);
1878 newreact(
"H2*,H+=>H+,H,H",
"rh2s_dis_h",1.,0.,0.);
1879 newreact(
"H2*,H3+=>H3+,H,H",
"rh2s_dis_h",1.,0.,0.);
1880 newreact(
"H2*,e-=>e-,H,H",
"rh2s_dis_h",1.,0.,0.);
1883 newreact(
"He,H,H=>H2,He",
"bh2g_dis_h",1.,0.,0.);
1884 newreact(
"H+,H,H=>H2,H+",
"bh2g_dis_h",1.,0.,0.);
1885 newreact(
"H3+,H,H=>H2,H3+",
"bh2g_dis_h",1.,0.,0.);
1886 newreact(
"e-,H,H=>H2,e-",
"bh2g_dis_h",1.,0.,0.);
1887 newreact(
"H,H,H=>H2*,H",
"bh2s_dis_h",1.,0.,0.);
1888 newreact(
"H2,H,H=>H2*,H2",
"bh2s_dis_h2",1.,0.,0.);
1889 newreact(
"H2,H,H=>H2*,H2*",
"bh2s_dis_h2",1.,0.,0.);
1890 newreact(
"H2*,H,H=>H2*,H2*",
"bh2s_dis_h2",1.,0.,0.);
1891 newreact(
"He,H,H=>H2*,He",
"bh2s_dis_h",1.,0.,0.);
1892 newreact(
"H+,H,H=>H2*,H+",
"bh2s_dis_h",1.,0.,0.);
1893 newreact(
"H3+,H,H=>H2*,H3+",
"bh2s_dis_h",1.,0.,0.);
1894 newreact(
"e-,H,H=>H2*,e-",
"bh2s_dis_h",1.,0.,0.);
1897 newreact(
"H2,PHOTON=>H2*",
"h2gexcit",1.,0.,0.);
1898 newreact(
"H2*,PHOTON=>H,H",
"h2sdissoc",1.,0.,0.);
1899 newreact(
"H2,PHOTON=>H,H",
"h2gdissoc",1.,0.,0.);
1900 newreact(
"HeH+,PHOTON=>H,He+",
"gamheh",1.,0.,0.);
1905 newreact(
"OHgrn,Hgrn=>H2Ogrn",
"grn_react",1.,0.,0.);
1910 fprintf(stderr,
"Finished testing against UDFA database\n");
1924 it->second->index = index;
1926 mole.reaction_rks.resize( index );
1928 memset( &
mole.reaction_rks[0], 0, (
unsigned long)(index)*
sizeof(
double) );
1937 DEBUG_ENTRY(
"mole_generate_isotopologue_reactions()");
1939 bool lgDebug =
false;
1942 vector<string> newReactionStrings;
1943 vector<mole_reaction*> oldReactions;
1948 bool lgFound =
false;
1950 for(
long i=0; i<it->second->nreactants; ++i )
1953 for(
nAtoms_i atom=it->second->reactants[i]->nAtom.begin(); atom != it->second->reactants[i]->nAtom.end(); ++atom )
1955 if( atom->first->label() == atom_old )
1969 fprintf(
ioQQQ,
"DEBUGGG mole_generate_isotopologue_reactions %s ..........\n", it->first.c_str() );
1971 for(
long i=0; i<it->second->nreactants; ++i )
1977 vector<string> react_iso_labels;
1979 vector< int > numAtoms;
1980 string embellishments;
1982 bool lgParseOK =
parse_species_label( it->second->reactants[i]->label.c_str(), atomsLeftToRight, numAtoms, embellishments );
1983 ASSERT( lgParseOK ==
true );
1992 for(
unsigned j=0; j<react_iso_labels.size(); ++j )
1994 for(
long k=0; k<it->second->nproducts; ++k )
2000 vector<string> prod_iso_labels;
2001 atomsLeftToRight.clear();
2003 embellishments.clear();
2005 lgParseOK =
parse_species_label( it->second->products[k]->label.c_str(), atomsLeftToRight, numAtoms, embellishments );
2006 ASSERT( lgParseOK ==
true );
2015 for(
unsigned l=0; l<prod_iso_labels.size(); ++l )
2017 string react_string;
2020 for(
long i1=0; i1<i; ++i1 )
2022 react_string += it->second->reactants[i1]->label;
2023 if( i1 != it->second->nreactants-1 )
2024 react_string +=
",";
2026 react_string += react_iso_labels[j];
2027 if( i != it->second->nreactants-1 )
2028 react_string +=
",";
2029 for(
long i2=i+1; i2<it->second->nreactants; ++i2 )
2031 react_string += it->second->reactants[i2]->label;
2032 if( i2 != it->second->nreactants-1 )
2033 react_string +=
",";
2036 react_string +=
"=>";
2038 for(
long k1=0; k1<k; ++k1 )
2040 react_string += it->second->products[k1]->label;
2041 if( k1 != it->second->nproducts-1 )
2042 react_string +=
",";
2044 react_string += prod_iso_labels[l];
2045 if( k != it->second->nproducts-1 )
2046 react_string +=
",";
2047 for(
long k2=k+1; k2<it->second->nproducts; ++k2 )
2049 react_string += it->second->products[k2]->label;
2050 if( k2 != it->second->nproducts-1 )
2051 react_string +=
",";
2055 fprintf(
ioQQQ,
"DEBUGGG mole_generate_isotopologue_reactions .................... %s\n",
2056 react_string.c_str() );
2060 newReactionStrings.push_back( canon_react_string );
2061 oldReactions.push_back( it->second.get_ptr() );
2068 ASSERT( oldReactions.size() == newReactionStrings.size() );
2071 vector<mole_reaction*>::const_iterator it2 = oldReactions.begin();
2072 for( vector<string>::const_iterator it1 = newReactionStrings.begin(); it1 != newReactionStrings.end(); ++it1, ++it2 )
2078 newreact( it1->c_str(), (*it2)->name(), (*it2)->a, (*it2)->b, (*it2)->c );
2088 char chLabel[50], chLabelSave[50];
2095 strcpy( chLabel, p->second->label.c_str() );
2096 strcpy( chLabelSave, p->second->label.c_str() );
2097 char *str = chLabel;
2098 const char *delim =
"=>";
2099 char *chNewProducts = strtok( str, delim );
2100 char *chNewReactants = strtok( NULL, delim );
2101 char chNewReaction[50];
2103 strcpy( chNewReaction, chNewReactants );
2104 strcat( chNewReaction,
"=>" );
2105 strcat( chNewReaction, chNewProducts );
2113 if(
trace.lgTraceMole )
2115 fprintf(
ioQQQ,
"Warning! No reverse reaction for %30s. ", p->second->label.c_str() );
2116 fprintf(
ioQQQ,
"\n" );
2137 double ln_result = 0.;
2143 ln_result += log(fac);
2151 ln_result -= log(fac);
2153 double result = exp( ln_result );
2163 if( sp->
label ==
"PHOTON" || sp->
label ==
"CRPHOT" )
2169 else if( sp->
label ==
"CRP" || sp->
label ==
"grn" )
2187 vector<long> numForm, numDest;
2194 for(
long i=0; i<rate->nreactants; ++i)
2196 ++numDest[ rate->reactants[i]->index ];
2199 for(
long i=0; i<rate->nproducts; ++i)
2201 ++numForm[ rate->products[i]->index ];
2205 for(
unsigned i=0; i<numForm.size(); ++i )
2207 if( numForm[i]==0 && numDest[i]==0 )
2209 if( numForm[i]>1 && numDest[i]>1 )
2213 fprintf(
ioQQQ,
"DEBUGGG mole_cmp_num_in_out_reactions %*s: in %4li out %4li\n",
CHARS_SPECIES,
mole_global.list[i]->label.c_str(), numForm[i], numDest[i] );
2222 char *label, *reactstr, *f;
2237STATIC void newreact(
const char label[],
const char fun[],
double a,
double b,
double c)
2241 ratefunc rate = findfunc(fun);
2242 if (rate.get_ptr() == NULL)
2244 fprintf(stderr,
"Rate function %s not known for reaction %s. Aborting. Sorry.\n",fun,label);
2248 rate->label = label;
2262 const char *rateLabelPtr = rate->label.c_str();
2266 rate->udfastate =
ABSENT;
2270 if (rate->udfastate ==
ABSENT)
2272 fprintf(stderr,
"Reaction %s not in UDFA\n",rateLabelPtr);
2277 if (rate->nreactants == 2 && rate->reactants[0]->mole_mass!=0.0 && rate->reactants[1]->mole_mass!=0.0 )
2279 rate->reduced_mass = 1./(1./rate->reactants[0]->mole_mass+1./rate->reactants[1]->mole_mass);
2283 rate->reduced_mass = 0.;
2293 fprintf(
ioQQQ,
"Warning: duplicate reaction %s -- using new version\n",rateLabelPtr);
2304 rate->reactants[i] = NULL;
2306 rate->nreactants = 0;
2310 rate->products[i] = NULL;
2312 rate->nproducts = 0;
2314 bool lgProd =
false;
2316 for(
int i=0;!i || label[i-1]!=
'\0';i++)
2318 if(label[i] ==
',' || label[i] ==
'=' || label[i] ==
'\0')
2323 if(
trace.lgTraceMole )
2324 fprintf(
ioQQQ,
"Mole_reactions: ignoring reaction %s (species %s not active)\n",label,buf.c_str());
2332 fprintf(stderr,
"Mole_reactions: Too many reactants in %s, only %d allowed\n",label,
MAXREACTANTS);
2335 rate->reactants[rate->nreactants] = sp;
2342 fprintf(stderr,
"Mole_reactions: Too many products in %s, only %d allowed\n",label,
MAXPRODUCTS);
2345 rate->products[rate->nproducts] = sp;
2351 if (label[i] !=
'>')
2353 fprintf(
ioQQQ,
"Format error in reaction %s\n",label);
2365 ASSERT( rate->nreactants );
2366 ASSERT( rate->nproducts );
2376 ratefunc rate = findfunc(
"null");
2377 rate->label = label;
2399 for(
long i=0; i<rate->nreactants; ++i )
2401 newLabel += rate->reactants[i]->label;
2402 if( i != rate->nreactants-1 )
2406 for(
long i=0; i<rate->nproducts; ++i )
2408 newLabel += rate->products[i]->label;
2409 if( i != rate->nproducts-1 )
2414 rate->label = newLabel;
2423 for (
long k=0;k<rate->nreactants;k++)
2425 rate->rvector[k] = NULL;
2426 rate->rvector_excit[k] = NULL;
2429 for (
long k=0;k<rate->nproducts;k++)
2431 rate->pvector[k] = NULL;
2432 rate->pvector_excit[k] = NULL;
2436 for (
long i=0;i<rate->nproducts;i++)
2438 if (rate->pvector[i] == NULL)
2440 for (
long k=0;k<rate->nreactants;k++)
2442 if (rate->rvector[k] == NULL)
2444 if (rate->products[i] == rate->reactants[k])
2446 rate->rvector[k] = rate->products[i];
2447 rate->pvector[i] = rate->reactants[k];
2456 for (
long i=0;i<rate->nproducts;i++)
2458 if (rate->pvector[i] == NULL)
2460 for (
long k=0;k<rate->nreactants;k++)
2462 if (rate->rvector[k] == NULL)
2464 if (rate->products[i]->groupnum != -1 &&
2465 rate->products[i]->groupnum ==
2466 rate->reactants[k]->groupnum)
2468 rate->rvector[k] = rate->products[i];
2469 rate->pvector[i] = rate->reactants[k];
2478 for (
long i=0;i<rate->nproducts;i++)
2480 if (rate->pvector[i] == NULL && rate->pvector_excit[i] == NULL)
2482 for (
long k=0;k<rate->nreactants;k++)
2484 if (rate->rvector[k] == NULL && rate->rvector_excit[k] == NULL )
2488 rate->rvector_excit[k] = rate->products[i];
2489 rate->pvector_excit[i] = rate->reactants[k];
2546 sparsefp = fopen(
"sparse.pbm",
"w");
2547 fprintf(sparsefp,
"P4\n%d %d\n",
2555 ch = (ch << 1) | (c[i][j] != 0.0);
2577 int dcharge, n,
sign;
2581 for (n=0;n<rate->nreactants;n++)
2583 for(
nAtoms_i it = rate->reactants[n]->nAtom.begin(); it != rate->reactants[n]->nAtom.end(); ++it )
2584 nel[it->first] += it->second;
2585 dcharge += rate->reactants[n]->charge;
2587 for (n=0;n<rate->nproducts;n++)
2589 for(
nAtoms_i it = rate->products[n]->nAtom.begin(); it != rate->products[n]->nAtom.end(); ++it )
2590 nel[it->first] -= it->second;
2591 dcharge -= rate->products[n]->charge;
2595 fprintf(stderr,
"Reaction %s charge out of balance by %d\n",
2596 rate->label.c_str(),dcharge);
2600 for(
nAtoms_i it = nel.begin(); it != nel.end(); ++it )
2608 fprintf(stderr,
"Error: reaction %s %s %d of element %s\n",
2609 rate->label.c_str(),
sign==1?
"destroys":
"creates",
2611 it->first->label().c_str() );
2619#define gzopen(file,mode) fopen(file,mode)
2620#define gzgets(fp,buf,siz) fgets(buf,siz,fp)
2621#define gzclose(fp) fclose(fp)
2629 class formula_species {
2634 bool operator< (
const formula_species &a,
const formula_species &b)
2639 if (a.reactants[i]<b.reactants[i])
2641 else if (a.reactants[i]>b.reactants[i])
2646 if (a.products[i]<b.products[i])
2648 else if (a.products[i]>b.products[i])
2654 struct udfa_reaction_s {
2658 double a, b, c, trange[2];
2662static map <formula_species,struct udfa_reaction_s *>
udfatab;
2672 fprintf(stderr,
"Error, could not read %s\n",file);
2685#define FLTEQ(a,b) (fabs((a)-(b)) <= 1e-6*fabs((a)+(b)))
2689 struct udfa_reaction_s *r;
2690 map <formula_species, struct udfa_reaction_s *>::iterator p;
2691 unsigned int havespecies = 1, i, n;
2696 r = (
struct udfa_reaction_s *)
MALLOC(
sizeof(
struct udfa_reaction_s));
2703 r->l.reactants[n] = NULL;
2715 if (!strncmp(f,
"CRPHOT",6))
2724 r->l.products[n] = NULL;
2762 r->source = f[0]?f[0]:
'?';
2765 r->trange[n] = atof(f);
2777 fprintf(stderr,
"Duplicate reaction\n");
2798 map<formula_species, struct udfa_reaction_s *>::iterator p;
2799 struct udfa_reaction_s *u;
2803 s.reactants[n] = rate->reactants[n];
2807 s.products[n] = rate->products[n];
2818 if (
FLTEQ(rate->a,u->a) &&
FLTEQ(rate->b,u->b) &&
FLTEQ(rate->c,u->c))
2834 ratefunc findfunc (
const char name[])
2842 enum { DEBUG_MOLE =
false };
2854 long index = rate.
index;
2857 mole.reaction_rks[index] = newrk;
2860 if (fabs(newrk-oldrk) > 0.1*newrk)
2861 fprintf(
ioQQQ,
"%s: %15.8g => %15.8g\n",
2862 rate.
label.c_str(),oldrk,newrk);
2868 enum { DEBUG_MOLE =
false };
2872 if (
mole.old_reaction_rks.size() == 0 )
2875 mole.old_reaction_rks.resize(
mole.reaction_rks.size());
2881 double bigchange = 0.;
2882 unsigned long bigindex = ULONG_MAX;
2883 for (
unsigned long index=0; index<
mole.reaction_rks.size(); ++index)
2885 double oldrk =
mole.old_reaction_rks[index],
2886 newrk =
mole.reaction_rks[index],
2887 sum = oldrk+newrk, diff = newrk-oldrk;
2890 double change = fabs(diff)/sum;
2891 if (change > bigchange)
2903 if (bigindex == (
unsigned) rate.
index)
2905 double oldrk =
mole.old_reaction_rks[bigindex],
2906 newrk =
mole.reaction_rks[bigindex];
2909 pct = 100.*(newrk-oldrk)/oldrk;
2910 fprintf(
ioQQQ,
"Zone %ld, big chemistry rate change %s:"
2911 " %15.8g => %15.8g (%6.2g%%)\n",
2919 for (
unsigned long index=0; index<
mole.reaction_rks.size(); ++index)
2921 mole.old_reaction_rks[index] =
mole.reaction_rks[index];
2927 double T_ortho_para_crit;
2937# ifndef IGNORE_QUANTUM_HEATING
2945 h2.rate_grain_op_conserve = 0.;
2947 h2.rate_grain_J1_to_J0 = 0.;
2950 for(
size_t nd=0; nd <
gv.bin.size(); nd++ )
2952# ifndef IGNORE_QUANTUM_HEATING
2954 double *qtemp, *qprob;
2955 bool lgUseQHeat =
gv.lgGrainPhysicsOn &&
gv.bin[nd]->lgQHeat;
2964 double sticking_probability_H = 1./(1. + 0.04*sqrt(
gv.bin[nd]->tedust+
phycon.te) +
2967# ifndef IGNORE_QUANTUM_HEATING
2971 qtemp = (
double*)
MALLOC((
size_t)(
NQGRID*
sizeof(
double)));
2972 qprob = (
double*)
MALLOC((
size_t)(
NQGRID*
sizeof(
double)));
2974 qheat(qtemp,qprob,&qnbin,nd);
2976 if(
gv.bin[nd]->lgUseQHeat )
2984 qtemp[0] =
gv.bin[nd]->tedust;
2987 gv.bin[nd]->rate_h2_form_grains_HM79 = 0.;
2989 for( k=0; k < qnbin; k++ )
2999 double conversion_efficiency_HM79 = 1/(1. + 1e4*
sexp(920./qtemp[k]));
3000 sticking_probability_H = 1./(1. + 0.04*sqrt(qtemp[k]+
phycon.te) +
3003 gv.bin[nd]->rate_h2_form_grains_HM79 += qprob[k] * sticking_probability_H *
3004 conversion_efficiency_HM79;
3010 gv.bin[nd]->rate_h2_form_grains_HM79 *= 0.5 * AveVelH *
3011 gv.bin[nd]->IntArea/4. *
gv.bin[nd]->cnv_H_pCM3;
3013 ASSERT(
gv.bin[nd]->rate_h2_form_grains_HM79 > 0. );
3026 double conversion_efficiency_HM79 = 1/(1. + 1e4*
sexp(920./
gv.bin[nd]->tedust));
3031 gv.bin[nd]->rate_h2_form_grains_HM79 = 0.5 * AveVelH *
gv.bin[nd]->IntArea/4. *
3033 gv.bin[nd]->cnv_H_pCM3 * sticking_probability_H * conversion_efficiency_HM79;
3034 ASSERT(
gv.bin[nd]->rate_h2_form_grains_HM79 > 0. );
3037# ifndef IGNORE_QUANTUM_HEATING
3046 double sqrt_term = 35.399494936611667;
3048 gv.bin[nd]->rate_h2_form_grains_CT02 = 0.;
3050 for( k=0; k < qnbin; k++ )
3052 double beta_alpha = 0.25 * sqrt_term *
sexp(200./qtemp[k] );
3054 double xi = 1./ (1. + 1.3e13*
sexp(1.5*1e4/qtemp[k])*sqrt_term/(2.*f) );
3056 double beta = 3e12 *
sexp( 320. / qtemp[k] );
3059 double recombination_efficiency_CT02 = xi / (1. + 0.005*f/2./
SDIV(beta) + beta_alpha );
3060 sticking_probability_H = 1./(1. + 0.04*sqrt(qtemp[k]+
phycon.te) +
3066 gv.bin[nd]->rate_h2_form_grains_CT02 += qprob[k] * sticking_probability_H *
3067 recombination_efficiency_CT02;
3074 gv.bin[nd]->rate_h2_form_grains_CT02 *= 0.5 * AveVelH *
3075 gv.bin[nd]->IntArea/4. *
gv.bin[nd]->cnv_H_pCM3;
3077 ASSERT(
gv.bin[nd]->rate_h2_form_grains_CT02 > 0. );
3091 double sqrt_term = 35.399494936611667;
3092 double beta_alpha = 0.25 * sqrt_term *
sexp(200./
gv.bin[nd]->tedust );
3094 double xi = 1./ (1. + 1.3e13*
sexp(1.5*1e4/
gv.bin[nd]->tedust )*sqrt_term/(2.*f) );
3096 double beta = 3e12 *
sexp( 320. /
gv.bin[nd]->tedust );
3099 double recombination_efficiency_CT02 = xi / (1. + 0.005*f/2./
SDIV(beta) + beta_alpha );
3107 double Td =
gv.bin[nd]->tedust;
3108 recombination_efficiency_CT02 = 1./(
3109 ( 1. + 4.609569629185726e24*
sexp(45000./Td) ) *
3110 ( 1. +
sexp(800./Td) / (0.5389970511202651 *
sexp(540./Td) + 5.6333909478365e-14*sqrt(Td)) )
3117 double Td =
gv.bin[nd]->tedust;
3118 recombination_efficiency_CT02 = 1./(
3119 ( 1. + 6.998161265697586e24*
sexp(45000./Td) ) *
3120 ( 1. +
sexp(450./Td) / (0.4266153643741715 *
sexp(340./Td) + 2.5335919594255e-14*sqrt(Td)) )
3130 gv.bin[nd]->rate_h2_form_grains_CT02 = 0.5 * AveVelH *
gv.bin[nd]->IntArea/4. *
3131 gv.bin[nd]->cnv_H_pCM3 * sticking_probability_H * recombination_efficiency_CT02;
3132 ASSERT(
gv.bin[nd]->rate_h2_form_grains_CT02 > 0. );
3135# ifndef IGNORE_QUANTUM_HEATING
3137 sticking_probability_H = 1./(1. + 0.04*sqrt(
gv.bin[nd]->tedust+
phycon.te) +
3151 h2.rate_grain_op_conserve += AveVelH2 *
gv.bin[nd]->IntArea/4. *
3152 gv.bin[nd]->cnv_H_pCM3 * sticking_probability_H;
3177 T_ortho_para_crit = 2. *
hmi.Tad / log(
POW2(60. *1.1e11)*
hmi.Tad);
3178 if(
gv.bin[nd]->tedust < T_ortho_para_crit )
3180 double efficiency_opr =
sexp(60.*1.1e11*sqrt(
hmi.Tad)*
sexp(
hmi.Tad/
gv.bin[nd]->tedust));
3183 h2.rate_grain_J1_to_J0 += AveVelH2 *
gv.bin[nd]->IntArea/4. *
3184 gv.bin[nd]->cnv_H_pCM3 * sticking_probability_H * efficiency_opr;
3196 h2.rate_grain_op_conserve *=
h2.lgH2_grain_deexcitation;
3197 h2.rate_grain_J1_to_J0 *=
h2.lgH2_grain_deexcitation;
3203 for(
size_t nd=0; nd <
gv.bin.size(); nd++ )
3205 gv.bin[nd]->rate_h2_form_grains_CT02 = 0.;
3206 gv.bin[nd]->rate_h2_form_grains_HM79 = 0.;
3209 h2.rate_grain_op_conserve = 0.;
3211 h2.rate_grain_J1_to_J0 = 0.;
3218 gv.rate_h2_form_grains_used_total = 0.;
3219 for(
size_t nd=0; nd <
gv.bin.size(); nd++ )
3221 if(
hmi.chJura ==
'C' )
3226 gv.bin[nd]->rate_h2_form_grains_used =
3227 gv.bin[nd]->rate_h2_form_grains_CT02*
hmi.ScaleJura;
3228 gv.rate_h2_form_grains_used_total +=
gv.bin[nd]->rate_h2_form_grains_used;
3230 else if(
hmi.chJura ==
'T' )
3233 gv.bin[nd]->rate_h2_form_grains_used =
3234 gv.bin[nd]->rate_h2_form_grains_HM79*
hmi.ScaleJura;
3235 gv.rate_h2_form_grains_used_total +=
gv.bin[nd]->rate_h2_form_grains_used;
3237 else if(
hmi.chJura ==
'S' )
3241 gv.bin[nd]->rate_h2_form_grains_used =
3244 gv.rate_h2_form_grains_used_total +=
gv.bin[nd]->rate_h2_form_grains_used;
3249 else if(
hmi.chJura ==
'F' )
3253 gv.bin[nd]->rate_h2_form_grains_used =
hmi.rate_h2_form_grains_set*
dense.gas_phase[
ipHYDROGEN] /
gv.bin.size();
3254 gv.rate_h2_form_grains_used_total +=
gv.bin[nd]->rate_h2_form_grains_used;
3257 ASSERT(
gv.rate_h2_form_grains_used_total >= 0. );
3259# ifndef IGNORE_QUANTUM_HEATING
3260 printf(
" fnzone %.2f H2 rate %.4e\n",
fnzone,
gv.rate_h2_form_grains_used_total );
3270 static double teused=-1;
3280 cr_H2dis_ELWERT_H2g,
3281 cr_H2dis_ELWERT_H2s;
3298 if(
hmi.exphmi > 0. )
3305 hmi.rel_pop_LTE_Hmin = 0.;
3315 hmi.rel_pop_LTE_H2p =
SAHA/(
phycon.te32*exphp)*(4./(2.*1.))*3.634e-5;
3319 hmi.rel_pop_LTE_H2p = 0.;
3329 hmi.rel_pop_LTE_H3p =
SAHA/(
phycon.te32*ex3hp)*(4./(2.*1.))*3.634e-5;
3333 hmi.rel_pop_LTE_H3p = 0.;
3342 if(
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
3345 hmi.rel_pop_LTE_H2g =
h2.rel_pop_LTE_g;
3346 hmi.rel_pop_LTE_H2s =
h2.rel_pop_LTE_s;
3359 hmi.rel_pop_LTE_H2g =
SAHA/(
phycon.te32*exph2)*(1./(2.*2.))*3.634e-5;
3363 hmi.rel_pop_LTE_H2g = 0.;
3375 hmi.rel_pop_LTE_H2s =
SAHA/(
phycon.te32*exph2s)*(1./(2.*2.))*3.634e-5;
3379 hmi.rel_pop_LTE_H2s = 0.;
3388 enum {DEBUG_LOC=
false};
3393 fprintf(
ioQQQ,
"ph2lteee\t%.2e\t%.1e\t%.1e\n",
3394 hmi.rel_pop_LTE_H2g,
3406 if( (*diatom)->lgEnabled &&
mole_global.lgStancil )
3407 (*diatom)->Mol_Photo_Diss_Rates();
3431 0.055502 , &
hmi.HMinus_induc_rec_rate , &
hmi.HMinus_induc_rec_cooling, &photoHeat );
3439 enum {DEBUG_LOC=
false};
3443 fprintf(
ioQQQ,
"hminphoto\t%li\t%li\t%.2e\n",
nzone,
conv.nPres2Ioniz ,
hmi.HMinus_photo_rate );
3448 hmi.HMinus_induc_rec_rate *=
hmi.rel_pop_LTE_Hmin*
dense.eden;
3457 enum {DEBUG_LOC=
false};
3459 if( DEBUG_LOC &&
nzone>400)
3467 hmi.HMinus_photo_rate,
3469 hmi.HMinus_photo_rate*0.05);
3483 hmi.UV_Cont_rel2_Habing_TH85_depth = 0.;
3489 for( i=
rfield.ipG0_TH85_lo; i <
rfield.ipG0_TH85_hi; ++i )
3491 hmi.UV_Cont_rel2_Habing_TH85_depth += ((
rfield.flux[0][i-1]) + (
rfield.ConInterOut[i-1])+
3499 hmi.UV_Cont_rel2_Habing_TH85_depth = (
realnum)(
hmi.UV_Cont_rel2_Habing_TH85_depth*
EN1RYD/1.6e-3);
3503 hmi.UV_Cont_rel2_Habing_TH85_face =
hmi.UV_Cont_rel2_Habing_TH85_depth;
3508 hmi.UV_Cont_rel2_Habing_spec_depth = 0.;
3509 for( i=
rfield.ipG0_spec_lo; i <
rfield.ipG0_spec_hi; ++i )
3511 hmi.UV_Cont_rel2_Habing_spec_depth += (
rfield.flux[0][i-1] +
rfield.ConInterOut[i-1]+
3514 hmi.UV_Cont_rel2_Habing_spec_depth = (
realnum)(
hmi.UV_Cont_rel2_Habing_spec_depth*
EN1RYD/1.6e-3);
3518 if( !
conv.nTotalIoniz )
3520 hmi.UV_Cont_rel2_Draine_DB96_face = 0.f;
3522 for( i=
rfield.ipG0_DB96_lo; i <
rfield.ipG0_DB96_hi; ++i )
3524 hmi.UV_Cont_rel2_Draine_DB96_face += (
rfield.flux[0][i-1] +
rfield.ConInterOut[i-1]+
3535 hmi.UV_Cont_rel2_Draine_DB96_face =
hmi.UV_Cont_rel2_Draine_DB96_face/(1.232e7f * 1.71f);
3553 hmi.H2_Solomon_dissoc_rate_TH85_H2g = 3.4e-11 *
hmi.UV_Cont_rel2_Habing_TH85_depth * h2esc;
3554 hmi.H2_Solomon_dissoc_rate_TH85_H2s = 3.4e-11 *
hmi.UV_Cont_rel2_Habing_TH85_depth * h2esc;
3555 hmi.H2_H2g_to_H2s_rate_TH85 =
hmi.H2_Solomon_dissoc_rate_TH85_H2g*9.;
3558 hmi.H2_Solomon_dissoc_rate_BHT90_H2g = 3.4e-11 *
hmi.UV_Cont_rel2_Habing_TH85_depth * h2esc;
3559 hmi.H2_Solomon_dissoc_rate_BHT90_H2s = 3.4e-11 *
hmi.UV_Cont_rel2_Habing_TH85_depth * h2esc;
3560 hmi.H2_H2g_to_H2s_rate_BHT90 =
hmi.H2_Solomon_dissoc_rate_BHT90_H2g*9.;
3568 double sqrtx = sqrt(1.+x);
3572 double fshield = 0.965/
POW2(1.+x/b5) + 0.035/sqrtx *
3580 hmi.UV_Cont_rel2_Draine_DB96_depth =
hmi.UV_Cont_rel2_Draine_DB96_face *
3589 hmi.H2_Solomon_dissoc_rate_BD96_H2g = 4.6e-11 *
hmi.UV_Cont_rel2_Draine_DB96_depth * fshield;
3590 hmi.H2_Solomon_dissoc_rate_BD96_H2s = 4.6e-11 *
hmi.UV_Cont_rel2_Draine_DB96_depth * fshield;
3595 hmi.H2_Solomon_dissoc_rate_BD96_H2g = 5.18e-11* (
hmi.UV_Cont_rel2_Habing_TH85_face/1.66f)
3596 *
sexp(3.02*
rfield.extin_mag_V_point)* fshield;
3597 hmi.H2_Solomon_dissoc_rate_BD96_H2s = 5.18e-11* (
hmi.UV_Cont_rel2_Habing_TH85_face/1.66f)
3598 *
sexp(3.02*
rfield.extin_mag_V_point)* fshield;
3605 hmi.H2_H2g_to_H2s_rate_BD96 = 5.67*
hmi.H2_Solomon_dissoc_rate_BD96_H2g;
3624 double x_H2g, x_H2s,
3625 fshield_H2g, fshield_H2s,
3627 static double a_H2g, a_H2s,
3644 if(
hmi.UV_Cont_rel2_Draine_DB96_face <= 1.)
3648 else if(
hmi.UV_Cont_rel2_Draine_DB96_face >= 1e7)
3654 log_G0_face = log10(
hmi.UV_Cont_rel2_Draine_DB96_face);
3660 a_H2g = 0.06 * log_G0_face + 1.32;
3661 a_H2s = 0.375 * log_G0_face + 2.125;
3663 e1_H2g = -0.05 * log_G0_face + 2.25;
3664 e1_H2s = -0.125 * log_G0_face + 2.625;
3666 e2_H2g = -0.005 * log_G0_face + 0.625;
3668 b_H2g = -4.0e-3 * log_G0_face + 3.2e-2;
3675 k_f_H2s =
MAX2(0.1,2.375 * log_G0_face - 1.875 );
3678 k_H2g_to_H2s =
MAX2(1.,-1.75 * log_G0_face + 11.25);
3686 f_H2s = k_f_H2s * pow((
double)
hmi.UV_Cont_rel2_Draine_DB96_depth, 0.2 );
3693 fshield_H2g = 0.965/pow(1.+x_H2g/b5,e1_H2g) + b_H2g/pow(1.+x_H2g/b5,e2_H2g);
3694 fshield_H2s = 0.965/pow(1.+x_H2s/b5,e1_H2s);
3697 hmi.H2_Solomon_dissoc_rate_ELWERT_H2g = a_H2g * 4.6e-11 * fshield_H2g *
hmi.UV_Cont_rel2_Draine_DB96_depth;
3698 hmi.H2_Solomon_dissoc_rate_ELWERT_H2s = a_H2s * 4.6e-11 * fshield_H2s * (
hmi.UV_Cont_rel2_Draine_DB96_depth + f_H2s);
3701 hmi.H2_H2g_to_H2s_rate_ELWERT = k_H2g_to_H2s *
hmi.H2_Solomon_dissoc_rate_ELWERT_H2g;
3706 hmi.H2_photodissoc_ELWERT_H2s =
hmi.UV_Cont_rel2_Draine_DB96_depth*1e-11;
3707 hmi.H2_photodissoc_ELWERT_H2g =
hmi.H2_photodissoc_ELWERT_H2s * 1.0e-10;
3711 hmi.H2_Solomon_dissoc_rate_ELWERT_H2g = 0.;
3712 hmi.H2_Solomon_dissoc_rate_ELWERT_H2s = 0.;
3713 hmi.H2_photodissoc_ELWERT_H2s = 0.;
3714 hmi.H2_photodissoc_ELWERT_H2g = 0.;
3718 hmi.H2_photodissoc_TH85 =
hmi.UV_Cont_rel2_Habing_TH85_depth*1e-11;
3721 hmi.H2_photodissoc_BHT90 =
hmi.UV_Cont_rel2_Habing_TH85_depth*1e-11;
3737 cr_H2dis_ELWERT_H2g =
secondaries.x12tot*5e-8 *
hmi.lgLeidenCRHack;
3738 cr_H2dis_ELWERT_H2s =
secondaries.x12tot*4e-2 *
hmi.lgLeidenCRHack;
3750 if(
h2.lgEnabled &&
h2.lgEvaluated &&
hmi.lgH2_Chemistry_BigH2 )
3758 hmi.H2_Solomon_dissoc_rate_used_H2g =
h2.Solomon_dissoc_rate_g;
3759 ASSERT(
hmi.H2_Solomon_dissoc_rate_used_H2g >= 0. );
3761 hmi.H2_Solomon_dissoc_rate_used_H2s =
h2.Solomon_dissoc_rate_s;
3762 ASSERT(
hmi.H2_Solomon_dissoc_rate_used_H2s >= 0. );
3765 hmi.H2_H2g_to_H2s_rate_used =
h2.gs_rate();
3766 ASSERT(
hmi.H2_H2g_to_H2s_rate_used >= 0. );
3772 hmi.H2_photodissoc_used_H2s =
h2.photodissoc_BigH2_H2s;
3775 hmi.H2_photodissoc_used_H2s =
h2.Cont_Dissoc_Rate_H2s;
3776 ASSERT(
hmi.H2_photodissoc_used_H2s >= 0. );
3780 hmi.H2_photodissoc_used_H2g =
h2.photodissoc_BigH2_H2g;
3783 hmi.H2_photodissoc_used_H2g =
h2.Cont_Dissoc_Rate_H2g;
3784 ASSERT(
hmi.H2_photodissoc_used_H2g >= 0. );
3786 else if(
hmi.chH2_small_model_type ==
'T' )
3790 hmi.H2_Solomon_dissoc_rate_used_H2g =
hmi.H2_Solomon_dissoc_rate_TH85_H2g + cr_H2dis;
3792 hmi.H2_Solomon_dissoc_rate_used_H2s =
hmi.H2_Solomon_dissoc_rate_TH85_H2s + cr_H2dis;
3793 hmi.H2_H2g_to_H2s_rate_used =
hmi.H2_H2g_to_H2s_rate_TH85 + cr_H2s;
3797 hmi.H2_photodissoc_used_H2s =
hmi.H2_photodissoc_TH85;
3800 hmi.H2_photodissoc_used_H2g =
hmi.H2_photodissoc_TH85*1.0e-10f;
3803 else if(
hmi.chH2_small_model_type ==
'H' )
3808 hmi.H2_Solomon_dissoc_rate_used_H2g =
hmi.H2_Solomon_dissoc_rate_BHT90_H2g + cr_H2dis;
3810 hmi.H2_Solomon_dissoc_rate_used_H2s =
hmi.H2_Solomon_dissoc_rate_BHT90_H2s + cr_H2dis;
3811 hmi.H2_H2g_to_H2s_rate_used =
hmi.H2_H2g_to_H2s_rate_BHT90 + cr_H2s;
3815 hmi.H2_photodissoc_used_H2s =
hmi.H2_photodissoc_BHT90;
3818 hmi.H2_photodissoc_used_H2g =
hmi.H2_photodissoc_BHT90*1.0e-10f;
3821 else if(
hmi.chH2_small_model_type ==
'B' )
3825 hmi.H2_Solomon_dissoc_rate_used_H2g =
hmi.H2_Solomon_dissoc_rate_BD96_H2g + cr_H2dis;
3827 hmi.H2_Solomon_dissoc_rate_used_H2s =
hmi.H2_Solomon_dissoc_rate_BD96_H2s + cr_H2dis;
3829 hmi.H2_H2g_to_H2s_rate_used =
hmi.H2_H2g_to_H2s_rate_BD96 + cr_H2s;
3834 hmi.H2_photodissoc_used_H2s =
hmi.H2_photodissoc_TH85;
3837 hmi.H2_photodissoc_used_H2g =
hmi.H2_photodissoc_TH85*1.0e-10f;
3839 else if(
hmi.chH2_small_model_type ==
'E' )
3843 hmi.H2_Solomon_dissoc_rate_used_H2g =
hmi.H2_Solomon_dissoc_rate_ELWERT_H2g + cr_H2dis_ELWERT_H2g;
3844 hmi.H2_Solomon_dissoc_rate_used_H2s =
hmi.H2_Solomon_dissoc_rate_ELWERT_H2s + cr_H2dis_ELWERT_H2s;
3845 hmi.H2_H2g_to_H2s_rate_used =
hmi.H2_H2g_to_H2s_rate_ELWERT + cr_H2s;
3850 hmi.H2_photodissoc_used_H2s =
hmi.H2_photodissoc_ELWERT_H2s;
3851 hmi.H2_photodissoc_used_H2g =
hmi.H2_photodissoc_ELWERT_H2g;
3858 enum {DEBUG_LOC=
false};
3860 if( DEBUG_LOC &&
h2.lgEnabled )
3862 fprintf(
ioQQQ,
" Solomon H2 dest rates: TH85 %.2e BD96 %.2e Big %.2e excit rates: TH85 %.2e Big %.2e\n",
3863 hmi.H2_Solomon_dissoc_rate_TH85_H2g,
3864 hmi.H2_Solomon_dissoc_rate_BD96_H2g,
3865 h2.Solomon_dissoc_rate_g,
3866 hmi.H2_H2g_to_H2s_rate_TH85 ,
3881 return &(*p->second);
3957 for(
int i=0;i<rate.
nreactants && ipthis == -1;i++)
3966 double ratevi = rate.
a * rate.
rk();
3990 fprintf(stderr,
"=>");
3994 fprintf(stderr,
"\n");
4028 int ipspin = 0, ipfreein = 0;
4036 int ipspout = 0, ipfreeout = 0;
4046 int newsp = ipspout-ipspin;
4051 int nbondsbroken = ipfreeout-ipfreein;
4052 if (nbondsbroken <= 0)
4055 double fracbroken = nbondsbroken/((double)ipfreeout);
4056 ASSERT( fracbroken <= 1.0 );
4067 double ratesp = ratevi*newsp;
4069 ratesp *= fracbroken;
4103 double ratevi = rate.
a * rate.
rk();
4108 ratev += ipthis*ratevi;
4130 double heating = 0.;
4131 map<double,string> heatMap;
4143 bool lgCanSkip =
false;
4171 realnum reaction_enthalpy = 0.;
4180 for(
long i=0; i < rate.
nproducts; ++i )
4190 double heat = reaction_enthalpy*rate_tot*(1e10/
AVOGADRO);
4191 heatMap[heat] = rate.
label;
4201 for( map<double,string>::reverse_iterator it = heatMap.rbegin(); it != heatMap.rend(); ++it, ++index )
4203 fprintf(
ioQQQ,
"DEBUGGG heat %li\t%li\t%.6e\t%s\n", index,
nzone, it->first, it->second.c_str() );
4208 for( map<double,string>::iterator it = heatMap.begin(); it != heatMap.end(); ++it, ++index )
4210 if( it->first >= 0. )
4212 fprintf(
ioQQQ,
"DEBUGGG cool %li\t%li\t%.6e\t%s\n", index,
nzone, it->first, it->second.c_str() );
4221void mole_punch(FILE *punit,
const char speciesname[],
const char args[],
bool lgHeader,
bool lgData,
double depth)
4233 fprintf (punit,
"#Depth");
4237 fprintf (punit,
"%.5e",depth);
4252 ( strcmp( args,
"CATA" )==0 && rate.
rvector[i]!=NULL ) ||
4253 strcmp( args,
"ALL " )==0 )
4264 ( strcmp( args,
"CATA" )==0 && rate.
pvector[i]!=NULL ) ||
4265 strcmp( args,
"ALL " )==0 )
4274 fprintf(punit,
"\t%s",rate.
label.c_str());
4278 ratevi =
mole.reaction_rks[ rate.
index ];
4283 fprintf(punit,
"\t%.3e",ratevi);
4287 fprintf(punit,
"\n");
sys_float sexp(sys_float x)
bool fp_equal(sys_float x, sys_float y, int n=3)
NORETURN void TotalInsanity(void)
sys_float SDIV(sys_float x)
#define DEBUG_ENTRY(funcname)
double powi(double, long int)
static t_version & Inst()
double Solomon_dissoc_rate_g
double photodissoc_BigH2_H2g
molecule * pvector_excit[MAXPRODUCTS]
virtual double rk() const =0
molecule * pvector[MAXPRODUCTS]
molecule * rvector_excit[MAXREACTANTS]
molecule * reactants[MAXREACTANTS]
molecule * rvector[MAXREACTANTS]
molecule * products[MAXPRODUCTS]
map< const count_ptr< chem_atom >, int, element_pointer_value_less > nAtomsMap
bool isMonatomic(void) const
double **** PhotoRate_Shell
static void sort(MoleculeList::iterator start, MoleculeList::iterator end)
double dissoc_rate(const char chSpecies[]) const
vector< double > reaction_rks
double sink_rate(const molecule *const sp, const mole_reaction &rate) const
double source_rate_tot(const char chSpecies[]) const
double findrate(const char buf[]) const
double findrk(const char buf[]) const
double chem_heat(void) const
valarray< class molezone > species
double sink_rate_tot(const char chSpecies[]) const
void GammaPrt(long int ipLoEnr, long int ipHiEnr, long int ipOpac, FILE *ioFILE, double total, double threshold)
double GammaK(long int ipLoEnr, long int ipHiEnr, long int ipOpac, double yield1, t_phoHeat *photoHeat)
double GammaBn(long int ipLoEnr, long int ipHiEnr, long int ipOpac, double thresh, double *ainduc, double *rcool, t_phoHeat *photoHeat)
bool operator<(const count_ptr< T > &a, const count_ptr< T > &b)
FILE * open_data(const char *fname, const char *mode, access_scheme scheme)
UNUSED const realnum BIGFLOAT
realnum GetDopplerWidth(realnum massAMU)
realnum GetAveVelocity(realnum massAMU)
void qheat(vector< double > &, vector< double > &, long *, size_t)
diatomics hd("hd", 4100., &hmi.HD_total, Yan_H2_CS)
vector< diatomics * > diatoms
diatomics h2("h2", 4100., &hmi.H2_total, Yan_H2_CS)
vector< diatomics * >::iterator diatom_iter
t_iso_sp iso_sp[NISO][LIMELM]
void iso_photo(long ipISO, long nelem)
t_mole_global mole_global
void create_isotopologues_one(ChemAtomList &atoms, vector< int > &numAtoms, string atom_old, string atom_new, string embellishments, vector< string > &newLabels)
molezone * findspecieslocal(const char buf[])
mole_reaction * mole_findrate_s(const char buf[])
ChemAtomList unresolved_atom_list
molecule::nAtomsMap::iterator nAtoms_i
bool lgDifferByExcitation(const molecule &mol1, const molecule &mol2)
bool parse_species_label(const char label[], ChemAtomList &atomsLeftToRight, vector< int > &numAtoms, string &embellishments)
molecule * findspecies(const char buf[])
vector< count_ptr< chem_atom > > ChemAtomList
map< string, count_ptr< mole_reaction > >::iterator mole_reaction_i
map< string, count_ptr< mole_reaction > >::const_iterator mole_reaction_ci
void mole_create_react(void)
STATIC double mole_get_equilibrium_constant(const char buf[])
STATIC void newreact(const char label[], const char fun[], double a, double b, double c)
STATIC void compare_udfa(const count_ptr< mole_reaction > &rate)
STATIC bool lgReactBalance(const count_ptr< mole_reaction > &rate)
STATIC void read_data(const char file[], void(*parse)(char *s))
STATIC string canonicalize_reaction_label(const char label[])
STATIC long parse_reaction(count_ptr< mole_reaction > &rate, const char label[])
STATIC void mole_check_reverse_reactions(void)
void mole_punch(FILE *punit, const char speciesname[], const char args[], bool lgHeader, bool lgData, double depth)
STATIC char * getcsvfield(char **s, char c)
void mole_rk_bigchange(void)
STATIC void mole_h_reactions(void)
STATIC void mole_generate_isotopologue_reactions(string atom_old, string atom_new)
STATIC void plot_sparsity(void)
STATIC void parse_base(char *s)
mole_reaction * mole_findrate_s(const char buf[])
STATIC double mole_partition_function(const molecule *const sp)
STATIC void mole_h2_grain_form(void)
double frac_H2star_hminus(void)
STATIC void canonicalize_reaction(count_ptr< mole_reaction > &rate)
static map< formula_species, struct udfa_reaction_s * > udfatab
STATIC void register_reaction_vectors(count_ptr< mole_reaction > rate)
void mole_update_rks(void)
void mole_cmp_num_in_out_reactions()
STATIC void parse_udfa(char *s)
map< string, count_ptr< mole_reaction > > functab
map< string, count_ptr< mole_reaction > > reactab
UNUSED const double BOLTZMANN
UNUSED const double ELECTRON_MASS
UNUSED const double EN1RYD
UNUSED const double EN1EV
UNUSED const double ATOMIC_MASS_UNIT
UNUSED const double HBAReV
UNUSED const double AVOGADRO
UNUSED const double HION_LTE_POP
double esc_PRD_1side(double tau, double a)
double esca0k2(double taume)
t_secondaries secondaries
double CharExcRecTo[NCX][LIMELM][LIMELM+1]
double CharExcIonOf[NCX][LIMELM][LIMELM+1]
double xIonDense[LIMELM][LIMELM+1]
realnum AtomicWeight[LIMELM]
bool lgH2_Chemistry_BigH2
realnum UV_Cont_rel2_Draine_DB96_depth
double H2_Solomon_dissoc_rate_used_H2g
double H2_H2g_to_H2s_rate_used
double HMinus_induc_rec_rate
double H2_Solomon_dissoc_rate_used_H2s
double H2_photodissoc_used_H2g
double H2_photodissoc_used_H2s