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)
109 integer,
intent(in) :: &
111 ncol, & !< Number of horizontal grid points
112 nlev, & !< Number of vertical layers
113 ico2, & !< Flag for co2 radiation scheme
117 logical,
intent(in) :: &
118 doswrad, & !< Call SW radiation?
120 real(kind_phys),
intent(in) :: &
121 fhswr, & !< Frequency of SW radiation call.
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
132 real(kind_phys),
dimension(:),
intent(in) :: &
135 tsfc, & !< Surface skin temperature (K)
136 coslat, & !< Cosine(latitude)
137 sinlat, & !< Sine(latitude)
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)
144 real(kind_phys),
dimension(:,:,:),
intent(in) :: &
146 character(len=*),
dimension(:),
intent(in),
optional :: &
150 character(len=*),
intent(out) :: &
152 integer,
intent(out) :: &
153 errflg, & !< Error flag
155 integer,
intent(inout) :: &
156 isfc, & !< Vertical index for surface
158 logical,
intent(inout) :: &
160 real(kind_phys),
intent(inout) :: &
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
168 integer,
dimension(:),
intent(inout) :: &
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
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
200 if (.not. (doswrad .or. dolwrad))
return
205 top_at_1 = (prsi(1,1) .lt. prsi(1, nlev))
221 q_lay(1:ncol,:) = qgrs(1:ncol,:,1)
222 where(q_lay .lt. 1.e-6) q_lay = 1.e-6
225 p_lev(1:ncol,:) = prsi(1:ncol,:)
228 p_lay(1:ncol,:) = prsl(1:ncol,:)
231 t_lay(1:ncol,:) = tgrs(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())
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())
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())
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())
252 call cmp_tlev(ncol,nlev,lw_gas_props%get_press_min(),p_lay,t_lay,p_lev,tsfc,t_lev)
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())
271 es = min( p_lay(icol,ilay), fpvs( t_lay(icol,ilay) ) )
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))
281 deltap = abs(p_lev(:,2:nlev+1)-p_lev(:,1:nlev))
282 con_rdog = con_rd/con_g
287 deltaz(icol,ilay) = con_rdog * abs(log(p_lev(icol,ilay+1)) - log(p_lev(icol,ilay))) * tv_lay(icol,ilay)
290 hgtb(nlev+1) = 0._kind_phys
292 hgtb(ilay)= hgtb(ilay+1) + deltaz(icol,ilay)
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))
301 do ilay = nlev-1, 1, -1
302 deltazc(icol,ilay) = hgtc(ilay) - hgtc(ilay+1)
304 deltazc(icol,nlev) = hgtc(nlev) - hgtb(nlev+1)
308 deltaz(icol,ilay) = con_rdog * abs(log(p_lev(icol,ilay)) - log(p_lev(icol,ilay+1))) * tv_lay(icol,ilay)
311 hgtb(1) = 0._kind_phys
313 hgtb(ilay+1)= hgtb(ilay) + deltaz(icol,ilay)
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))
323 deltazc(icol,ilay) = hgtc(ilay) - hgtc(ilay-1)
325 deltazc(icol,1) = hgtc(1) - hgtb(1)
336 o3_lay(icol,ilay) = max( con_epsqs, qgrs(icol,ilay,i_o3) )
341 call ozphys%run_o3clim(xlat, prslk, con_pi, o3_lay)
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)
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.)
361 raddt = min(fhswr, fhlwr)
367 if (top_at_1) isfc_ilev = isfc + 1
369 tsfg(1:ncol) = t_lev(1:ncol,isfc_ilev)
370 tsfa(1:ncol) = t_lay(1:ncol,isfc)
376 call coszmn (xlon, sinlat, coslat, solhr, ncol, me, coszen, coszdg)
381 if (coszen(icol) >= 0.0001)
then
394 do iband=1,lw_gas_props%get_nband()
395 sfc_emiss_byband(iband,:) = semis
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.