Skip to content

Commit 7b20346

Browse files
committed
BUG: inelastic collisions, etc. bugfixes in cet.
now it compiles
1 parent 293fb6c commit 7b20346

File tree

1 file changed

+18
-10
lines changed

1 file changed

+18
-10
lines changed

src/calc_electron_temperature.cpp

Lines changed: 18 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -95,6 +95,10 @@ void Ions::calc_electron_temperature(Neutrals neutrals, Grid grid) {
9595
std::vector<arma_cube> Qenc;
9696
arma_cube Qencm, Qencp, Qenc_v;
9797

98+
// Inelastic electron-neutral collisions:
99+
std::vector<arma_cube> Qenc_inelastic;
100+
arma_cube Qrotm, Qrotp, Qf, Qexc, Qvib_O2, Qvib_N2;
101+
98102
// Initialize everything to zero!
99103

100104
epsilon.set_size(grid.get_nLons(), grid.get_nLats(), grid.get_nAlts());
@@ -152,10 +156,13 @@ void Ions::calc_electron_temperature(Neutrals neutrals, Grid grid) {
152156

153157
// electron-neutral inelastic collisions
154158
if (input.get_do_electron_neutral_inelastic_collisional_heating()) {
155-
Qenc = calc_electron_neutral_inelastic_collisions(*this, neutrals);
156-
Qencp = Qenc[0];
157-
Qencm = Qenc[1];
158-
Qenc_v = Qenc[2]; // Friction
159+
Qenc_inelastic = calc_electron_neutral_inelastic_collisions(*this, neutrals);
160+
Qrotm = Qenc_inelastic[0];
161+
Qrotp = Qenc_inelastic[1];
162+
Qf = Qenc_inelastic[2];
163+
Qexc = Qenc_inelastic[3];
164+
Qvib_O2 = Qenc_inelastic[4];
165+
Qvib_N2 = Qenc_inelastic[5];
159166
}
160167

161168
electron_temperature_scgc = neutrals.temperature_scgc;
@@ -485,11 +492,11 @@ std::vector<arma_cube> calc_electron_neutral_inelastic_collisions(Ions &ions, Ne
485492
// GITM clamped e- temp before calculating, but that makles things really hard here.
486493
// Instead we will use all Te's, and then limit the outputs after.
487494

488-
arma_cube logQ = (5.0148e-31*pow(Te,9) - 1.5346e-26*pow(Te,8)
489-
+ 2.0127e-22*pow(Te,7) - 1.4791e-18*pow(Te,6)
490-
+ 6.6865e-15*pow(Te,5) - 1.9228e-11*pow(Te,4)
491-
+ 3.5187e-8*pow(Te,3) - 3.996e-5*pow(Te,2)
492-
+ 0.0267*Te - 19.9171);
495+
logQ = (5.0148e-31*pow(Te,9) - 1.5346e-26*pow(Te,8)
496+
+ 2.0127e-22*pow(Te,7) - 1.4791e-18*pow(Te,6)
497+
+ 6.6865e-15*pow(Te,5) - 1.9228e-11*pow(Te,4)
498+
+ 3.5187e-8*pow(Te,3) - 3.996e-5*pow(Te,2)
499+
+ 0.0267*Te - 19.9171);
493500
// GITM's Te_6000 was from 300 - 6000, which corresponds to ~-15.9 & 198.3 for logQ
494501
// Mask the values outside of this range...
495502
// TODO: Should we do it this way or just use the clamp?
@@ -565,5 +572,6 @@ std::vector<arma_cube> calc_electron_neutral_inelastic_collisions(Ions &ions, Ne
565572
}
566573

567574
Qvib_N2 = -ne % nn2 * 1.e-12 % Qvib_N2 * 1.6e-13;
568-
575+
576+
return std::vector<arma_cube> {Qrotm, Qrotp, Qf, Qexc, Qvib_O2, Qvib_N2};
569577
}

0 commit comments

Comments
 (0)