Cloudy
Spectral Synthesis Code for Astrophysics
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
h2_priv.h
Go to the documentation of this file.
1 /* This file is part of Cloudy and is copyright (C)1978-2023 by Gary J. Ferland and
2  * others. For conditions of distribution and use see copyright notice in license.txt */
3 
4 #ifndef H2_PRIV_H_
5 #define H2_PRIV_H_
6 
7 #include "transition.h"
8 #include "container_classes.h"
9 
10 typedef void (*linefunc)(const TransitionProxy& t,
11  bool lgShield_this_zone,
12  realnum pestrk,
13  realnum DopplerWidth);
14 class Parser;
15 class molecule;
16 class CollRateCoeffArray;
17 
20 const int N_X_COLLIDER = 5;
22 const int chN_X_COLLIDER = 10;
23 // colliders are (see h2.cpp diatoms_init)
24 // 0 H
25 // 1 He
26 // 2 H2 ortho
27 // 3 H2 para
28 // 4 H+
29 
31 const int nTE_HMINUS = 7;
32 
34 const int N_ELEC = 7;
35 
37 const realnum H2_logte_hminus[nTE_HMINUS] = {1.,1.47712,2.,2.47712,3.,3.47712,4.};
38 
39 struct diss_level
40 {
41  long n, v, j;
42 };
43 
44 class diss_tran
45 {
46 public:
47  explicit diss_tran( diss_level a, diss_level b )
48  {
49  initial = a;
50  final = b;
51  energies.clear();
52  xsections.clear();
53  rate_coeff = 0.;
54  };
55  diss_level initial, final;
56  vector<double> energies;
57  vector<double> xsections;
58  double rate_coeff;
59 };
60 
62 {
64  {
65  magic = 0;
66  filename = "";
67  };
68  long magic;
69  string filename;
70 };
71 
72 class diatomics
73 {
74 public:
75  double Abund() const
76  {
77  return *dense_total;
78  }
79  void GetIndices( long& ipHi, long& ipLo, const char* chLine, long& i ) const;
80  void CalcPhotoionizationRate(void);
81  double (*photoion_opacity_fun)( double energy );
82  long OpacityCreate( vector<double>& stack );
83  double GetHeatRate( const diss_tran& tran );
84  double GetDissociationRate( const diss_tran& tran );
85  double MolDissocOpacity( const diss_tran& tran, const double& Mol_Ene );
86  double Cont_Diss_Heat_Rate( void );
87  void Mol_Photo_Diss_Rates( void );
91  double GetExcitedElecDensity(void);
92  realnum GetXColden( long iVib, long iRot );
93 
94  long int getLine( long iElecHi, long iVibHi, long iRotHi, long iElecLo, long iVibLo, long iRotLo, double *relint, double *absint );
95 
96  /* compute rate coefficient for a single quenching collision */
97  realnum H2_CollidRateEvalOne( long iVibHi, long iRotHi, long iVibLo, long iRotLo, long ipHi, long ipLo, long nColl, double temp_K );
98 
99  void H2_Calc_Average_Rates( void );
100 
101  void H2_X_sink_and_source( void );
102  /*H2_X_coll_rate_evaluate find collisional rates within X -
103  * this is one time upon entry into H2_LevelPops */
104  void H2_X_coll_rate_evaluate( void );
105 
106  /*H2_Level_low_matrix evaluate lower populations within X */
107  /* total abundance within matrix */
108  void H2_Level_low_matrix(realnum abundance );
109 
110  /*H2_Read_Cosmicray_distribution read distribution function for H2 population following cosmic ray collisional excitation
111  void H2_Read_Cosmicray_distribution(void); */
112 
113  // read energies for all electronic levels
114  void H2_ReadEnergies();
115  void H2_ReadEnergies( long int nelec, vector<int>& n, vector<int>& v, vector<int>&J, vector<double>& eWN );
116 
120  void H2_ReadDissprob( long int nelec );
121 
123  void H2_CollidRateEvalAll( void );
124 
128  void H2_CollidRateRead( long int nColl );
129 
133  void H2_ReadTransprob( long int nelec, TransitionList &trans );
134 
136  void H2_Read_hminus_distribution(void);
137 
139  void mole_H2_form( void );
140 
142  void mole_H2_LTE( void );
143 
146  void H2_Solomon_rate( void );
147 
149  double gs_rate( void );
150 
153  void H2_zero_pops_too_low( void );
154 
156  void init(void);
157 
159  void H2_ContPoint( void );
160 
162  double H2_DR(void);
163 
165  double H2_Accel(void);
166 
168  void H2_RT_OTS( void );
169 
171  double H2_RadPress(void);
172 
174  void H2_LinesAdd(void);
175 
177  void H2_Reset( void );
178 
180  double H2_InterEnergy(void);
181 
185  void H2_Colden( const char *chLabel );
186 
188  void H2_Cooling(void);
189 
192  double LTE_Cooling_per_H2();
193 
198  void H2_Punch_line_data(
199  FILE* ioPUN ,
200  bool lgDoAll );
201 
207  void H2_PunchLineStuff( FILE * io , realnum xLimit , long index);
208 
210  void H2_RT_diffuse(void);
211 
214  void H2_RTMake( linefunc line_one );
215 
217  void H2_RT_tau_inc(void);
218 
220  void H2_Prt_Zone(void);
221 
222  // print departure coefficients for all X levels
223  void H2_PrtDepartCoef(void);
224 
226  void H2_LineZero( void );
227 
229  void H2_RT_tau_reset( void );
230 
232  void H2_LevelPops( bool &lgPopsConverged, double &old_value, double &new_value );
233 
240  void H2_PunchDo( FILE* io , char chJOB[] , const char chTime[] , long int ipPun );
241 
243  void H2_ParseSave( Parser &p, ostringstream& chHeader );
244 
246  double H2_itrzn( void );
247 
252  void H2_Prt_column_density( FILE *ioMEAN );
253 
254  void set_numLevelsMatrix( long numLevels );
255 
256  void H2_ReadDissocEnergies( void );
257 
260 
266 
267  /* total spontaneous dissociation rate [s-1],
268  * summed over excited electronic states, weighted by pops */
270 
273 
277 
282 
287 
291 
295 
296  double HeatDiss;
297  double HeatDexc;
298  double HeatDexc_old;
301 
303  double Average_A;
314 
318 
322 
324 
328  para_density;
329 
330  // single precision versions of the above
333 
335  double ortho_colden ,
336  para_colden;
337 
338  /* old and older ortho - para ratios, used to determine whether soln is converged */
340 
341  // these remember the largest and smallest factors needed to
342  // renormalize the H2 chemistry
343  double renorm_max ,
344  renorm_min;
345 
346  // this will say how many times the large H2 molecule has been called in this zone -
347  // if not called (due to low H2 abundance) then not need to update its line arrays
348  long int nCall_this_zone;
349 
350  // flag saying whether to bother with the large molecule at all,
351  // default is false, set true with atom h2 on command
352  bool lgEnabled;
353 
354  // this is the number of electronic levels to include in the output - default is 1,
355  // only X. changed with PRINT LINES H2 ELECTRONIC and option on PUNCH H2 LINES commands
357 
361 
364 
367 
368  // include collision rates that come from real calculations,
369  // off with atom h2 collisions off command
371 
373  bool lgLTE;
374 
377 
378  // which set of He - H2 collisions to use? default is ORNL, other
379  // is Le BOURlet
381 
382  // flag saying whether (true) or not to use ORNL H2 - H2 collisions
385 
390 
391  long int loop_h2_oscil;
392  long int nzoneEval;
393 
396 
399 
401 
403  int nTRACE;
404 
408  n_trace_full,
410 
411  // the number of electronic quantum states to include.
412  // To do both Lyman and Werner bands want nelec = 3
413  long int n_elec_states;
414 
415  /* this is fraction of population that is within levels done with matrix */
416  double frac_matrix;
417 
418  /* used to recall the temperature used for last set of Boltzmann factors */
419  double TeUsedBoltz;
420  double TeUsedColl;
421 
422  explicit diatomics( const string& a, const double& e_star, const double* const abund, double (*fun)(double) ) ;
423 
429  vector< diss_tran > Diss_Trans;
430 
431 private:
432  string label;
433  string shortlabel;
434  string path;
435 
438  /* >> chng 05 jul 15, TE, H2g = sum (v=0, J=0,1) */
439  /* >>chng 05 jul 29, to 0.5 eV, this goes up to J=8 for v=0 */
440  /* >>chng 05 aug 03, slight upward change in energy to include the J=8 level,
441  * also give energy in waveumbers for simplicity (save h2 levels give energy in ryd) */
442  /*#define ENERGY_H2_STAR (0.5/EVRYD/WAVNRYD)*/
443  /* energy of v=0, J=8 is 4051.73, J=9 is 5001.97
444  * v=1, J=0 is 4161.14 */
445 public:
446  const double ENERGY_H2_STAR;
447 
448 private:
449  // pointer to the density of the species
450  const double* const dense_total;
451 
453 
454  /* these vars are private for H2 but uses same style as all other header files -
455  * the extern is extern in all except cddefines */
456 
458  long int nEner_H2_ground;
459 
462 
465 
468 
469  //int H2_nRot_add_ortho_para[N_ELEC];
472  long int nVib_hi[N_ELEC];
474  valarray<long> nRot_hi[N_ELEC];
477  long int Jlowest[N_ELEC];
484  vector<CollRateCoeffArray> RateCoefTable;
485 
486  // quantities dealing with chemistry
487 #if 1
488 #endif
489 
490  // these are quantities for each state
491 #if 1
492 
496 
502 #endif
503 
505 
506  // these are quantities for states with X
507 #if 1
527 
530 #endif
531 
532  valarray<realnum> H2_X_source;
533  valarray<realnum> H2_X_sink;
534 
538 
540  double H2_den_s , H2_den_g;
541 
544 
545  valarray<long> ipVib_H2_energy_sort;
546  valarray<long> ipElec_H2_energy_sort;
547  valarray<long> ipRot_H2_energy_sort;
550 
553  long int nXLevelsMatrix;
554  long int ndim_allocated;
556  AulDest,
557  AulPump,
559  vector<double> pops, create, destroy, depart, stat_levn, excit;
560 
561  long int levelAsEval;
562  bool lgFirst;
563  long int nzone_eval;
565 
568 
572 
574  long int nH2_pops;
575  long int nH2_zone;
576 
579 
584 
585  /* LTE cooling per molecule */
586  vector<double> LTE_Temp, LTE_cool;
587 
592 
593 public:
594  /* Read LTE cooling per molecule */
596 
597  /* interpolate over LTE cooling array */
598  double interpolate_LTE_Cooling( double Temp );
599 };
600 
601 /* compute H2 continuum dissociation cross sections */
602 double MolDissocCrossSection( const diss_tran& tran, const double& Mol_Ene );
603 
604 double Yan_H2_CS( double energy_ryd /* photon energy in ryd */);
605 
606 /* compute H2 continuum dissoication opacities */
607 //double MolDissocOpacity( const diss_tran& tran, const double& Mol_Ene );
608 
609 #endif /* H2_PRIV_H_ */
610 
multi_arr< double, 2 > H2_rad_rate_in
Definition: h2_priv.h:511
int nTRACE
Definition: h2_priv.h:403
char chH2ColliderLabels[N_X_COLLIDER][chN_X_COLLIDER]
Definition: h2_priv.h:452
multi_arr< double, 2 > H2_col_rate_out
Definition: h2_priv.h:510
const int N_ELEC
Definition: h2_priv.h:34
multi_arr< realnum, 3 > H2_dissprob
Definition: h2_priv.h:493
const double ENERGY_H2_STAR
Definition: h2_priv.h:446
double Average_A
Definition: h2_priv.h:303
double gs_rate(void)
Definition: mole_h2_etc.cpp:111
realnum GetXColden(long iVib, long iRot)
Definition: mole_h2.cpp:2351
double renorm_min
Definition: h2_priv.h:343
double rel_pop_LTE_s
Definition: h2_priv.h:290
double rel_pop_LTE_g
Definition: h2_priv.h:289
STATIC long int ipPun
Definition: save_do.cpp:368
double H2_DissocEnergies[N_ELEC]
Definition: h2_priv.h:470
void SolveExcitedElectronicLevels(void)
Definition: mole_h2.cpp:1937
void H2_LineZero(void)
Definition: mole_h2.cpp:440
void H2_Punch_line_data(FILE *ioPUN, bool lgDoAll)
Definition: mole_h2_io.cpp:1065
double TeUsedBoltz
Definition: h2_priv.h:419
double Average_collH2_excit
Definition: h2_priv.h:307
long int n_elec_states
Definition: h2_priv.h:413
multi_arr< double, 2 > AulDest
Definition: h2_priv.h:555
bool lgImgMatrix
Definition: h2_priv.h:591
double HeatChangeOld
Definition: h2_priv.h:300
long ip_photo_opac_thresh
Definition: h2_priv.h:320
realnum mass_amu
Definition: h2_priv.h:400
double spon_diss_tot
Definition: h2_priv.h:269
vector< double > LTE_cool
Definition: h2_priv.h:586
double Abund() const
Definition: h2_priv.h:75
void H2_LevelPops(bool &lgPopsConverged, double &old_value, double &new_value)
Definition: mole_h2.cpp:896
int n_trace_full
Definition: h2_priv.h:406
double ortho_para_older
Definition: h2_priv.h:339
void H2_ParseSave(Parser &p, ostringstream &chHeader)
Definition: mole_h2_io.cpp:74
double average_energy_s
Definition: h2_priv.h:294
multi_arr< double, 2 > pops_per_vib
Definition: h2_priv.h:461
t_coll_source coll_source[N_X_COLLIDER]
Definition: h2_priv.h:323
double HeatDexc_deriv
Definition: h2_priv.h:299
multi_arr< realnum, 2 > H2_coll_dissoc_rate_coef_H2
Definition: h2_priv.h:529
double Cont_Dissoc_Rate_H2g
Definition: h2_priv.h:285
Definition: h2_priv.h:44
bool lgPrtMatrix
Definition: h2_priv.h:589
molecule * sp_star
Definition: h2_priv.h:425
bool lgREAD_DATA
Definition: h2_priv.h:259
void H2_RTMake(linefunc line_one)
Definition: mole_h2.cpp:388
multi_arr< double, 2 > CollRate_levn
Definition: h2_priv.h:555
valarray< long > ipVib_H2_energy_sort
Definition: h2_priv.h:545
long int Jlowest[N_ELEC]
Definition: h2_priv.h:477
double ortho_colden
Definition: h2_priv.h:335
realnum H2_CollidRateEvalOne(long iVibHi, long iRotHi, long iVibLo, long iRotLo, long ipHi, long ipLo, long nColl, double temp_K)
Definition: mole_h2_coll.cpp:99
const int nTE_HMINUS
Definition: h2_priv.h:31
double xSTDNoise
Definition: h2_priv.h:395
double H2_to_H_limit
Definition: h2_priv.h:398
bool lgEvaluated
Definition: h2_priv.h:317
void H2_CollidRateRead(long int nColl)
Definition: mole_h2_coll.cpp:163
vector< CollRateCoeffArray > RateCoefTable
Definition: h2_priv.h:484
void H2_Prt_column_density(FILE *ioMEAN)
Definition: mole_h2_io.cpp:369
double rate_grain_op_conserve
Definition: h2_priv.h:280
vector< double > LTE_Temp
Definition: h2_priv.h:586
int nElecLevelOutput
Definition: h2_priv.h:356
double H2_den_g
Definition: h2_priv.h:540
multi_arr< realnum, 2 > H2_coll_dissoc_rate_coef
Definition: h2_priv.h:526
double HeatChange
Definition: h2_priv.h:300
Definition: h2_priv.h:72
multi_arr< double, 2 > AulEscp
Definition: h2_priv.h:555
void H2_Read_hminus_distribution(void)
Definition: mole_h2_io.cpp:942
TransitionList trans
Definition: h2_priv.h:427
Definition: mole.h:145
Definition: parser.h:41
double Average_collH2_dissoc_g
Definition: h2_priv.h:312
void GetIndices(long &ipHi, long &ipLo, const char *chLine, long &i) const
Definition: mole_h2_coll.cpp:202
long int nEner_H2_ground
Definition: h2_priv.h:458
void H2_X_sink_and_source(void)
Definition: mole_h2.cpp:52
long int ndim_allocated
Definition: h2_priv.h:554
int n_trace_matrix
Definition: h2_priv.h:406
double Average_collH2_deexcit
Definition: h2_priv.h:305
multi_arr< realnum, 3 > CollRateErrFac
Definition: h2_priv.h:483
Definition: h2_priv.h:39
multi_arr< realnum, 2 > H2_X_formation
Definition: h2_priv.h:514
const realnum H2_logte_hminus[nTE_HMINUS]
Definition: h2_priv.h:37
double Yan_H2_CS(double energy_ryd)
Definition: mole_h2_etc.cpp:358
double Solomon_dissoc_rate_g
Definition: h2_priv.h:271
double TeUsedColl
Definition: h2_priv.h:420
Definition: h2_priv.h:61
multi_arr< realnum, 3 > CollRateCoeff
Definition: h2_priv.h:482
const double *const dense_total
Definition: h2_priv.h:450
valarray< realnum > H2_X_sink
Definition: h2_priv.h:533
long int levelAsEval
Definition: h2_priv.h:561
vector< diss_tran > Diss_Trans
Definition: h2_priv.h:429
long int iteration_evaluated
Definition: h2_priv.h:564
multi_arr< double, 3 > Cont_Dissoc_Rate
Definition: h2_priv.h:286
long n
Definition: h2_priv.h:41
multi_arr< long int, 3 > ipEnergySort
Definition: h2_priv.h:548
multi_arr< double, 3 > H2_old_populations
Definition: h2_priv.h:498
double LTE_Cooling_per_H2()
Definition: mole_h2_etc.cpp:438
double GetExcitedElecDensity(void)
Definition: mole_h2.cpp:2530
bool lgH2_ortho_para_coll_on
Definition: h2_priv.h:376
molecule * sp
Definition: h2_priv.h:424
double Average_collH2_dissoc_s
Definition: h2_priv.h:313
void H2_CollidRateEvalAll(void)
Definition: mole_h2_coll.cpp:13
long int iterationAsEval
Definition: h2_priv.h:504
multi_arr< realnum, 3 > H2_disske
Definition: h2_priv.h:494
double ortho_density
Definition: h2_priv.h:326
diss_tran(diss_level a, diss_level b)
Definition: h2_priv.h:47
double HeatDexc_old
Definition: h2_priv.h:298
int n_trace_iterations
Definition: h2_priv.h:406
t_abund abund
Definition: abund.cpp:5
realnum para_density_f
Definition: h2_priv.h:331
void H2_ReadEnergies()
Definition: mole_h2_io.cpp:634
void mole_H2_LTE(void)
Definition: mole_h2_etc.cpp:225
long OpacityCreate(vector< double > &stack)
Definition: mole_h2_etc.cpp:167
void H2_RT_diffuse(void)
Definition: mole_h2.cpp:369
void mole_H2_form(void)
Definition: mole_h2_form.cpp:14
double para_density
Definition: h2_priv.h:326
vector< double > energies
Definition: h2_priv.h:56
Definition: quantumstate.h:35
void H2_RT_tau_reset(void)
Definition: mole_h2.cpp:456
multi_arr< realnum, 3 > H2_X_hminus_formation_distribution
Definition: h2_priv.h:543
double ortho_para_old
Definition: h2_priv.h:339
bool lgLTE
Definition: h2_priv.h:373
double photodissoc_BigH2_H2s
Definition: h2_priv.h:264
bool lgColl_dissoc_coll
Definition: h2_priv.h:366
double energy(const genericState &gs)
Definition: generic_state.cpp:51
string filename
Definition: h2_priv.h:69
long int nzoneAsEval
Definition: h2_priv.h:504
long int nLevels_per_elec[N_ELEC]
Definition: h2_priv.h:479
double xMeanNoise
Definition: h2_priv.h:395
vector< double > destroy
Definition: h2_priv.h:559
bool lgEnabled
Definition: h2_priv.h:352
double H2_den_s
Definition: h2_priv.h:540
Definition: proxy_iterator.h:58
long int nzoneEval
Definition: h2_priv.h:392
multi_arr< double, 3 > H2_populations_LTE
Definition: h2_priv.h:499
string label
Definition: h2_priv.h:432
long int nCall_this_zone
Definition: h2_priv.h:348
string path
Definition: h2_priv.h:434
void Mol_Photo_Diss_Rates(void)
Definition: mole_dissociate.cpp:131
void H2_RT_tau_inc(void)
Definition: mole_h2.cpp:410
void(* linefunc)(const TransitionProxy &t, bool lgShield_this_zone, realnum pestrk, realnum DopplerWidth)
Definition: h2_priv.h:10
void H2_Prt_Zone(void)
Definition: mole_h2_io.cpp:307
bool lgH2_PAH2_ORNL
Definition: h2_priv.h:384
bool lgH2_NOISE
Definition: h2_priv.h:387
double Average_collH_deexcit
Definition: h2_priv.h:306
const int N_X_COLLIDER
Definition: h2_priv.h:20
bool lgFirst
Definition: h2_priv.h:562
long int nCall_this_iteration
Definition: h2_priv.h:583
Definition: atmdat.h:12
double H2_InterEnergy(void)
void H2_RT_OTS(void)
Definition: mole_h2.cpp:2421
long magic
Definition: h2_priv.h:67
float realnum
Definition: cddefines.h:127
multi_arr< realnum, 2 > H2_X_colden_LTE
Definition: h2_priv.h:520
double photodissoc_BigH2_H2g
Definition: h2_priv.h:265
multi_arr< bool, 2 > lgH2_radiative
Definition: h2_priv.h:571
long int nH2_pops
Definition: h2_priv.h:574
long int nXLevelsMatrix
Definition: h2_priv.h:553
double MolDissocOpacity(const diss_tran &tran, const double &Mol_Ene)
Definition: mole_dissociate.cpp:95
bool lgH2_He_ORNL
Definition: h2_priv.h:380
multi_arr< realnum, 2 > H2_X_Hmin_back
Definition: h2_priv.h:516
void H2_PrtDepartCoef(void)
Definition: mole_h2_io.cpp:336
void H2_LinesAdd(void)
Definition: mole_h2_io.cpp:44
double(* photoion_opacity_fun)(double energy)
Definition: h2_priv.h:81
void H2_Cooling(void)
Definition: mole_h2.cpp:2185
long ip_photo_opac_offset
Definition: h2_priv.h:321
multi_arr< realnum, 6 > H2_SaveLine
Definition: h2_priv.h:567
long v
Definition: h2_priv.h:41
multi_arr< double, 3 > H2_rad_rate_out
Definition: h2_priv.h:495
bool lgColl_deexec_Calc
Definition: h2_priv.h:363
double average_energy_g
Definition: h2_priv.h:293
double H2_DR(void)
Definition: mole_h2.cpp:2415
valarray< long > ipElec_H2_energy_sort
Definition: h2_priv.h:546
void H2_ContPoint(void)
Definition: mole_h2.cpp:277
long int loop_h2_oscil
Definition: h2_priv.h:391
void Read_Mol_Diss_cross_sections(void)
Definition: mole_dissociate.cpp:12
multi_arr< realnum, 3 > H2_X_grain_formation_distribution
Definition: h2_priv.h:537
realnum ortho_density_f
Definition: h2_priv.h:331
double rate_grain_J1_to_J0
Definition: h2_priv.h:281
double GetHeatRate(const diss_tran &tran)
Definition: mole_dissociate.cpp:227
double HeatDexc
Definition: h2_priv.h:297
void H2_Calc_Average_Rates(void)
Definition: mole_h2.cpp:2445
double ortho_para_current
Definition: h2_priv.h:339
double frac_matrix
Definition: h2_priv.h:416
long int nzone_eval
Definition: h2_priv.h:563
void H2_Solomon_rate(void)
Definition: mole_h2_etc.cpp:24
diatomics(const string &a, const double &e_star, const double *const abund, double(*fun)(double))
Definition: h2.cpp:13
double HeatDiss
Definition: h2_priv.h:296
vector< double > stat_levn
Definition: h2_priv.h:559
void H2_ReadTransprob(long int nelec, TransitionList &trans)
Definition: mole_h2_io.cpp:404
double photo_heat_soft
Definition: h2_priv.h:262
bool lgH2_grain_deexcitation
Definition: h2_priv.h:370
multi_arr< double, 2 > H2_X_rate_to_elec_excited
Definition: h2_priv.h:524
double photo_heat_hard
Definition: h2_priv.h:263
string shortlabel
Definition: h2_priv.h:433
void H2_zero_pops_too_low(void)
Definition: mole_h2_etc.cpp:184
multi_arr< double, 2 > H2_X_rate_from_elec_excited
Definition: h2_priv.h:522
long int getLine(long iElecHi, long iVibHi, long iRotHi, long iElecLo, long iVibLo, long iRotLo, double *relint, double *absint)
Definition: mole_h2_io.cpp:1865
double renorm_max
Definition: h2_priv.h:343
qList states
Definition: h2_priv.h:426
void H2_Colden(const char *chLabel)
Definition: mole_h2.cpp:2370
double Average_collH_excit
Definition: h2_priv.h:308
double H2_RadPress(void)
Definition: mole_h2.cpp:318
void H2_X_coll_rate_evaluate(void)
Definition: mole_h2.cpp:199
void H2_PunchLineStuff(FILE *io, realnum xLimit, long index)
Definition: mole_h2_io.cpp:1106
void H2_Read_LTE_cooling_per_H2()
Definition: mole_h2_io.cpp:1934
double Solomon_dissoc_rate_s
Definition: h2_priv.h:272
double photoionize_rate
Definition: h2_priv.h:261
bool lgH2_ORH2_ORNL
Definition: h2_priv.h:383
void H2_Level_low_matrix(realnum abundance)
Definition: mole_h2.cpp:472
multi_arr< double, 2 > AulPump
Definition: h2_priv.h:555
multi_arr< long int, 2 > ipTransitionSort
Definition: h2_priv.h:549
void H2_PunchDo(FILE *io, char chJOB[], const char chTime[], long int ipPun)
Definition: mole_h2_io.cpp:1143
double Cont_Dissoc_Rate_H2s
Definition: h2_priv.h:284
bool lgColl_gbar
Definition: h2_priv.h:360
bool lgH2_NOISECOSMIC
Definition: h2_priv.h:389
Definition: transition.h:297
double Solomon_elec_decay_g
Definition: h2_priv.h:275
vector< double > pops
Definition: h2_priv.h:559
void H2_Reset(void)
Definition: mole_h2_etc.cpp:298
TransitionList::iterator rad_end
Definition: h2_priv.h:428
void H2_ReadDissocEnergies(void)
Definition: mole_h2_io.cpp:793
valarray< long > nRot_hi[N_ELEC]
Definition: h2_priv.h:474
double pops_per_elec[N_ELEC]
Definition: h2_priv.h:481
Definition: transition.h:23
double Solomon_elec_decay_s
Definition: h2_priv.h:276
double para_colden
Definition: h2_priv.h:335
double H2_renorm_chemistry
Definition: h2_priv.h:464
vector< double > create
Definition: h2_priv.h:559
double Average_collH_dissoc_g
Definition: h2_priv.h:310
t_coll_source()
Definition: h2_priv.h:63
valarray< long > ipRot_H2_energy_sort
Definition: h2_priv.h:547
valarray< realnum > H2_X_source
Definition: h2_priv.h:532
double H2_Accel(void)
Definition: mole_h2.cpp:295
double rate_coeff
Definition: h2_priv.h:58
long int nzone_nlevel_set
Definition: h2_priv.h:578
long j
Definition: h2_priv.h:41
double Average_collH_dissoc_s
Definition: h2_priv.h:311
void CalcPhotoionizationRate(void)
Definition: mole_h2_etc.cpp:407
void init(void)
Definition: mole_h2_create.cpp:107
double interpolate_LTE_Cooling(double Temp)
Definition: mole_h2_etc.cpp:497
void set_numLevelsMatrix(long numLevels)
Definition: mole_h2_io.cpp:1927
void H2_ReadDissprob(long int nelec)
Definition: mole_h2_io.cpp:854
double MolDissocCrossSection(const diss_tran &tran, const double &Mol_Ene)
Definition: mole_dissociate.cpp:104
double GetDissociationRate(const diss_tran &tran)
Definition: mole_dissociate.cpp:194
vector< double > depart
Definition: h2_priv.h:559
vector< double > excit
Definition: h2_priv.h:559
multi_arr< realnum, 2 > H2_X_coll_rate
Definition: h2_priv.h:467
double H2_itrzn(void)
Definition: mole_h2.cpp:264
long int nVib_hi[N_ELEC]
Definition: h2_priv.h:472
vector< double > xsections
Definition: h2_priv.h:57
void SolveSomeGroundElectronicLevels(void)
Definition: mole_h2.cpp:2048
multi_arr< double, 2 > H2_col_rate_in
Definition: h2_priv.h:509
const int chN_X_COLLIDER
Definition: h2_priv.h:22
int n_trace_final
Definition: h2_priv.h:406
multi_arr< realnum, 2 > H2_X_colden
Definition: h2_priv.h:518
long int nH2_zone
Definition: h2_priv.h:575
multi_arr< int, 2 > H2_ipPhoto
Definition: h2_priv.h:508
diss_level initial
Definition: h2_priv.h:54
double Cont_Diss_Heat_Rate(void)
Definition: mole_dissociate.cpp:206
multi_arr< bool, 3 > H2_lgOrtho
Definition: h2_priv.h:501