Line data Source code
1 : module steady_state_tei
2 :
3 : use shr_kind_mod, only: r8 => shr_kind_r8 ! Real kind to declare variables
4 : use physics_buffer, only: pbuf_get_index, & !
5 : physics_buffer_desc, & !
6 : pbuf_get_field, & ! Needed to access physics buffer
7 : pbuf_set_field
8 : use physics_types, only: physics_state ! Structures containing physics state variables
9 : use ppgrid, only: pcols, pver, pverp, begchunk, endchunk
10 : use cam_abortutils, only: endrun
11 :
12 : implicit none
13 :
14 : private
15 :
16 : !------------------------
17 : ! PUBLIC: interfaces
18 : !------------------------
19 : public :: steady_state_tei_init
20 : public :: steady_state_tei_tend
21 :
22 : integer :: indxTe = -1
23 : integer :: indxTi = -1
24 : integer :: indxQt = -1
25 : integer :: indxAR = -1
26 : integer :: indxO1 = -1
27 : integer :: indxO2 = -1
28 :
29 : real(r8) :: op_mass
30 :
31 : contains
32 :
33 : !==============================================================================
34 :
35 0 : subroutine steady_state_tei_init(pbuf2d)
36 :
37 : !-----------------------------------------------------------------------
38 : ! Time independent initialization for ionosphere simulation.
39 : !-----------------------------------------------------------------------
40 :
41 : use constituents, only: cnst_get_ind
42 : use mo_chem_utls, only: get_spc_ndx ! Routine to get index of adv_mass array for short lived species
43 : use chem_mods, only: adv_mass ! Array holding mass values for short lived species
44 : use infnan, only: nan, assignment(=)
45 :
46 : type(physics_buffer_desc), pointer :: pbuf2d(:,:)
47 : real(r8) :: nanval
48 :
49 0 : call cnst_get_ind( 'O2', indxO2, abort=.true. )
50 0 : call cnst_get_ind( 'O', indxO1, abort=.true. )
51 :
52 : !-------------------------------------------------------------------------------------------------------------------
53 : ! Get physics buffer indices
54 : !-------------------------------------------------------------------------------------------------------------------
55 0 : indxTe = pbuf_get_index( 'TElec' )
56 0 : indxTi = pbuf_get_index( 'TIon' )
57 0 : indxQt = pbuf_get_index( 'QTeAur' )
58 0 : indxAR = pbuf_get_index( 'AurIPRateSum' )
59 :
60 0 : nanval=nan
61 0 : call pbuf_set_field(pbuf2d, indxQt, nanval)
62 0 : call pbuf_set_field(pbuf2d, indxAR, nanval)
63 :
64 0 : op_mass = adv_mass(get_spc_ndx('Op'))
65 :
66 0 : end subroutine steady_state_tei_init
67 :
68 :
69 : !==============================================================================
70 :
71 0 : subroutine steady_state_tei_tend(state,istate, dse_tend, pbuf)
72 : use tei_mod, only: settei
73 : use air_composition, only: mbarv ! Constituent dependent mbar
74 : use solar_parms_data, only: f107=>solar_parms_f107 ! 10.7 cm solar flux
75 : use mo_apex, only: alatm
76 : use perf_mod, only: t_startf, t_stopf ! timing utils
77 : use ionos_state_mod, only: ionos_state
78 : !-------------------------------------------------------------------------------------
79 : ! Calculate dry static energy and O+ tendency for extended ionosphere simulation
80 : !-------------------------------------------------------------------------------------
81 :
82 : !------------------------------Arguments--------------------------------
83 :
84 : type(physics_state), intent(in) :: state ! physics state structure
85 : type(ionos_state),target,intent(in) :: istate ! ionosphere state structure
86 : real(r8), intent(inout) :: dse_tend(:,:) ! dry static energy tendency
87 : type(physics_buffer_desc),pointer :: pbuf(:) ! physics buffer
88 :
89 : !---------------------------Local variables-------------------------------
90 :
91 : real(r8),dimension(pver,state%ncol) :: &
92 0 : tn, & ! neutral temperature (deg K)
93 0 : o2, & ! molecular oxygen (mmr)
94 0 : o1, & ! atomic oxygen (mmr)
95 0 : n2, & ! molecular nitrogen (mmr)
96 0 : ne, & ! electron density (cm3)
97 0 : te, & ! electron temperature (from previous time step) (K)
98 0 : ti, & ! ion temperature (from previous time step) (K)
99 0 : op, & ! O+ number dens (/cm^3)
100 0 : o2p, & ! O2+ number dens (/cm^3)
101 0 : nop, & ! NO+ number dens (/cm^3)
102 0 : barm, & ! mean molecular weight (g/mole)
103 0 : qji_ti ! joule heating from qjoule_ti (used ui,vi) ev/g/s
104 :
105 0 : real(r8) :: chi(state%ncol) ! solar zenith angle (radians)
106 0 : real(r8) :: qtot(pver,state%ncol) ! total ionization rate (s-1 cm-3)
107 0 : real(r8) :: pmid(pver,state%ncol) ! mid-level press (dyne/cm^2) ! 10._r8*pmid(:ncol,k)
108 0 : real(r8) :: pint(pverp,state%ncol) ! interface press (dyne/cm^2)
109 :
110 0 : real(r8) :: te_out(pver,state%ncol) ! output electron temperature (deg K)
111 0 : real(r8) :: ti_out(pver,state%ncol) ! output ion temperature (deg K)
112 0 : real(r8) :: qtotal_out(pver,state%ncol) ! heating rate of neutrals
113 :
114 : integer :: lchnk ! Chunk number
115 : integer :: ncol ! Number of columns in chunk
116 : integer :: i, k
117 :
118 : real(r8), parameter :: evergs = 1.602e-12_r8 ! 1 eV = 1.602e-12 ergs (ergs/eV)
119 : real(r8), parameter :: sToQConv = 6.24E15_r8 ! Conversion from J/kg/s to ev/g/s
120 : real(r8), parameter :: n2_mass = 28._r8 ! N2 molecular weight kg/kmol
121 :
122 0 : real(r8), pointer :: te_ptr(:,:) ! Pointer to electron temperature in pbuf (K)
123 0 : real(r8), pointer :: ti_ptr(:,:) ! Pointer to ion temperature in pbuf (K)
124 0 : real(r8), pointer :: qteaur(:)
125 0 : real(r8), pointer :: aurIPRateSum(:,:) ! Auroral ion production sum for O2+,O+,N2+
126 : ! (s-1 cm-3 from module mo_aurora)
127 :
128 0 : call t_startf ('steady_state_tei_tend')
129 :
130 : !----------------------------------------------------------------
131 : ! Get number of this chunk
132 : !----------------------------------------------------------------
133 0 : lchnk = state%lchnk
134 0 : ncol = state%ncol
135 :
136 0 : call pbuf_get_field(pbuf, indxTe, te_ptr)
137 0 : call pbuf_get_field(pbuf, indxTi, ti_ptr)
138 0 : call pbuf_get_field(pbuf, indxQt, qteaur)
139 0 : call pbuf_get_field(pbuf, indxAR, aurIPRateSum)
140 :
141 0 : do i =1,ncol
142 0 : chi(i) = acos(istate%cosZenAngR(i))
143 0 : qji_ti(1:pver,i) = dse_tend(i,pver:1:-1)*sToQConv*evergs ! J/kg/s --> ev/g/s --> ergs/s/g
144 : ! take the part going into heating the neutrals
145 0 : qji_ti(1:pver,i) = qji_ti(1:pver,i) * mbarv(i,pver:1:-1,lchnk) / (mbarv(i,pver:1:-1,lchnk)+op_mass)
146 0 : qtot(1:pver,i) = istate%sumIonPRates(i,pver:1:-1) + aurIPRateSum(i,pver:1:-1)
147 0 : te(1:pver,i) = te_ptr(i,pver:1:-1)
148 0 : ti(1:pver,i) = ti_ptr(i,pver:1:-1)
149 0 : tn(1:pver,i) = state%t(i,pver:1:-1)
150 0 : pmid(1:pver,i) = state%pmid(i,pver:1:-1) *10._r8 ! dynes/cm2
151 0 : pint(1:pverp,i) = state%pint(i,pverp:1:-1)*10._r8
152 0 : o2(1:pver,i) = state%q(i,pver:1:-1,indxO2)
153 0 : o1(1:pver,i) = state%q(i,pver:1:-1,indxO1)
154 0 : n2(1:pver,i) = istate%n2_mmr(i,pver:1:-1)
155 0 : barm(1:pver,i) = mbarv(i,pver:1:-1,lchnk)
156 0 : ne(1:pver,i) = istate%ndensE(i,pver:1:-1)
157 0 : op(1:pver,i) = istate%ndensOp(i,pver:1:-1)
158 0 : o2p(1:pver,i) = istate%ndensO2p(i,pver:1:-1)
159 0 : nop(1:pver,i) = istate%ndensNOp(i,pver:1:-1)
160 : end do
161 :
162 : call settei(tn,o2,o1,n2,ne,te,ti,op,o2p,nop,barm,qji_ti, &
163 0 : f107, chi, qtot, qteaur(:ncol), alatm(:ncol,lchnk), &
164 0 : istate%dipmag(:ncol,1), pmid, pint, 1,pver,ncol, &
165 0 : te_out, ti_out, qtotal_out )
166 :
167 0 : do i =1,ncol
168 0 : te_ptr(i,pver:1:-1) = te_out(1:pver,i)
169 0 : ti_ptr(i,pver:1:-1) = ti_out(1:pver,i)
170 0 : dse_tend(i,pver:1:-1) = qtotal_out(1:pver,i) / evergs / sToQConv ! ergs/g/s --> eV/g/s --> J/kg/sec
171 : end do
172 :
173 0 : call t_stopf ('steady_state_tei_tend')
174 :
175 0 : end subroutine steady_state_tei_tend
176 :
177 : end module steady_state_tei
|