CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_rrtmgp_pre.F90
1
4
6 use machine, only: kind_phys
7 use funcphys, only: fpvs
10 use module_ozphys, only: ty_ozphys
11 use mo_gas_concentrations, only: ty_gas_concs
12 use radiation_tools, only: check_error_msg,cmp_tlev
13 use rrtmgp_lw_gas_optics, only: lw_gas_props
14
15 implicit none
16
17 real(kind_phys), parameter :: &
18 amd = 28.9644_kind_phys, &
19 amw = 18.0154_kind_phys, &
20 amo3 = 47.9982_kind_phys, &
21 amdw = amd/amw, &
22 amdo3 = amd/amo3
23 real(kind_phys), parameter :: eps = 1.0e-6_kind_phys
24 real(kind_phys), parameter :: oneminus = 1.0_kind_phys - eps
25 real(kind_phys), parameter :: ftiny = 1.0e-12_kind_phys
26
27 ! Save trace gas indices.
28 integer :: istr_h2o, istr_co2, istr_o3, istr_n2o, istr_ch4, istr_o2, istr_ccl4, &
29 istr_cfc11, istr_cfc12, istr_cfc22
30
31 public gfs_rrtmgp_pre_run,gfs_rrtmgp_pre_init
32contains
33
37 subroutine gfs_rrtmgp_pre_init(nGases, active_gases, active_gases_array, errmsg, errflg)
38 ! Inputs
39 integer, intent(in) :: &
40 ngases
41 character(len=*), intent(in) :: &
42 active_gases
43 character(len=*), dimension(:), intent(out), optional :: &
44 active_gases_array
45
46 ! Outputs
47 character(len=*), intent(out) :: &
48 errmsg
49 integer, intent(out) :: &
50 errflg
51
52 ! Local variables
53 character(len=1) :: tempstr
54 integer :: ij, count
55 integer,dimension(nGases,2) :: gasindices
56
57 ! Initialize
58 errmsg = ''
59 errflg = 0
60
61 if (len(active_gases) .eq. 0) return
62
63 ! Which gases are active? Provided via physics namelist.
64
65 ! Pull out gas names from list...
66 ! First grab indices in character array corresponding to start:end of gas name.
67 gasindices(1,1)=1
68 count=1
69 do ij=1,len(active_gases)
70 tempstr=trim(active_gases(ij:ij))
71 if (tempstr .eq. '_') then
72 gasindices(count,2)=ij-1
73 gasindices(count+1,1)=ij+1
74 count=count+1
75 endif
76 enddo
77 gasindices(ngases,2)=len(trim(active_gases))
78
79 ! Now extract the gas names
80 do ij=1,ngases
81 active_gases_array(ij) = active_gases(gasindices(ij,1):gasindices(ij,2))
82 if(trim(active_gases_array(ij)) .eq. 'h2o') istr_h2o = ij
83 if(trim(active_gases_array(ij)) .eq. 'co2') istr_co2 = ij
84 if(trim(active_gases_array(ij)) .eq. 'o3') istr_o3 = ij
85 if(trim(active_gases_array(ij)) .eq. 'n2o') istr_n2o = ij
86 if(trim(active_gases_array(ij)) .eq. 'ch4') istr_ch4 = ij
87 if(trim(active_gases_array(ij)) .eq. 'o2') istr_o2 = ij
88 if(trim(active_gases_array(ij)) .eq. 'ccl4') istr_ccl4 = ij
89 if(trim(active_gases_array(ij)) .eq. 'cfc11') istr_cfc11 = ij
90 if(trim(active_gases_array(ij)) .eq. 'cfc12') istr_cfc12 = ij
91 if(trim(active_gases_array(ij)) .eq. 'cfc22') istr_cfc22 = ij
92 enddo
93
94 end subroutine gfs_rrtmgp_pre_init
95
99 subroutine gfs_rrtmgp_pre_run(me, nCol, nLev, i_o3, doSWrad, doLWrad, fhswr, fhlwr, &
100 xlat, xlon, prsl, tgrs, prslk, prsi, qgrs, tsfc, coslat, sinlat, con_g, con_rd, &
101 con_eps, con_epsm1, con_fvirt, con_epsqs, solhr, raddt, p_lay, t_lay, p_lev, t_lev, &
102 vmr_o2, vmr_h2o, vmr_o3, vmr_ch4, &
103 vmr_n2o, vmr_co2, tsfg, tsfa, qs_lay, q_lay, tv_lay, &
104 relhum, deltaZ, deltaZc, deltaP, active_gases_array, &
105 tsfc_radtime, coszen, coszdg, top_at_1, iSFC, iTOA, nDay, idxday, semis, &
106 sfc_emiss_byband, ico2, ozphys, con_pi, errmsg, errflg)
107
108 ! Inputs
109 integer, intent(in) :: &
110 me, & !< MPI rank
111 ncol, & !< Number of horizontal grid points
112 nlev, & !< Number of vertical layers
113 ico2, & !< Flag for co2 radiation scheme
114 i_o3
115 type(ty_ozphys),intent(in) :: &
116 ozphys
117 logical, intent(in) :: &
118 doswrad, & !< Call SW radiation?
119 dolwrad
120 real(kind_phys), intent(in) :: &
121 fhswr, & !< Frequency of SW radiation call.
122 fhlwr
123 real(kind_phys), intent(in) :: &
124 con_g, & !< Physical constant: gravitational constant
125 con_rd, & !< Physical constant: gas-constant for dry air
126 con_eps, & !< Physical constant: Epsilon (Rd/Rv)
127 con_epsm1, & !< Physical constant: Epsilon (Rd/Rv) minus one
128 con_fvirt, & !< Physical constant: Inverse of epsilon minus one
129 con_epsqs, & !< Physical constant: Minimum saturation mixing-ratio (kg/kg)
130 con_pi, & !< Physical constant: Pi
131 solhr
132 real(kind_phys), dimension(:), intent(in) :: &
133 xlon, & !< Longitude
134 xlat, & !< Latitude
135 tsfc, & !< Surface skin temperature (K)
136 coslat, & !< Cosine(latitude)
137 sinlat, & !< Sine(latitude)
138 semis
139 real(kind_phys), dimension(:,:), intent(in) :: &
140 prsl, & !< Pressure at model-layer centers (Pa)
141 tgrs, & !< Temperature at model-layer centers (K)
142 prslk, & !< Exner function at model layer centers (1)
143 prsi
144 real(kind_phys), dimension(:,:,:), intent(in) :: &
145 qgrs
146 character(len=*), dimension(:), intent(in), optional :: &
147 active_gases_array
148
149 ! Outputs
150 character(len=*), intent(out) :: &
151 errmsg
152 integer, intent(out) :: &
153 errflg, & !< Error flag
154 nday
155 integer, intent(inout) :: &
156 isfc, & !< Vertical index for surface
157 itoa
158 logical, intent(inout) :: &
159 top_at_1
160 real(kind_phys), intent(inout) :: &
161 raddt
162 real(kind_phys), dimension(:), intent(inout) :: &
163 tsfg, & !< Ground temperature
164 tsfa, & !< Skin temperature
165 tsfc_radtime, & !< Surface temperature at radiation timestep
166 coszen, & !< Cosine of SZA
167 coszdg
168 integer, dimension(:), intent(inout) :: &
169 idxday
170 real(kind_phys), dimension(:,:), intent(inout), optional :: &
171 p_lay, & !< Pressure at model-layer
172 t_lay, & !< Temperature at model layer
173 q_lay, & !< Water-vapor mixing ratio (kg/kg)
174 tv_lay, & !< Virtual temperature at model-layers
175 relhum, & !< Relative-humidity at model-layers
176 qs_lay, & !< Saturation vapor pressure at model-layers
177 deltaz, & !< Layer thickness (m)
178 deltazc, & !< Layer thickness (m) (between layer centers)
179 deltap, & !< Layer thickness (Pa)
180 p_lev, & !< Pressure at model-interface
181 sfc_emiss_byband, & !<
182 t_lev, & !< Temperature at model-interface
183 vmr_o2, vmr_h2o, vmr_o3, vmr_ch4, vmr_n2o, vmr_co2
184
185 ! Local variables
186 integer :: i, j, icol, iband, ilay, ilev, isfc_ilev
187 real(kind_phys) :: es, tem1, tem2, pfac
188 real(kind_phys), dimension(nLev+1) :: hgtb
189 real(kind_phys), dimension(nLev) :: hgtc
190 real(kind_phys), dimension(nCol,nLev) :: o3_lay
191 real(kind_phys), dimension(nCol,nLev, NF_VGAS) :: gas_vmr
192 real(kind_phys) :: con_rdog
193
194 ! Initialize CCPP error handling variables
195 errmsg = ''
196 errflg = 0
197
198 nday = 0
199 idxday = 0
200 if (.not. (doswrad .or. dolwrad)) return
201
202 ! #######################################################################################
203 ! What is vertical ordering?
204 ! #######################################################################################
205 top_at_1 = (prsi(1,1) .lt. prsi(1, nlev))
206 if (top_at_1) then
207 isfc = nlev
208 itoa = 1
209 isfc_ilev = isfc + 1
210 else
211 isfc = 1
212 itoa = nlev
213 isfc_ilev = 1
214 endif
215
216 ! #######################################################################################
217 ! Compute some fields needed by RRTMGP
218 ! #######################################################################################
219
220 ! Water-vapor mixing-ratio
221 q_lay(1:ncol,:) = qgrs(1:ncol,:,1)
222 where(q_lay .lt. 1.e-6) q_lay = 1.e-6
223
224 ! Pressure at layer-interface
225 p_lev(1:ncol,:) = prsi(1:ncol,:)
226
227 ! Pressure at layer-center
228 p_lay(1:ncol,:) = prsl(1:ncol,:)
229
230 ! Temperature at layer-center
231 t_lay(1:ncol,:) = tgrs(1:ncol,:)
232
233 ! Bound temperature/pressure at layer centers.
234 do ilay=1,nlev
235 do icol=1,ncol
236 if (t_lay(icol,ilay) .le. lw_gas_props%get_temp_min()) then
237 t_lay(icol,ilay) = lw_gas_props%get_temp_min() + epsilon(lw_gas_props%get_temp_min())
238 endif
239 if (p_lay(icol,ilay) .le. lw_gas_props%get_press_min()) then
240 p_lay(icol,ilay) = lw_gas_props%get_press_min() + epsilon(lw_gas_props%get_press_min())
241 endif
242 if (t_lay(icol,ilay) .ge. lw_gas_props%get_temp_max()) then
243 t_lay(icol,ilay) = lw_gas_props%get_temp_max() - epsilon(lw_gas_props%get_temp_max())
244 endif
245 if (p_lay(icol,ilay) .ge. lw_gas_props%get_press_max()) then
246 p_lay(icol,ilay) = lw_gas_props%get_press_max() - epsilon(lw_gas_props%get_press_max())
247 endif
248 enddo
249 enddo
250
251 ! Temperature at layer-interfaces
252 call cmp_tlev(ncol,nlev,lw_gas_props%get_press_min(),p_lay,t_lay,p_lev,tsfc,t_lev)
253 do ilev=1,nlev+1
254 do icol=1,ncol
255 if (t_lev(icol,ilev) .le. lw_gas_props%get_temp_min()) t_lev(icol,ilev) = &
256 lw_gas_props%get_temp_min() + epsilon(lw_gas_props%get_temp_min())
257 if (t_lev(icol,ilev) .ge. lw_gas_props%get_temp_max()) t_lev(icol,ilev) = &
258 lw_gas_props%get_temp_max() - epsilon(lw_gas_props%get_temp_max())
259 enddo
260 enddo
261
262 ! Save surface temperature at radiation time-step, used for LW flux adjustment betwen
263 ! radiation calls.
264 tsfc_radtime = tsfc
265
266 ! Compute a bunch of thermodynamic fields needed by the cloud microphysics schemes.
267 ! Relative humidity, saturation mixing-ratio, vapor mixing-ratio, virtual temperature,
268 ! layer thickness,...
269 do ilay=1,nlev
270 do icol=1,ncol
271 es = min( p_lay(icol,ilay), fpvs( t_lay(icol,ilay) ) ) ! fpvs and prsl in pa
272 qs_lay(icol,ilay) = max( con_epsqs, con_eps * es / (p_lay(icol,ilay) + con_epsm1*es) )
273 relhum(icol,ilay) = max( 0._kind_phys, min( 1._kind_phys, max(con_epsqs, q_lay(icol,ilay))/qs_lay(icol,ilay) ) )
274 tv_lay(icol,ilay) = t_lay(icol,ilay) * (1._kind_phys + con_fvirt*q_lay(icol,ilay))
275 enddo
276 enddo
277
278 !
279 ! Compute layer-thickness between layer boundaries (deltaZ) and layer centers (deltaZc)
280 !
281 deltap = abs(p_lev(:,2:nlev+1)-p_lev(:,1:nlev))
282 con_rdog = con_rd/con_g
283 do icol=1,ncol
284 if (top_at_1) then
285 ! Layer thickness (m)
286 do ilay=1,nlev
287 deltaz(icol,ilay) = con_rdog * abs(log(p_lev(icol,ilay+1)) - log(p_lev(icol,ilay))) * tv_lay(icol,ilay)
288 enddo
289 ! Height at layer boundaries
290 hgtb(nlev+1) = 0._kind_phys
291 do ilay=nlev,1,-1
292 hgtb(ilay)= hgtb(ilay+1) + deltaz(icol,ilay)
293 enddo
294 ! Height at layer centers
295 do ilay = nlev, 1, -1
296 pfac = abs(log(p_lev(icol,ilay+1)) - log(p_lay(icol,ilay))) / &
297 abs(log(p_lev(icol,ilay+1)) - log(p_lev(icol,ilay)))
298 hgtc(ilay) = hgtb(ilay+1) + pfac * (hgtb(ilay) - hgtb(ilay+1))
299 enddo
300 ! Layer thickness between centers
301 do ilay = nlev-1, 1, -1
302 deltazc(icol,ilay) = hgtc(ilay) - hgtc(ilay+1)
303 enddo
304 deltazc(icol,nlev) = hgtc(nlev) - hgtb(nlev+1)
305 else
306 ! Layer thickness (m)
307 do ilay=nlev,1,-1
308 deltaz(icol,ilay) = con_rdog * abs(log(p_lev(icol,ilay)) - log(p_lev(icol,ilay+1))) * tv_lay(icol,ilay)
309 enddo
310 ! Height at layer boundaries
311 hgtb(1) = 0._kind_phys
312 do ilay=1,nlev
313 hgtb(ilay+1)= hgtb(ilay) + deltaz(icol,ilay)
314 enddo
315 ! Height at layer centers
316 do ilay = 1, nlev
317 pfac = abs(log(p_lev(icol,ilay)) - log(p_lay(icol,ilay) )) / &
318 abs(log(p_lev(icol,ilay)) - log(p_lev(icol,ilay+1)))
319 hgtc(ilay) = hgtb(ilay) + pfac * (hgtb(ilay+1) - hgtb(ilay))
320 enddo
321 ! Layer thickness between centers
322 do ilay = 2, nlev
323 deltazc(icol,ilay) = hgtc(ilay) - hgtc(ilay-1)
324 enddo
325 deltazc(icol,1) = hgtc(1) - hgtb(1)
326 endif
327 enddo
328
329 ! #######################################################################################
330 ! Get layer ozone mass mixing ratio
331 ! #######################################################################################
332
333 if (i_o3 > 0) then
334 do ilay=1,nlev
335 do icol=1,ncol
336 o3_lay(icol,ilay) = max( con_epsqs, qgrs(icol,ilay,i_o3) )
337 enddo
338 enddo
339 ! OR Use climatological ozone data
340 else
341 call ozphys%run_o3clim(xlat, prslk, con_pi, o3_lay)
342 endif
343
344 ! #######################################################################################
345 ! Set gas concentrations for RRTMGP
346 ! #######################################################################################
347 ! Call getgases(), to set up non-prognostic gas volume mixing ratios (gas_vmr).
348 call getgases (p_lev/100., xlon, xlat, ncol, nlev, ico2, top_at_1, con_pi, gas_vmr)
349 vmr_o2 = gas_vmr(:,:,4)
350 vmr_ch4 = gas_vmr(:,:,3)
351 vmr_n2o = gas_vmr(:,:,2)
352 vmr_co2 = gas_vmr(:,:,1)
353
354 ! Compute volume mixing-ratios for ozone (mmr) and specific-humidity.
355 vmr_h2o = merge((q_lay/(1-q_lay))*amdw, 0., q_lay .ne. 1.)
356 vmr_o3 = merge(o3_lay*amdo3, 0., o3_lay .gt. 0.)
357
358 ! #######################################################################################
359 ! Radiation time step (output) (Is this really needed?) (Used by some diagnostics)
360 ! #######################################################################################
361 raddt = min(fhswr, fhlwr)
362
363 ! #######################################################################################
364 ! Setup surface ground temperature and ground/air skin temperature if required.
365 ! #######################################################################################
366 isfc_ilev = 1
367 if (top_at_1) isfc_ilev = isfc + 1
368
369 tsfg(1:ncol) = t_lev(1:ncol,isfc_ilev)
370 tsfa(1:ncol) = t_lay(1:ncol,isfc)
371
372 ! #######################################################################################
373 ! Compute cosine of zenith angle (only when SW is called)
374 ! #######################################################################################
375 if (doswrad) then
376 call coszmn (xlon, sinlat, coslat, solhr, ncol, me, coszen, coszdg)
377 ! For SW gather daylit points
378 nday = 0
379 idxday = 0
380 do icol = 1, ncol
381 if (coszen(icol) >= 0.0001) then
382 nday = nday + 1
383 idxday(nday) = icol
384 endif
385 enddo
386 else
387 nday = 0
388 idxday = 0
389 endif
390
391 ! #######################################################################################
392 ! Surface emissivity
393 ! #######################################################################################
394 do iband=1,lw_gas_props%get_nband()
395 sfc_emiss_byband(iband,:) = semis
396 enddo
397
398 end subroutine gfs_rrtmgp_pre_run
399
400end module gfs_rrtmgp_pre
subroutine merge(len, lsoil, iy, im, id, ih, fh, deltsfc, slmskl, slmskw, sihfcs, sicfcs, vmnfcs, vmxfcs, slpfcs, absfcs, tsffcs, wetfcs, snofcs, zorfcs, albfcs, aisfcs, cvfcs, cvbfcs, cvtfcs, cnpfcs, smcfcs, stcfcs, slifcs, vegfcs, vetfcs, sotfcs, socfcs, alffcs, sihanl, sicanl, vmnanl, vmxanl, slpanl, absanl, tsfanl, tsfan2, wetanl, snoanl, zoranl, albanl, aisanl, cvanl, cvbanl, cvtanl, cnpanl, smcanl, stcanl, slianl, veganl, vetanl, sotanl, socanl, alfanl, ctsfl, calbl, caisl, csnol, csmcl, czorl, cstcl, cvegl, ctsfs, calbs, caiss, csnos, csmcs, czors, cstcs, cvegs, ccv, ccvb, ccvt, ccnp, cvetl, cvets, csotl, csots, csocl, csocs, calfl, calfs, csihl, csihs, csicl, csics, cvmnl, cvmns, cvmxl, cvmxs, cslpl, cslps, cabsl, cabss, irttsf, irtwet, irtsno, irtzor, irtalb, irtais, irttg3, irtscv, irtacn, irtsmc, irtstc, irtveg, irtvmn, irtvmx, irtslp, irtabs, irtvet, irtsot, irtsoc, irtalf, landice, me)
This subroutine merges analysis and forecast.
Definition sfcsub.F:4790
subroutine count(slimsk, sno, ijmax)
This subroutine counts number of points for the four surface conditions.
Definition sfcsub.F:2631
subroutine, public coszmn(xlon, sinlat, coslat, solhr, im, me, coszen, coszdg)
This subroutine computes mean cos solar zenith angle over SW calling interval.
subroutine, public getgases(plvl, xlon, xlat, imax, lmax, ico2flg, top_at_1, con_pi, gasdat)
This subroutine sets up global distribution of radiation absorbing gases in volume mixing ratio....
integer, parameter, public nf_vgas
This module provides an Application Program Interface (API) for computing basic thermodynamic physics...
Definition funcphys.f90:26
The operational GFS currently parameterizes ozone production and destruction based on monthly mean co...
This module sets up astronomy quantities for solar radiation calculations.
This module sets up constant gas rofiles, such as co2, ch4, n2o, o2, and those of cfc gases.
Derived type containing data and procedures needed by ozone photochemistry parameterization Note All ...