Line data Source code
1 : module radiation
2 :
3 : !---------------------------------------------------------------------------------
4 : !
5 : ! CAM interface to the legacy 'camrt' radiation code
6 : !
7 : !---------------------------------------------------------------------------------
8 :
9 : use shr_kind_mod, only: r8=>shr_kind_r8, cl=>shr_kind_cl
10 : use spmd_utils, only: masterproc
11 : use ppgrid, only: pcols, pver, pverp, begchunk, endchunk
12 : use physics_types, only: physics_state, physics_ptend
13 : use phys_grid, only: get_ncols_p
14 : use camsrfexch, only: cam_out_t, cam_in_t
15 : use physconst, only: cpair, cappa
16 : use time_manager, only: get_nstep, is_first_restart_step, &
17 : get_curr_calday, get_step_size
18 : use cam_control_mod, only: lambm0, obliqr, mvelpp, eccen
19 :
20 : use radae, only: abstot_3d, absnxt_3d, emstot_3d, initialize_radbuffer, ntoplw
21 :
22 : use scamMod, only: scm_crm_mode, single_column,have_cld,cldobs,&
23 : have_clwp,clwpobs,have_tg,tground
24 :
25 : use cam_grid_support, only: cam_grid_write_attr, cam_grid_id, &
26 : cam_grid_header_info_t, cam_grid_dimensions, &
27 : cam_grid_write_dist_array, cam_grid_read_dist_array
28 :
29 : use cam_history, only: outfld, hist_fld_active
30 : use cam_history_support, only: fillvalue
31 :
32 : use cam_pio_utils, only: cam_pio_def_dim
33 :
34 : use pio, only: file_desc_t, var_desc_t, &
35 : pio_double, pio_int, pio_noerr, &
36 : pio_seterrorhandling, pio_bcast_error, &
37 : pio_inq_varid, &
38 : pio_def_var, pio_def_dim, &
39 : pio_put_var, pio_get_var, pio_put_att
40 :
41 : use cam_abortutils, only: endrun
42 : use error_messages, only: handle_err
43 : use perf_mod, only: t_startf, t_stopf
44 : use cam_logfile, only: iulog
45 :
46 : implicit none
47 : private
48 : save
49 :
50 : public :: &
51 : radiation_readnl, &! read namelist variables
52 : radiation_register, &! registers radiation physics buffer fields
53 : radiation_nextsw_cday, &! calendar day of next radiation calculation
54 : radiation_do, &! query which radiation calcs are done this timestep
55 : radiation_init, &! initialization
56 : radiation_define_restart, &!
57 : radiation_write_restart, &!
58 : radiation_read_restart, &!
59 : radiation_tend, &! compute heating rates and fluxes
60 : rad_out_t ! type for diagnostic outputs
61 :
62 : real(r8), public, protected :: nextsw_cday ! future radiation calday for surface models
63 :
64 : type rad_out_t
65 : real(r8) :: solin(pcols) ! Solar incident flux
66 : real(r8) :: fsntoa(pcols) ! Net solar flux at TOA
67 : real(r8) :: fsutoa(pcols) ! upwelling solar flux at TOA
68 : real(r8) :: fsntoac(pcols) ! Clear sky net solar flux at TOA
69 : real(r8) :: fsnirt(pcols) ! Near-IR flux absorbed at toa
70 : real(r8) :: fsnrtc(pcols) ! Clear sky near-IR flux absorbed at toa
71 : real(r8) :: fsnirtsq(pcols) ! Near-IR flux absorbed at toa >= 0.7 microns
72 : real(r8) :: fsntc(pcols) ! Clear sky total column abs solar flux
73 : real(r8) :: fsnsc(pcols) ! Clear sky surface abs solar flux
74 : real(r8) :: fsdsc(pcols) ! Clear sky surface downwelling solar flux
75 : real(r8) :: flut(pcols) ! Upward flux at top of model
76 : real(r8) :: flutc(pcols) ! Upward Clear Sky flux at top of model
77 : real(r8) :: flntc(pcols) ! Clear sky lw flux at model top
78 : real(r8) :: flnsc(pcols) ! Clear sky lw flux at srf (up-down)
79 : real(r8) :: fldsc(pcols) ! Clear sky lw flux at srf (down)
80 : real(r8) :: flwds(pcols) ! Down longwave flux at surface
81 : real(r8) :: fsnr(pcols)
82 : real(r8) :: flnr(pcols)
83 : real(r8) :: fsds(pcols) ! Surface solar down flux
84 : real(r8) :: fln200(pcols) ! net longwave flux interpolated to 200 mb
85 : real(r8) :: fln200c(pcols) ! net clearsky longwave flux interpolated to 200 mb
86 : real(r8) :: fsn200(pcols) ! fns interpolated to 200 mb
87 : real(r8) :: fsn200c(pcols) ! fcns interpolated to 200 mb
88 : real(r8) :: sols(pcols) ! Solar downward visible direct to surface
89 : real(r8) :: soll(pcols) ! Solar downward near infrared direct to surface
90 : real(r8) :: solsd(pcols) ! Solar downward visible diffuse to surface
91 : real(r8) :: solld(pcols) ! Solar downward near infrared diffuse to surface
92 : real(r8) :: qrsc(pcols,pver) ! clearsky shortwave radiative heating rate
93 : real(r8) :: qrlc(pcols,pver) ! clearsky longwave radiative heating rate
94 : real(r8) :: fsdtoa(pcols) ! Solar input = Flux Solar Downward Top of Atmosphere
95 : real(r8) :: swcf(pcols) ! shortwave cloud forcing
96 : real(r8) :: lwcf(pcols) ! longwave cloud forcing
97 :
98 : real(r8) :: tot_cld_vistau(pcols,pver)
99 : real(r8) :: tot_icld_vistau(pcols,pver)
100 : real(r8) :: liq_icld_vistau(pcols,pver) ! in-cld liq cloud optical depth (only during day, night = fillvalue)
101 : real(r8) :: ice_icld_vistau(pcols,pver) ! in-cld ice cloud optical depth (only during day, night = fillvalue)
102 : end type rad_out_t
103 :
104 : ! Namelist variables
105 :
106 : character(len=cl) :: absems_data
107 : integer :: iradsw = -1 ! freq. of shortwave radiation calc in time steps (positive)
108 : ! or hours (negative).
109 : integer :: iradlw = -1 ! frequency of longwave rad. calc. in time steps (positive)
110 : ! or hours (negative).
111 : integer :: iradae = -12 ! frequency of absorp/emis calc in time steps (positive)
112 : ! or hours (negative).
113 : integer :: irad_always = 0 ! Specifies length of time in timesteps (positive)
114 : ! or hours (negative) SW/LW radiation will be
115 : ! run continuously from the start of an
116 : ! initial or restart run
117 : logical :: use_rad_dt_cosz = .false. ! if true use zenith angle averaged over
118 : ! interval between radiation calculations
119 :
120 : ! Physics buffer indices
121 : integer :: qrs_idx = 0
122 : integer :: qrl_idx = 0
123 : integer :: fsds_idx = 0
124 : integer :: fsns_idx = 0
125 : integer :: fsnt_idx = 0
126 : integer :: flns_idx = 0
127 : integer :: flnt_idx = 0
128 : integer :: cld_idx = 0
129 : integer :: rel_idx = 0
130 : integer :: rei_idx = 0
131 : integer :: cicewp_idx = -1
132 : integer :: cliqwp_idx = -1
133 : integer :: cldemis_idx = -1
134 : integer :: cldtau_idx = -1
135 : integer :: nmxrgn_idx = -1
136 : integer :: pmxrgn_idx = -1
137 :
138 : ! averaging time interval for zenith angle
139 : real(r8) :: dt_avg = 0._r8
140 :
141 : real(r8), parameter :: cgs2mks = 1.e-3_r8
142 :
143 : ! PIO descriptors (for restarts)
144 :
145 : type(var_desc_t), allocatable :: abstot_desc(:)
146 : type(var_desc_t) :: emstot_desc, absnxt_desc(4)
147 : type(var_desc_t) :: nextsw_cday_desc
148 :
149 : logical :: use_rad_uniform_angle = .false. ! if true, use the namelist rad_uniform_angle for the zenith calculation
150 : real(r8) :: rad_uniform_angle = -99._r8
151 :
152 : !===============================================================================
153 : contains
154 : !===============================================================================
155 :
156 1536 : subroutine radiation_readnl(nlfile)
157 :
158 : ! Read radiation_nl namelist group.
159 :
160 : use namelist_utils, only: find_group_name
161 : use units, only: getunit, freeunit
162 : use spmd_utils, only: mpicom, mstrid=>masterprocid, mpi_integer, mpi_logical, &
163 : mpi_character, mpi_real8
164 :
165 : character(len=*), intent(in) :: nlfile ! filepath for file containing namelist input
166 :
167 : ! Local variables
168 : integer :: unitn, ierr
169 : integer :: dtime ! timestep size
170 : character(len=*), parameter :: sub = 'radiation_readnl'
171 :
172 : namelist /radiation_nl/ absems_data, iradsw, iradlw, iradae, irad_always, &
173 : use_rad_dt_cosz, use_rad_uniform_angle, rad_uniform_angle
174 : !-----------------------------------------------------------------------------
175 :
176 1536 : if (masterproc) then
177 2 : unitn = getunit()
178 2 : open( unitn, file=trim(nlfile), status='old' )
179 2 : call find_group_name(unitn, 'radiation_nl', status=ierr)
180 2 : if (ierr == 0) then
181 2 : read(unitn, radiation_nl, iostat=ierr)
182 2 : if (ierr /= 0) then
183 0 : call endrun(sub // ':: ERROR reading namelist')
184 : end if
185 : end if
186 2 : close(unitn)
187 2 : call freeunit(unitn)
188 : end if
189 :
190 : ! Broadcast namelist variables
191 1536 : call mpi_bcast(absems_data, len(absems_data), mpi_character, mstrid, mpicom, ierr)
192 1536 : if (ierr /= 0) call endrun(sub//": FATAL: mpi_bcast: absems_data")
193 1536 : call mpi_bcast(iradsw, 1, mpi_integer, mstrid, mpicom, ierr)
194 1536 : if (ierr /= 0) call endrun(sub//": FATAL: mpi_bcast: iradsw")
195 1536 : call mpi_bcast(iradlw, 1, mpi_integer, mstrid, mpicom, ierr)
196 1536 : if (ierr /= 0) call endrun(sub//": FATAL: mpi_bcast: iradlw")
197 1536 : call mpi_bcast(iradae, 1, mpi_integer, mstrid, mpicom, ierr)
198 1536 : if (ierr /= 0) call endrun(sub//": FATAL: mpi_bcast: iradae")
199 1536 : call mpi_bcast(irad_always, 1, mpi_integer, mstrid, mpicom, ierr)
200 1536 : if (ierr /= 0) call endrun(sub//": FATAL: mpi_bcast: irad_always")
201 1536 : call mpi_bcast(use_rad_dt_cosz, 1, mpi_logical, mstrid, mpicom, ierr)
202 1536 : if (ierr /= 0) call endrun(sub//": FATAL: mpi_bcast: use_rad_dt_cosz")
203 1536 : call mpi_bcast(use_rad_uniform_angle, 1, mpi_logical, mstrid, mpicom, ierr)
204 1536 : if (ierr /= 0) call endrun(sub//": FATAL: mpi_bcast: use_rad_uniform_angle")
205 1536 : call mpi_bcast(rad_uniform_angle, 1, mpi_real8, mstrid, mpicom, ierr)
206 1536 : if (ierr /= 0) call endrun(sub//": FATAL: mpi_bcast: rad_uniform_angle")
207 :
208 1536 : if (use_rad_uniform_angle .and. rad_uniform_angle == -99._r8) then
209 0 : call endrun(sub // ' ERROR - use_rad_uniform_angle is set to .true, but rad_uniform_angle is not set ')
210 : end if
211 :
212 : ! Convert iradsw, iradlw and irad_always from hours to timesteps if necessary
213 1536 : dtime = get_step_size()
214 1536 : if (iradsw < 0) iradsw = nint((-iradsw *3600._r8)/dtime)
215 1536 : if (iradlw < 0) iradlw = nint((-iradlw *3600._r8)/dtime)
216 1536 : if (irad_always < 0) irad_always = nint((-irad_always*3600._r8)/dtime)
217 :
218 : ! Convert iradae from hours to timesteps if necessary and check that
219 : ! iradae must be an even multiple of iradlw
220 1536 : if (iradae < 0) iradae = nint((-iradae*3600._r8)/dtime)
221 1536 : if (mod(iradae,iradlw)/=0) then
222 0 : write(iulog,*) sub//': iradae must be an even multiple of iradlw.'
223 0 : write(iulog,*)' iradae = ',iradae,', iradlw = ',iradlw
224 0 : call endrun(sub//': iradae must be an even multiple of iradlw.')
225 : end if
226 :
227 : !-----------------------------------------------------------------------
228 : ! Print runtime options to log.
229 : !-----------------------------------------------------------------------
230 :
231 1536 : if (masterproc) then
232 2 : write(iulog,*) 'CAMRT radiation scheme parameters:'
233 2 : write(iulog,10) iradsw, iradlw, iradae, irad_always, use_rad_dt_cosz
234 2 : write(iulog,*) ' Abs/Emis dataset: ', trim(absems_data)
235 : end if
236 :
237 : 10 format(' Frequency (timesteps) of Shortwave Radiation calc: ',i5/, &
238 : ' Frequency (timesteps) of Longwave Radiation calc: ',i5/, &
239 : ' Frequency (timesteps) of Absorptivity/Emissivity calc:',i5/, &
240 : ' SW/LW calc done every timestep for first N steps. N= ',i5/, &
241 : ' Use average zenith angle: ',l5)
242 :
243 :
244 1536 : end subroutine radiation_readnl
245 :
246 : !================================================================================================
247 :
248 1536 : subroutine radiation_register
249 :
250 : ! Register radiation fields in the physics buffer
251 :
252 : use physics_buffer, only: pbuf_add_field, dtype_r8
253 : use radiation_data, only: rad_data_register
254 :
255 1536 : call pbuf_add_field('QRS' , 'global',dtype_r8,(/pcols,pver/), qrs_idx) ! shortwave radiative heating rate
256 1536 : call pbuf_add_field('QRL' , 'global',dtype_r8,(/pcols,pver/), qrl_idx) ! longwave radiative heating rate
257 :
258 1536 : call pbuf_add_field('FSDS' , 'global',dtype_r8,(/pcols/), fsds_idx) ! Surface solar downward flux
259 :
260 1536 : call pbuf_add_field('FSNS' , 'global',dtype_r8,(/pcols/), fsns_idx) ! Surface net shortwave flux
261 1536 : call pbuf_add_field('FSNT' , 'global',dtype_r8,(/pcols/), fsnt_idx) ! Top-of-model net shortwave flux
262 1536 : call pbuf_add_field('FLNS' , 'global',dtype_r8,(/pcols/), flns_idx) ! Surface net longwave flux
263 1536 : call pbuf_add_field('FLNT' , 'global',dtype_r8,(/pcols/), flnt_idx) ! Top-of-model net longwave flux
264 :
265 1536 : call rad_data_register()
266 :
267 1536 : end subroutine radiation_register
268 :
269 : !================================================================================================
270 :
271 9731472 : function radiation_do(op, timestep)
272 :
273 : ! Returns true if the specified operation is done this timestep.
274 :
275 : character(len=*), intent(in) :: op ! name of operation
276 : integer, intent(in), optional:: timestep
277 : logical :: radiation_do ! return value
278 :
279 : ! Local variables
280 : integer :: nstep ! current timestep number
281 : !-----------------------------------------------------------------------
282 :
283 9731472 : if (present(timestep)) then
284 2250792 : nstep = timestep
285 : else
286 7480680 : nstep = get_nstep()
287 : end if
288 :
289 5241528 : select case (op)
290 :
291 : case ('sw') ! do a shortwave heating calc this timestep?
292 : radiation_do = nstep == 0 .or. iradsw == 1 &
293 : .or. (mod(nstep-1,iradsw) == 0 .and. nstep /= 1) &
294 5241528 : .or. nstep <= irad_always
295 :
296 : case ('lw') ! do a longwave heating calc this timestep?
297 : radiation_do = nstep == 0 .or. iradlw == 1 &
298 : .or. (mod(nstep-1,iradlw) == 0 .and. nstep /= 1) &
299 2990736 : .or. nstep <= irad_always
300 :
301 : case ('absems') ! do an absorptivity/emissivity calculation this timestep?
302 : radiation_do = nstep == 0 .or. iradae == 1 &
303 1495368 : .or. (mod(nstep-1,iradae) == 0 .and. nstep /= 1)
304 :
305 : case ('aeres') ! write absorptivity/emissivity to restart file this timestep?
306 3840 : radiation_do = mod(nstep,iradae) /= 0
307 :
308 : case default
309 9731472 : call endrun('radiation_do: unknown operation:'//op)
310 :
311 : end select
312 1536 : end function radiation_do
313 :
314 : !================================================================================================
315 :
316 1495368 : real(r8) function radiation_nextsw_cday()
317 :
318 : ! Returns calendar day of next sw radiation calculation
319 :
320 : ! Local variables
321 : integer :: nstep ! timestep counter
322 : logical :: dosw ! true => do shosrtwave calc
323 : integer :: offset ! offset for calendar day calculation
324 : integer :: dtime ! integer timestep size
325 : real(r8):: calday ! calendar day of
326 : real(r8):: caldayp1 ! calendar day of next time-step
327 : !-----------------------------------------------------------------------
328 :
329 1495368 : radiation_nextsw_cday = -1._r8
330 1495368 : dosw = .false.
331 1495368 : nstep = get_nstep()
332 1495368 : dtime = get_step_size()
333 1495368 : offset = 0
334 : do while (.not. dosw)
335 2250792 : nstep = nstep + 1
336 2250792 : offset = offset + dtime
337 2250792 : if (radiation_do('sw', nstep)) then
338 1495368 : radiation_nextsw_cday = get_curr_calday(offset=offset)
339 : dosw = .true.
340 : end if
341 : end do
342 1495368 : if(radiation_nextsw_cday == -1._r8) then
343 0 : call endrun('error in radiation_nextsw_cday')
344 : end if
345 :
346 : ! determine if next radiation time-step not equal to next time-step
347 1495368 : if (get_nstep() >= 1) then
348 1492272 : caldayp1 = get_curr_calday(offset=int(dtime))
349 1492272 : if (caldayp1 /= radiation_nextsw_cday) radiation_nextsw_cday = -1._r8
350 : end if
351 :
352 1495368 : end function radiation_nextsw_cday
353 :
354 : !================================================================================================
355 :
356 1536 : subroutine radiation_init(pbuf2d)
357 :
358 : ! Initialize the radiation parameterization, add fields to the history buffer
359 :
360 : use cam_history, only: addfld, add_default, horiz_only
361 : use physconst, only: gravit, cpair, epsilo, stebol, &
362 : pstd, mwdry, mwco2, mwo3
363 :
364 : use physics_buffer, only: physics_buffer_desc, pbuf_get_index
365 : use radsw, only: radsw_init
366 : use radlw, only: radlw_init
367 : use radae, only: radae_init
368 : use radconstants, only: radconstants_init
369 : use rad_solar_var, only: rad_solar_var_init
370 : use radiation_data, only: rad_data_init
371 : use phys_control, only: phys_getopts
372 : use time_manager, only: get_step_size
373 :
374 : ! args
375 : type(physics_buffer_desc), pointer :: pbuf2d(:,:)
376 :
377 :
378 : ! Local variables
379 : integer :: nstep ! current timestep number
380 : logical :: history_amwg ! output the variables used by the AMWG diag package
381 : logical :: history_vdiag ! output the variables used by the AMWG variability diag package
382 : logical :: history_budget ! output tendencies and state variables for CAM4
383 : ! temperature, water vapor, cloud ice and cloud
384 : ! liquid budgets.
385 : integer :: history_budget_histfile_num ! output history file number for budget fields
386 :
387 : integer :: dtime
388 : !-----------------------------------------------------------------------
389 :
390 1536 : call radconstants_init()
391 1536 : call rad_solar_var_init()
392 :
393 1536 : call radsw_init(gravit)
394 1536 : call radlw_init(gravit, stebol)
395 : call radae_init( &
396 : gravit, epsilo, stebol, pstd, mwdry, &
397 1536 : mwco2, mwo3, absems_data)
398 :
399 1536 : call rad_data_init(pbuf2d)
400 :
401 : ! Set the radiation timestep for cosz calculations if requested using the adjusted iradsw value from radiation
402 1536 : if (use_rad_dt_cosz) then
403 0 : dtime = get_step_size()
404 0 : dt_avg = real(iradsw*dtime, r8)
405 : end if
406 :
407 : ! Surface components to get radiation computed today
408 1536 : if (.not. is_first_restart_step()) then
409 768 : nextsw_cday = get_curr_calday()
410 : end if
411 :
412 : ! Get physics buffer indices
413 1536 : cld_idx = pbuf_get_index('CLD')
414 1536 : rel_idx = pbuf_get_index('REL')
415 1536 : rei_idx = pbuf_get_index('REI')
416 :
417 : ! "irad_always" is number of time steps to execute radiation continuously from start of
418 : ! initial OR restart run
419 1536 : nstep = get_nstep()
420 1536 : if ( irad_always > 0) then
421 0 : nstep = get_nstep()
422 0 : irad_always = irad_always + nstep
423 : end if
424 :
425 : ! Shortwave radiation
426 1536 : call addfld ('SOLIN', horiz_only, 'A','W/m2','Solar insolation', sampling_seq='rad_lwsw')
427 : call addfld ('SOLL', horiz_only, 'A','W/m2','Solar downward near infrared direct to surface', &
428 1536 : sampling_seq='rad_lwsw')
429 1536 : call addfld ('SOLS', horiz_only, 'A','W/m2','Solar downward visible direct to surface', sampling_seq='rad_lwsw')
430 : call addfld ('SOLLD', horiz_only, 'A','W/m2','Solar downward near infrared diffuse to surface', &
431 1536 : sampling_seq='rad_lwsw')
432 1536 : call addfld ('SOLSD', horiz_only, 'A','W/m2','Solar downward visible diffuse to surface', sampling_seq='rad_lwsw')
433 3072 : call addfld ('QRS', (/ 'lev' /), 'A','K/s', 'Solar heating rate', sampling_seq='rad_lwsw')
434 3072 : call addfld ('QRSC', (/ 'lev' /), 'A','K/s', 'Clearsky solar heating rate', sampling_seq='rad_lwsw')
435 1536 : call addfld ('FSNS', horiz_only, 'A','W/m2','Net solar flux at surface', sampling_seq='rad_lwsw')
436 1536 : call addfld ('FSNT', horiz_only, 'A','W/m2','Net solar flux at top of model', sampling_seq='rad_lwsw')
437 1536 : call addfld ('FSNTOA', horiz_only, 'A','W/m2','Net solar flux at top of atmosphere', sampling_seq='rad_lwsw')
438 1536 : call addfld ('FSUTOA', horiz_only, 'A','W/m2','Upwelling solar flux at top of atmosphere', sampling_seq='rad_lwsw')
439 : call addfld ('FSNTOAC', horiz_only, 'A','W/m2','Clearsky net solar flux at top of atmosphere', &
440 1536 : sampling_seq='rad_lwsw')
441 1536 : call addfld ('FSDTOA', horiz_only, 'A','W/m2','Downwelling solar flux at top of atmosphere', sampling_seq='rad_lwsw')
442 1536 : call addfld ('FSN200', horiz_only, 'A','W/m2','Net shortwave flux at 200 mb', sampling_seq='rad_lwsw')
443 1536 : call addfld ('FSN200C', horiz_only, 'A','W/m2','Clearsky net shortwave flux at 200 mb', sampling_seq='rad_lwsw')
444 1536 : call addfld ('FSNTC', horiz_only, 'A','W/m2','Clearsky net solar flux at top of model', sampling_seq='rad_lwsw')
445 1536 : call addfld ('FSNSC', horiz_only, 'A','W/m2','Clearsky net solar flux at surface', sampling_seq='rad_lwsw')
446 : call addfld ('FSDSC', horiz_only, 'A','W/m2','Clearsky downwelling solar flux at surface', &
447 1536 : sampling_seq='rad_lwsw')
448 1536 : call addfld ('FSDS', horiz_only, 'A','W/m2','Downwelling solar flux at surface', sampling_seq='rad_lwsw')
449 3072 : call addfld ('FUS', (/ 'ilev' /), 'I','W/m2','Shortwave upward flux')
450 3072 : call addfld ('FDS', (/ 'ilev' /), 'I','W/m2','Shortwave downward flux')
451 3072 : call addfld ('FUSC', (/ 'ilev' /), 'I','W/m2','Shortwave clear-sky upward flux')
452 3072 : call addfld ('FDSC', (/ 'ilev' /), 'I','W/m2','Shortwave clear-sky downward flux')
453 : call addfld ('FSNIRTOA', horiz_only, 'A','W/m2','Net near-infrared flux (Nimbus-7 WFOV) at top of atmosphere', &
454 1536 : sampling_seq='rad_lwsw')
455 : call addfld ('FSNRTOAC', horiz_only, 'A','W/m2', &
456 1536 : 'Clearsky net near-infrared flux (Nimbus-7 WFOV) at top of atmosphere', sampling_seq='rad_lwsw')
457 : call addfld ('FSNRTOAS', horiz_only, 'A','W/m2','Net near-infrared flux (>= 0.7 microns) at top of atmosphere', &
458 1536 : sampling_seq='rad_lwsw')
459 1536 : call addfld ('FSNR', horiz_only, 'A','W/m2','Net solar flux at tropopause', sampling_seq='rad_lwsw')
460 1536 : call addfld ('SWCF', horiz_only, 'A','W/m2','Shortwave cloud forcing', sampling_seq='rad_lwsw')
461 :
462 : call addfld ('TOT_CLD_VISTAU', (/ 'lev' /), 'A','1', 'Total gbx cloud visible sw optical depth', &
463 3072 : sampling_seq='rad_lwsw',flag_xyfill=.true.)
464 : call addfld ('TOT_ICLD_VISTAU', (/ 'lev' /), 'A','1', 'Total in-cloud visible sw optical depth', &
465 3072 : sampling_seq='rad_lwsw',flag_xyfill=.true.)
466 : call addfld ('LIQ_ICLD_VISTAU', (/ 'lev' /), 'A','1', 'Liquid in-cloud visible sw optical depth', &
467 3072 : sampling_seq='rad_lwsw',flag_xyfill=.true.)
468 : call addfld ('ICE_ICLD_VISTAU', (/ 'lev' /), 'A','1', 'Ice in-cloud visible sw optical depth', &
469 3072 : sampling_seq='rad_lwsw',flag_xyfill=.true.)
470 :
471 : ! Longwave radiation
472 3072 : call addfld ('QRL', (/ 'lev' /), 'A','K/s', 'Longwave heating rate', sampling_seq='rad_lwsw')
473 3072 : call addfld ('QRLC', (/ 'lev' /), 'A','K/s', 'Clearsky longwave heating rate', sampling_seq='rad_lwsw')
474 1536 : call addfld ('FLNS', horiz_only, 'A','W/m2','Net longwave flux at surface', sampling_seq='rad_lwsw')
475 1536 : call addfld ('FLDS', horiz_only, 'A','W/m2','Downwelling longwave flux at surface', sampling_seq='rad_lwsw')
476 1536 : call addfld ('FLNT', horiz_only, 'A','W/m2','Net longwave flux at top of model', sampling_seq='rad_lwsw')
477 1536 : call addfld ('FLUT', horiz_only, 'A','W/m2','Upwelling longwave flux at top of model', sampling_seq='rad_lwsw')
478 : call addfld ('FLUTC', horiz_only, 'A','W/m2','Clearsky upwelling longwave flux at top of model', &
479 1536 : sampling_seq='rad_lwsw')
480 1536 : call addfld ('FLNTC', horiz_only, 'A','W/m2','Clearsky net longwave flux at top of model', sampling_seq='rad_lwsw')
481 1536 : call addfld ('FLN200', horiz_only, 'A','W/m2','Net longwave flux at 200 mb', sampling_seq='rad_lwsw')
482 1536 : call addfld ('FLN200C', horiz_only, 'A','W/m2','Clearsky net longwave flux at 200 mb', sampling_seq='rad_lwsw')
483 1536 : call addfld ('FLNR', horiz_only, 'A','W/m2','Net longwave flux at tropopause', sampling_seq='rad_lwsw')
484 1536 : call addfld ('FLNSC', horiz_only, 'A','W/m2','Clearsky net longwave flux at surface', sampling_seq='rad_lwsw')
485 : call addfld ('FLDSC', horiz_only, 'A','W/m2','Clearsky downwelling longwave flux at surface', &
486 1536 : sampling_seq='rad_lwsw')
487 1536 : call addfld ('LWCF', horiz_only, 'A','W/m2','Longwave cloud forcing', sampling_seq='rad_lwsw')
488 3072 : call addfld ('FUL', (/ 'ilev' /), 'I','W/m2','Longwave upward flux')
489 3072 : call addfld ('FDL', (/ 'ilev' /), 'I','W/m2','Longwave downward flux')
490 3072 : call addfld ('FULC', (/ 'ilev' /), 'I','W/m2','Longwave clear-sky upward flux')
491 3072 : call addfld ('FDLC', (/ 'ilev' /), 'I','W/m2','Longwave clear-sky downward flux')
492 :
493 : ! Heating rate needed for d(theta)/dt computation
494 3072 : call addfld ('HR', (/ 'lev' /), 'A','K/s', 'Heating rate needed for d(theta)/dt computation')
495 :
496 : ! determine default variables
497 : call phys_getopts(history_amwg_out = history_amwg, &
498 : history_vdiag_out = history_vdiag, &
499 : history_budget_out = history_budget, &
500 1536 : history_budget_histfile_num_out = history_budget_histfile_num)
501 :
502 1536 : if (history_amwg) then
503 : ! Shortwave variables
504 1536 : call add_default ('SOLIN ', 1, ' ')
505 1536 : call add_default ('QRS ', 1, ' ')
506 1536 : call add_default ('FSNS ', 1, ' ')
507 1536 : call add_default ('FSNT ', 1, ' ')
508 1536 : call add_default ('FSDTOA ', 1, ' ')
509 1536 : call add_default ('FSNTOA ', 1, ' ')
510 1536 : call add_default ('FSUTOA ', 1, ' ')
511 1536 : call add_default ('FSNTOAC ', 1, ' ')
512 1536 : call add_default ('FSNTC ', 1, ' ')
513 1536 : call add_default ('FSNSC ', 1, ' ')
514 1536 : call add_default ('FSDSC ', 1, ' ')
515 1536 : call add_default ('FSDS ', 1, ' ')
516 1536 : call add_default ('SWCF ', 1, ' ')
517 : ! Longwave variables
518 1536 : call add_default ('QRL ', 1, ' ')
519 1536 : call add_default ('FLNS ', 1, ' ')
520 1536 : call add_default ('FLDS ', 1, ' ')
521 1536 : call add_default ('FLNT ', 1, ' ')
522 1536 : call add_default ('FLUT ', 1, ' ')
523 1536 : call add_default ('FLUTC ', 1, ' ')
524 1536 : call add_default ('FLNTC ', 1, ' ')
525 1536 : call add_default ('FLNSC ', 1, ' ')
526 1536 : call add_default ('FLDSC ', 1, ' ')
527 1536 : call add_default ('LWCF ', 1, ' ')
528 : endif
529 1536 : if (single_column.and.scm_crm_mode) then
530 : ! Shortwave variables
531 0 : call add_default ('FUS ', 1, ' ')
532 0 : call add_default ('FUSC ', 1, ' ')
533 0 : call add_default ('FDS ', 1, ' ')
534 0 : call add_default ('FDSC ', 1, ' ')
535 : ! Longwave variables
536 0 : call add_default ('FUL ', 1, ' ')
537 0 : call add_default ('FULC ', 1, ' ')
538 0 : call add_default ('FDL ', 1, ' ')
539 0 : call add_default ('FDLC ', 1, ' ')
540 : endif
541 :
542 1536 : if ( history_budget .and. history_budget_histfile_num > 1 ) then
543 0 : call add_default ('QRL ', history_budget_histfile_num, ' ')
544 0 : call add_default ('QRS ', history_budget_histfile_num, ' ')
545 : end if
546 :
547 1536 : if (history_vdiag) then
548 0 : call add_default('FLUT',2,' ')
549 0 : call add_default('FLUT',3,' ')
550 : end if
551 :
552 1536 : cicewp_idx = pbuf_get_index('CICEWP')
553 1536 : cliqwp_idx = pbuf_get_index('CLIQWP')
554 1536 : cldemis_idx= pbuf_get_index('CLDEMIS')
555 1536 : cldtau_idx = pbuf_get_index('CLDTAU')
556 1536 : nmxrgn_idx = pbuf_get_index('NMXRGN')
557 1536 : pmxrgn_idx = pbuf_get_index('PMXRGN')
558 :
559 1536 : end subroutine radiation_init
560 :
561 : !===============================================================================
562 :
563 1536 : subroutine radiation_define_restart(file)
564 :
565 : ! define variables to be written to restart file
566 :
567 : ! arguments
568 : type(file_desc_t), intent(inout) :: file
569 :
570 : ! local variables
571 : integer :: i, ierr
572 : integer :: grid_id
573 : integer :: hdimcnt
574 : integer :: pver_id, pverp_id
575 : integer :: vsize
576 : integer :: dimids(4)
577 :
578 1536 : type(cam_grid_header_info_t) :: info
579 :
580 : character(len=16) :: pname
581 : !----------------------------------------------------------------------------
582 :
583 1536 : call pio_seterrorhandling(File, PIO_BCAST_ERROR)
584 :
585 1536 : ierr = pio_def_var(File, 'nextsw_cday', pio_double, nextsw_cday_desc)
586 1536 : ierr = pio_put_att(File, nextsw_cday_desc, 'long_name', 'future radiation calday for surface models')
587 :
588 1536 : if (radiation_do('aeres')) then
589 :
590 0 : grid_id = cam_grid_id('physgrid')
591 0 : call cam_grid_write_attr(File, grid_id, info)
592 0 : hdimcnt = info%num_hdims()
593 0 : do i = 1, hdimcnt
594 0 : dimids(i) = info%get_hdimid(i)
595 : end do
596 :
597 0 : call cam_pio_def_dim(File, 'lev', pver, pver_id, existOK=.true.)
598 0 : call cam_pio_def_dim(File, 'ilev', pverp, pverp_id, existOK=.true.)
599 :
600 0 : vsize = pverp - ntoplw + 1
601 0 : if (vsize /= pverp) then
602 0 : ierr = pio_def_dim(File, 'lwcols', vsize, dimids(hdimcnt+1))
603 : else
604 0 : dimids(hdimcnt+1) = pverp_id
605 : end if
606 :
607 : ! split into vsize variables to avoid excessive memory usage in IO
608 :
609 0 : allocate(abstot_desc(ntoplw:pverp))
610 :
611 0 : do i = ntoplw, pverp
612 0 : write(pname,'(a,i3.3)') 'NAL_absorp', i
613 0 : ierr = pio_def_var(File, trim(pname), pio_double, dimids(1:hdimcnt+1), abstot_desc(i))
614 : end do
615 :
616 0 : dimids(hdimcnt+1) = pverp_id
617 0 : ierr = pio_def_var(File, 'Emissivity', pio_double, dimids(1:hdimcnt+1), emstot_desc)
618 :
619 0 : dimids(hdimcnt+1) = pver_id
620 0 : do i=1,4
621 0 : write(pname,'(a,i3.3)') 'NN_absorp',i
622 0 : ierr = pio_def_var(File, pname, pio_double, dimids(1:hdimcnt+1), absnxt_desc(i))
623 : end do
624 :
625 : end if
626 :
627 3072 : end subroutine radiation_define_restart
628 :
629 : !===============================================================================
630 :
631 1536 : subroutine radiation_write_restart(file)
632 :
633 : ! write variables to restart file
634 :
635 : ! arguments
636 : type(file_desc_t), intent(inout) :: file
637 :
638 : ! local variables
639 : integer :: i, ierr
640 : integer :: physgrid
641 : integer :: dims(3), gdims(3)
642 : integer :: nhdims
643 : integer :: ncol
644 : !----------------------------------------------------------------------------
645 :
646 3072 : ierr = pio_put_var(File, nextsw_cday_desc, (/ nextsw_cday /))
647 :
648 1536 : if ( radiation_do('aeres') ) then
649 :
650 0 : physgrid = cam_grid_id('physgrid')
651 0 : call cam_grid_dimensions(physgrid, gdims(1:2), nhdims)
652 :
653 0 : do i = begchunk, endchunk
654 0 : ncol = get_ncols_p(i)
655 0 : if (ncol < pcols) then
656 0 : abstot_3d(ncol+1:pcols,:,:,i) = fillvalue
657 0 : absnxt_3d(ncol+1:pcols,:,:,i) = fillvalue
658 0 : emstot_3d(ncol+1:pcols,:,i) = fillvalue
659 : end if
660 : end do
661 :
662 : ! abstot_3d is written as a series of 3D variables
663 :
664 0 : dims(1) = size(abstot_3d, 1) ! Should be pcols
665 0 : dims(2) = size(abstot_3d, 2) ! Should be (pverp-ntoplw+1)
666 0 : dims(3) = size(abstot_3d, 4) ! Should be endchunk - begchunk + 1
667 0 : gdims(nhdims+1) = dims(2)
668 0 : do i = ntoplw, pverp
669 : call cam_grid_write_dist_array(File, physgrid, dims(1:3), &
670 0 : gdims(1:nhdims+1), abstot_3d(:,:,i,:), abstot_desc(i))
671 : end do
672 :
673 0 : dims(1) = size(emstot_3d, 1) ! Should be pcols
674 0 : dims(2) = size(emstot_3d, 2) ! Should be pverp
675 0 : dims(3) = size(emstot_3d, 3) ! Should be endchunk - begchunk + 1
676 0 : gdims(nhdims+1) = dims(2)
677 : call cam_grid_write_dist_array(File, physgrid, dims(1:3), &
678 0 : gdims(1:nhdims+1), emstot_3d, emstot_desc)
679 :
680 0 : dims(1) = size(absnxt_3d, 1) ! Should be pcols
681 0 : dims(2) = size(absnxt_3d, 2) ! Should be pver
682 0 : dims(3) = size(absnxt_3d, 4) ! Should be endchunk - begchunk + 1
683 0 : gdims(nhdims+1) = dims(2)
684 0 : do i = 1, 4
685 : call cam_grid_write_dist_array(File, physgrid, dims(1:3), &
686 0 : gdims(1:nhdims+1), absnxt_3d(:,:,i,:), absnxt_desc(i))
687 : end do
688 :
689 : ! module data was allocated in radiation_define_restart
690 0 : deallocate(abstot_desc)
691 : end if
692 :
693 1536 : end subroutine radiation_write_restart
694 :
695 : !===============================================================================
696 :
697 768 : subroutine radiation_read_restart(file)
698 :
699 : ! read variables from restart file
700 :
701 : ! arguments
702 : type(file_desc_t), intent(inout) :: file
703 :
704 : ! local variables
705 :
706 : integer :: err_handling
707 : integer :: ierr
708 : integer :: physgrid
709 : integer :: dims(3), gdims(3), nhdims
710 : integer :: vsize
711 : integer :: i
712 : real(r8) :: temp_var
713 :
714 : type(var_desc_t) :: vardesc
715 : character(len=16) :: pname
716 : !----------------------------------------------------------------------------
717 :
718 : ! Put this call here for now. It should move to an init method when the
719 : ! initialization and restart sequencing is unified.
720 768 : call initialize_radbuffer()
721 :
722 768 : if ( radiation_do('aeres') ) then
723 :
724 768 : call pio_seterrorhandling(File, PIO_BCAST_ERROR, err_handling)
725 768 : ierr = pio_inq_varid(File, 'Emissivity', vardesc)
726 768 : call pio_seterrorhandling(File, err_handling)
727 768 : if (ierr /= PIO_NOERR) then
728 768 : if (masterproc) write(iulog,*) 'Warning: Emissivity variable not found on restart file.'
729 768 : return
730 : end if
731 :
732 0 : physgrid = cam_grid_id('physgrid')
733 0 : call cam_grid_dimensions(physgrid, gdims(1:2), nhdims)
734 :
735 0 : dims(1) = pcols
736 0 : dims(2) = pverp
737 0 : dims(3) = endchunk - begchunk + 1
738 0 : gdims(nhdims+1) = dims(2)
739 :
740 : call cam_grid_read_dist_array(File, physgrid, dims(1:3), &
741 0 : gdims(1:nhdims+1), emstot_3d, vardesc)
742 :
743 0 : vsize = pverp - ntoplw + 1
744 0 : dims(2) = vsize
745 0 : gdims(nhdims+1) = dims(2)
746 :
747 0 : do i = ntoplw, pverp
748 0 : write(pname,'(a,i3.3)') 'NAL_absorp', i
749 0 : ierr = pio_inq_varid(File, trim(pname), vardesc)
750 : call cam_grid_read_dist_array(File, physgrid, dims(1:3), &
751 0 : gdims(1:nhdims+1), abstot_3d(:,:,i,:), vardesc)
752 : end do
753 :
754 0 : dims(2) = pver
755 0 : gdims(nhdims+1) = dims(2)
756 0 : do i = 1, 4
757 0 : write(pname,'(a,i3.3)') 'NN_absorp', i
758 0 : ierr = pio_inq_varid(File, trim(pname), vardesc)
759 : call cam_grid_read_dist_array(File, physgrid, dims(1:3), &
760 0 : gdims(1:nhdims+1), absnxt_3d(:,:,i,:), vardesc)
761 : end do
762 : end if
763 :
764 0 : ierr = pio_inq_varid(File, 'nextsw_cday', vardesc)
765 0 : ierr = pio_get_var(File, vardesc, temp_var)
766 0 : nextsw_cday = temp_var
767 :
768 : end subroutine radiation_read_restart
769 :
770 : !===============================================================================
771 :
772 5981472 : subroutine radiation_tend( &
773 : state, ptend, pbuf, cam_out, cam_in, net_flx, rd_out)
774 :
775 : !-----------------------------------------------------------------------
776 : ! Driver for radiation computation.
777 : !
778 : ! NOTE: Radiation uses cgs units, so conversions must be done from
779 : ! model fields to radiation fields.
780 : !-----------------------------------------------------------------------
781 :
782 : use physics_buffer, only: physics_buffer_desc, pbuf_get_field, pbuf_old_tim_idx
783 : use phys_grid, only: get_rlat_all_p, get_rlon_all_p
784 : use physics_types, only: physics_state, physics_ptend
785 : use time_manager, only: get_curr_calday
786 : use radheat, only: radheat_tend
787 : use physconst, only: cpair, stebol
788 : use radconstants, only: nlwbands, nswbands
789 : use radsw, only: radcswmx
790 : use radlw, only: radclwmx
791 : use rad_constituents, only: rad_cnst_get_gas, rad_cnst_out
792 : use aer_rad_props, only: aer_rad_props_sw, aer_rad_props_lw
793 : use interpolate_data, only: vertinterp
794 : use radiation_data, only: rad_data_write
795 : use cloud_cover_diags, only: cloud_cover_diags_out
796 : use tropopause, only: tropopause_find_cam, TROP_ALG_HYBSTOB, TROP_ALG_CLIMATE
797 : use orbit, only: zenith
798 :
799 : ! Arguments
800 : type(physics_state), target, intent(in) :: state
801 : type(physics_ptend), intent(out) :: ptend
802 : type(physics_buffer_desc), pointer :: pbuf(:)
803 : type(cam_out_t), intent(inout) :: cam_out
804 : type(cam_in_t), intent(in) :: cam_in
805 : real(r8), intent(out) :: net_flx(pcols)
806 : type(rad_out_t), target, optional, intent(out) :: rd_out
807 :
808 : ! Local variables
809 : type(rad_out_t), pointer :: rd ! allow rd_out to be optional by allocating a local object
810 : ! if the argument is not present
811 :
812 : integer :: i, k
813 : integer :: lchnk, ncol
814 :
815 : logical :: dosw, dolw, doabsems
816 1495368 : integer, pointer :: nmxrgn(:) ! pbuf pointer to Number of maximally overlapped regions
817 1495368 : real(r8),pointer :: pmxrgn(:,:) ! Maximum values of pressure for each
818 : ! maximally overlapped region.
819 : ! 0->pmxrgn(i,1) is range of pressure for
820 : ! 1st region,pmxrgn(i,1)->pmxrgn(i,2) for
821 : ! 2nd region, etc
822 :
823 1495368 : real(r8),pointer :: emis(:,:) ! Cloud longwave emissivity
824 1495368 : real(r8),pointer :: cldtau(:,:) ! Cloud longwave optical depth
825 1495368 : real(r8),pointer :: cicewp(:,:) ! in-cloud cloud ice water path
826 1495368 : real(r8),pointer :: cliqwp(:,:) ! in-cloud cloud liquid water path
827 :
828 : real(r8) :: cltot(pcols) ! Diagnostic total cloud cover
829 : real(r8) :: cllow(pcols) ! " low cloud cover
830 : real(r8) :: clmed(pcols) ! " mid cloud cover
831 : real(r8) :: clhgh(pcols) ! " hgh cloud cover
832 :
833 : real(r8) :: ftem(pcols,pver) ! Temporary workspace for outfld variables
834 :
835 : integer :: itim_old
836 1495368 : real(r8), pointer, dimension(:,:) :: rel ! liquid effective drop radius (microns)
837 1495368 : real(r8), pointer, dimension(:,:) :: rei ! ice effective drop size (microns)
838 1495368 : real(r8), pointer, dimension(:,:) :: cld ! cloud fraction
839 1495368 : real(r8), pointer, dimension(:,:) :: qrs ! shortwave radiative heating rate
840 1495368 : real(r8), pointer, dimension(:,:) :: qrl ! longwave radiative heating rate
841 :
842 : real(r8) :: calday ! current calendar day
843 : real(r8) :: clat(pcols) ! current latitudes(radians)
844 : real(r8) :: clon(pcols) ! current longitudes(radians)
845 : real(r8) :: coszrs(pcols) ! Cosine solar zenith angle
846 :
847 : real(r8) :: fns(pcols,pverp) ! net shortwave flux
848 : real(r8) :: fcns(pcols,pverp) ! net clear-sky shortwave flux
849 : real(r8) :: fnl(pcols,pverp) ! net longwave flux
850 : real(r8) :: fcnl(pcols,pverp) ! net clear-sky longwave flux
851 :
852 : ! This is used by the chemistry.
853 1495368 : real(r8), pointer :: fsds(:) ! Surface solar down flux
854 :
855 : ! This is used for the energy checker and the Eulerian dycore.
856 1495368 : real(r8), pointer :: fsns(:) ! Surface solar absorbed flux
857 1495368 : real(r8), pointer :: fsnt(:) ! Net column abs solar flux at model top
858 1495368 : real(r8), pointer :: flns(:) ! Srf longwave cooling (up-down) flux
859 1495368 : real(r8), pointer :: flnt(:) ! Net outgoing lw flux at model top
860 :
861 : real(r8) :: pbr(pcols,pver) ! Model mid-level pressures (dynes/cm2)
862 : real(r8) :: pnm(pcols,pverp) ! Model interface pressures (dynes/cm2)
863 : real(r8) :: eccf ! Earth/sun distance factor
864 : real(r8) :: lwupcgs(pcols) ! Upward longwave flux in cgs units
865 :
866 1495368 : real(r8), pointer, dimension(:,:) :: n2o ! nitrous oxide mass mixing ratio
867 1495368 : real(r8), pointer, dimension(:,:) :: ch4 ! methane mass mixing ratio
868 1495368 : real(r8), pointer, dimension(:,:) :: cfc11 ! cfc11 mass mixing ratio
869 1495368 : real(r8), pointer, dimension(:,:) :: cfc12 ! cfc12 mass mixing ratio
870 1495368 : real(r8), pointer, dimension(:,:) :: o3 ! Ozone mass mixing ratio
871 1495368 : real(r8), pointer, dimension(:,:) :: o2 ! Oxygen mass mixing ratio
872 : real(r8), dimension(pcols) :: o2_col ! column oxygen mmr
873 1495368 : real(r8), pointer, dimension(:,:) :: co2 ! co2 mass mixing ratio
874 : real(r8), dimension(pcols) :: co2_col_mean ! co2 column mean mmr
875 1495368 : real(r8), pointer, dimension(:,:) :: sp_hum ! specific humidity
876 :
877 : ! Aerosol shortwave radiative properties
878 : real(r8) :: aer_tau (pcols,0:pver,nswbands) ! aerosol extinction optical depth
879 : real(r8) :: aer_tau_w (pcols,0:pver,nswbands) ! aerosol single scattering albedo * tau
880 : real(r8) :: aer_tau_w_g(pcols,0:pver,nswbands) ! aerosol asymmetry parameter * w * tau
881 : real(r8) :: aer_tau_w_f(pcols,0:pver,nswbands) ! aerosol forward scattered fraction * w * tau
882 :
883 : ! Aerosol longwave absorption optical depth
884 : real(r8) :: odap_aer(pcols,pver,nlwbands)
885 :
886 : ! Gathered indicies of day and night columns
887 : ! chunk_column_index = IdxDay(daylight_column_index)
888 : integer :: Nday ! Number of daylight columns
889 : integer :: Nnite ! Number of night columns
890 : integer, dimension(pcols) :: IdxDay ! Indicies of daylight coumns
891 : integer, dimension(pcols) :: IdxNite ! Indicies of night coumns
892 :
893 : character(*), parameter :: name = 'radiation_tend'
894 :
895 : ! tropopause diagnostic
896 : integer :: troplev(pcols)
897 : real(r8):: p_trop(pcols)
898 :
899 : logical :: write_output ! switch for outfld calls
900 : !----------------------------------------------------------------------
901 :
902 1495368 : lchnk = state%lchnk
903 1495368 : ncol = state%ncol
904 :
905 1495368 : calday = get_curr_calday()
906 :
907 1495368 : if (present(rd_out)) then
908 : rd => rd_out
909 : write_output = .false.
910 : else
911 1495368 : allocate(rd)
912 : write_output=.true.
913 : end if
914 :
915 1495368 : itim_old = pbuf_old_tim_idx()
916 5981472 : call pbuf_get_field(pbuf, cld_idx, cld, start=(/1,1,itim_old/), kount=(/pcols,pver,1/) )
917 :
918 1495368 : call pbuf_get_field(pbuf, qrs_idx,qrs)
919 1495368 : call pbuf_get_field(pbuf, qrl_idx,qrl)
920 :
921 1495368 : call pbuf_get_field(pbuf, fsds_idx, fsds)
922 :
923 1495368 : call pbuf_get_field(pbuf, fsns_idx, fsns)
924 1495368 : call pbuf_get_field(pbuf, fsnt_idx, fsnt)
925 1495368 : call pbuf_get_field(pbuf, flns_idx, flns)
926 1495368 : call pbuf_get_field(pbuf, flnt_idx, flnt)
927 :
928 1495368 : call pbuf_get_field(pbuf, rel_idx, rel)
929 1495368 : call pbuf_get_field(pbuf, rei_idx, rei)
930 :
931 : ! For CRM, make cloud equal to input observations:
932 1495368 : if (single_column.and.scm_crm_mode.and.have_cld) then
933 0 : do k = 1,pver
934 0 : cld(:ncol,k)= cldobs(k)
935 : enddo
936 : endif
937 :
938 : ! Cosine solar zenith angle for current time step
939 1495368 : call get_rlat_all_p(lchnk, ncol, clat)
940 1495368 : call get_rlon_all_p(lchnk, ncol, clon)
941 1495368 : if (use_rad_uniform_angle) then
942 0 : call zenith (calday, clat, clon, coszrs, ncol, dt_avg, uniform_angle=rad_uniform_angle)
943 : else
944 1495368 : call zenith (calday, clat, clon, coszrs, ncol, dt_avg)
945 : end if
946 :
947 : ! Gather night/day column indices.
948 1495368 : Nday = 0
949 1495368 : Nnite = 0
950 24969168 : do i = 1, ncol
951 24969168 : if ( coszrs(i) > 0.0_r8 ) then
952 11736900 : Nday = Nday + 1
953 11736900 : IdxDay(Nday) = i
954 : else
955 11736900 : Nnite = Nnite + 1
956 11736900 : IdxNite(Nnite) = i
957 : end if
958 : end do
959 :
960 1495368 : dosw = radiation_do('sw') ! do shortwave heating calc this timestep?
961 1495368 : dolw = radiation_do('lw') ! do longwave heating calc this timestep?
962 :
963 1495368 : doabsems = radiation_do('absems') ! do absorptivity/emissivity calc this timestep?
964 :
965 : ! Get time of next radiation calculation - albedos will need to be
966 : ! calculated by each surface model at this time
967 1495368 : nextsw_cday = radiation_nextsw_cday()
968 :
969 1495368 : if (dosw .or. dolw) then
970 :
971 : ! pbuf cloud properties set in cloud_diagnostics
972 749232 : call pbuf_get_field(pbuf, cicewp_idx, cicewp)
973 749232 : call pbuf_get_field(pbuf, cliqwp_idx, cliqwp)
974 749232 : call pbuf_get_field(pbuf, cldemis_idx, emis)
975 :
976 749232 : call pbuf_get_field(pbuf, cldtau_idx, cldtau)
977 :
978 749232 : call pbuf_get_field(pbuf, pmxrgn_idx, pmxrgn)
979 749232 : call pbuf_get_field(pbuf, nmxrgn_idx, nmxrgn)
980 :
981 : ! For CRM, make cloud liquid water path equal to input observations
982 749232 : if(single_column.and.scm_crm_mode.and.have_clwp)then
983 0 : do k=1,pver
984 0 : cliqwp(:ncol,k) = clwpobs(k)
985 : end do
986 : endif
987 :
988 : ! Get specific humidity
989 749232 : call rad_cnst_get_gas(0,'H2O', state, pbuf, sp_hum)
990 :
991 : ! Get ozone mass mixing ratio.
992 749232 : call rad_cnst_get_gas(0,'O3', state, pbuf, o3)
993 :
994 : ! Get CO2 mass mixing ratio and compute column mean values
995 749232 : call rad_cnst_get_gas(0,'CO2', state, pbuf, co2)
996 749232 : call calc_col_mean(state, co2, co2_col_mean)
997 :
998 : ! construct cgs unit reps of pmid and pint and get "eccf" - earthsundistancefactor
999 749232 : call radinp(ncol, state%pmid, state%pint, pbr, pnm, eccf)
1000 :
1001 : ! Solar radiation computation
1002 :
1003 749232 : if (hist_fld_active('FSNR') .or. hist_fld_active('FLNR')) then
1004 0 : call tropopause_find_cam(state, troplev, tropP=p_trop, primary=TROP_ALG_HYBSTOB, backup=TROP_ALG_CLIMATE)
1005 : endif
1006 :
1007 749232 : if (dosw) then
1008 :
1009 749232 : call t_startf('rad_sw')
1010 :
1011 : ! Get Oxygen mass mixing ratio.
1012 749232 : call rad_cnst_get_gas(0,'O2', state, pbuf, o2)
1013 749232 : call calc_col_mean(state, o2, o2_col)
1014 :
1015 : ! Get aerosol radiative properties.
1016 749232 : call t_startf('aero_optics_sw')
1017 : call aer_rad_props_sw(0, state, pbuf, nnite, idxnite, &
1018 749232 : aer_tau, aer_tau_w, aer_tau_w_g, aer_tau_w_f)
1019 749232 : call t_stopf('aero_optics_sw')
1020 :
1021 : call radcswmx(lchnk, &
1022 : ncol, pnm, pbr, sp_hum, o3, &
1023 : o2_col, cld, cicewp, cliqwp, rel, &
1024 : rei, eccf, coszrs, rd%solin, &
1025 : cam_in%asdir, cam_in%asdif, cam_in%aldir, cam_in%aldif, nmxrgn, &
1026 : pmxrgn, qrs, rd%qrsc, fsnt, rd%fsntc, rd%fsdtoa, &
1027 : rd%fsntoa, rd%fsutoa, rd%fsntoac, rd%fsnirt, rd%fsnrtc, rd%fsnirtsq, &
1028 : fsns, rd%fsnsc, rd%fsdsc, fsds, cam_out%sols, &
1029 : cam_out%soll, cam_out%solsd, cam_out%solld, fns, fcns, &
1030 : Nday, Nnite, IdxDay, IdxNite, co2_col_mean, &
1031 749232 : aer_tau, aer_tau_w, aer_tau_w_g, aer_tau_w_f , rd%liq_icld_vistau, rd%ice_icld_vistau )
1032 :
1033 749232 : call t_stopf('rad_sw')
1034 :
1035 : ! Output net fluxes at 200 mb
1036 749232 : call vertinterp(ncol, pcols, pverp, state%pint, 20000._r8, fcns, rd%fsn200c)
1037 749232 : call vertinterp(ncol, pcols, pverp, state%pint, 20000._r8, fns, rd%fsn200)
1038 749232 : if (hist_fld_active('FSNR')) then
1039 0 : do i = 1,ncol
1040 0 : call vertinterp(1, 1, pverp, state%pint(i,:), p_trop(i), fns(i,:), rd%fsnr(i))
1041 0 : rd%fsnr(i) = rd%fsnr(i)*cgs2mks
1042 : enddo
1043 : else
1044 12736944 : rd%fsnr(:) = 0._r8
1045 : endif
1046 :
1047 : ! Convert units of shortwave fields needed by rest of model from CGS to MKS
1048 :
1049 12510432 : do i=1,ncol
1050 11761200 : rd%solin(i) = rd%solin(i) *cgs2mks
1051 11761200 : fsds(i) = fsds(i) *cgs2mks
1052 11761200 : rd%fsnirt(i) = rd%fsnirt(i) *cgs2mks
1053 11761200 : rd%fsnrtc(i) = rd%fsnrtc(i) *cgs2mks
1054 11761200 : rd%fsnirtsq(i)= rd%fsnirtsq(i)*cgs2mks
1055 11761200 : fsnt(i) = fsnt(i) *cgs2mks
1056 11761200 : rd%fsdtoa(i) = rd%fsdtoa(i) *cgs2mks
1057 11761200 : fsns(i) = fsns(i) *cgs2mks
1058 11761200 : rd%fsntc(i) = rd%fsntc(i) *cgs2mks
1059 11761200 : rd%fsnsc(i) = rd%fsnsc(i) *cgs2mks
1060 11761200 : rd%fsdsc(i) = rd%fsdsc(i) *cgs2mks
1061 11761200 : rd%fsntoa(i) = rd%fsntoa(i) *cgs2mks
1062 11761200 : rd%fsutoa(i) = rd%fsutoa(i) *cgs2mks
1063 11761200 : rd%fsntoac(i) = rd%fsntoac(i) *cgs2mks
1064 11761200 : rd%fsn200(i) = rd%fsn200(i) *cgs2mks
1065 11761200 : rd%fsn200c(i) = rd%fsn200c(i) *cgs2mks
1066 12510432 : rd%swcf(i) = rd%fsntoa(i) - rd%fsntoac(i)
1067 : end do
1068 :
1069 : ! initialize tau_cld_vistau and tau_icld_vistau as fillvalue, they will stay fillvalue for night columns
1070 331909776 : rd%tot_icld_vistau(1:pcols,1:pver) = fillvalue
1071 331909776 : rd%tot_cld_vistau(1:pcols,1:pver) = fillvalue
1072 :
1073 : ! only do calcs for tot_cld_vistau and tot_icld_vistau on daytime columns
1074 6629832 : do i=1,Nday
1075 : ! sum the water and ice optical depths to get total in-cloud cloud optical depth
1076 11761200 : rd%tot_icld_vistau(IdxDay(i),1:pver) = rd%liq_icld_vistau(IdxDay(i),1:pver) + &
1077 170537400 : rd%ice_icld_vistau(IdxDay(i),1:pver)
1078 :
1079 : ! sum wat and ice, multiply by cloud fraction to get grid-box value
1080 : rd%tot_cld_vistau(IdxDay(i),1:pver) = (rd%liq_icld_vistau(IdxDay(i),1:pver) + &
1081 312421032 : rd%ice_icld_vistau(IdxDay(i),1:pver))*cld(IdxDay(i),1:pver)
1082 : end do
1083 :
1084 : ! add fillvalue for night columns
1085 6629832 : do i = 1, Nnite
1086 158776200 : rd%liq_icld_vistau(IdxNite(i),:) = fillvalue
1087 159525432 : rd%ice_icld_vistau(IdxNite(i),:) = fillvalue
1088 : end do
1089 :
1090 749232 : if (write_output) call radiation_output_sw(state, rd, cam_out, fsns, fsnt, fsds, qrs)
1091 :
1092 : end if ! dosw
1093 :
1094 : ! Longwave radiation computation
1095 :
1096 749232 : if (dolw) then
1097 :
1098 749232 : call t_startf("rad_lw")
1099 :
1100 : ! Convert upward longwave flux units to CGS
1101 :
1102 12510432 : do i=1,ncol
1103 11761200 : lwupcgs(i) = cam_in%lwup(i)*1000._r8
1104 11761200 : if (single_column .and. scm_crm_mode .and. have_tg) &
1105 749232 : lwupcgs(i) = 1000*stebol*tground(1)**4
1106 : end do
1107 :
1108 : ! Get gas phase constituents.
1109 749232 : call rad_cnst_get_gas(0,'N2O', state, pbuf, n2o)
1110 749232 : call rad_cnst_get_gas(0,'CH4', state, pbuf, ch4)
1111 749232 : call rad_cnst_get_gas(0,'CFC11', state, pbuf, cfc11)
1112 749232 : call rad_cnst_get_gas(0,'CFC12', state, pbuf, cfc12)
1113 :
1114 : ! absems requires lw absorption optical depth and transmission through aerosols
1115 749232 : call t_startf('aero_optics_lw')
1116 749232 : if (doabsems) call aer_rad_props_lw(0, state, pbuf, odap_aer)
1117 749232 : call t_stopf('aero_optics_lw')
1118 :
1119 : call radclwmx(lchnk, ncol, doabsems, &
1120 : lwupcgs, state%t, sp_hum, o3, pbr, &
1121 : pnm, state%lnpmid, state%lnpint, n2o, ch4, &
1122 : cfc11, cfc12, cld, emis, pmxrgn, &
1123 : nmxrgn, qrl, rd%qrlc, flns, flnt, rd%flnsc, &
1124 : rd%flntc, cam_out%flwds, rd%fldsc, rd%flut, rd%flutc, &
1125 749232 : fnl, fcnl, co2_col_mean, odap_aer)
1126 :
1127 749232 : call t_stopf("rad_lw")
1128 :
1129 : ! Output fluxes at 200 mb
1130 749232 : call vertinterp(ncol, pcols, pverp, state%pint, 20000._r8, fnl, rd%fln200)
1131 749232 : call vertinterp(ncol, pcols, pverp, state%pint, 20000._r8, fcnl, rd%fln200c)
1132 749232 : if (hist_fld_active('FLNR')) then
1133 0 : do i = 1,ncol
1134 0 : call vertinterp(1, 1, pverp, state%pint(i,:), p_trop(i), fnl(i,:), rd%flnr(i))
1135 : enddo
1136 : else
1137 12736944 : rd%flnr(:) = 0._r8
1138 : endif
1139 :
1140 : ! Convert units of longwave fields needed by rest of model from CGS to MKS
1141 :
1142 12510432 : do i = 1, ncol
1143 11761200 : flnt(i) = flnt(i) *cgs2mks
1144 11761200 : rd%flut(i) = rd%flut(i) *cgs2mks
1145 11761200 : rd%flutc(i) = rd%flutc(i) *cgs2mks
1146 11761200 : rd%lwcf(i) = rd%flutc(i) - rd%flut(i)
1147 11761200 : flns(i) = flns(i) *cgs2mks
1148 11761200 : rd%fldsc(i) = rd%fldsc(i) *cgs2mks
1149 11761200 : rd%flntc(i) = rd%flntc(i) *cgs2mks
1150 11761200 : rd%fln200(i) = rd%fln200(i) *cgs2mks
1151 11761200 : rd%fln200c(i) = rd%fln200c(i) *cgs2mks
1152 11761200 : rd%flnsc(i) = rd%flnsc(i) *cgs2mks
1153 11761200 : cam_out%flwds(i) = cam_out%flwds(i) *cgs2mks
1154 12510432 : rd%flnr(i) = rd%flnr(i) *cgs2mks
1155 : end do
1156 :
1157 749232 : if (write_output) call radiation_output_lw(state, rd, cam_out, flns, flnt, qrl)
1158 :
1159 : end if ! dolw
1160 :
1161 : ! Output aerosol mmr
1162 749232 : if (write_output) call rad_cnst_out(0, state, pbuf)
1163 :
1164 : ! Cloud cover diagnostics
1165 : ! radsw can change pmxrgn and nmxrgn so cldsav needs to follow radsw
1166 749232 : if (write_output) call cloud_cover_diags_out(lchnk, ncol, cld, state%pmid, nmxrgn, pmxrgn )
1167 :
1168 : else ! if (dosw .or. dolw) then
1169 :
1170 : ! convert radiative heating rates from Q*dp to Q for energy conservation
1171 20145672 : do k =1 , pver
1172 324673272 : do i = 1, ncol
1173 304527600 : qrs(i,k) = qrs(i,k)/state%pdel(i,k)
1174 323927136 : qrl(i,k) = qrl(i,k)/state%pdel(i,k)
1175 : end do
1176 : end do
1177 :
1178 : end if ! if (dosw .or. dolw) then
1179 :
1180 : ! output rad inputs and resulting heating rates
1181 1495368 : call rad_data_write( pbuf, state, cam_in, coszrs )
1182 :
1183 : ! Compute net radiative heating tendency
1184 : call radheat_tend(state, pbuf, ptend, qrl, qrs, fsns, &
1185 1495368 : fsnt, flns, flnt, cam_in%asdir, net_flx)
1186 :
1187 1495368 : if (write_output) then
1188 : ! Compute heating rate for dtheta/dt
1189 40374936 : do k=1,pver
1190 650693736 : do i=1,ncol
1191 649198368 : ftem(i,k) = (qrs(i,k) + qrl(i,k))/cpair * (1.e5_r8/state%pmid(i,k))**cappa
1192 : end do
1193 : end do
1194 1495368 : call outfld('HR ',ftem ,pcols ,lchnk )
1195 : end if
1196 :
1197 : ! convert radiative heating rates to Q*dp for energy conservation
1198 40374936 : do k =1 , pver
1199 650693736 : do i = 1, ncol
1200 610318800 : qrs(i,k) = qrs(i,k)*state%pdel(i,k)
1201 649198368 : qrl(i,k) = qrl(i,k)*state%pdel(i,k)
1202 : end do
1203 : end do
1204 :
1205 24969168 : cam_out%netsw(:ncol) = fsns(:ncol)
1206 :
1207 1495368 : if (.not. present(rd_out)) then
1208 1495368 : deallocate(rd)
1209 : end if
1210 2990736 : end subroutine radiation_tend
1211 :
1212 : !===============================================================================
1213 :
1214 749232 : subroutine radiation_output_sw(state, rd, cam_out, fsns, fsnt, fsds, qrs)
1215 :
1216 : ! Dump shortwave radiation information to history buffer (diagnostics)
1217 :
1218 : type(physics_state), intent(in) :: state
1219 : type(rad_out_t), intent(in) :: rd
1220 : type(cam_out_t), intent(in) :: cam_out
1221 : real(r8), intent(in) :: fsns(pcols) ! Surface solar absorbed flux
1222 : real(r8), intent(in) :: fsnt(pcols) ! Net column abs solar flux at model top
1223 : real(r8), intent(in) :: fsds(pcols) ! Surface solar down flux
1224 : real(r8), pointer :: qrs(:,:) ! shortwave radiative heating rate
1225 :
1226 : ! Local variables
1227 : integer :: lchnk, ncol
1228 : real(r8) :: ftem(pcols,pver) ! Temporary workspace for outfld variables
1229 : !----------------------------------------------------------------------
1230 :
1231 749232 : lchnk = state%lchnk
1232 749232 : ncol = state%ncol
1233 :
1234 326020464 : ftem(:ncol,:pver) = qrs(:ncol,:pver)/cpair
1235 749232 : call outfld('QRS ',ftem ,pcols,lchnk)
1236 :
1237 326020464 : ftem(:ncol,:pver) = rd%qrsc(:ncol,:pver)/cpair
1238 749232 : call outfld('QRSC ',ftem ,pcols,lchnk)
1239 :
1240 749232 : call outfld('SOLIN ',rd%solin ,pcols,lchnk)
1241 749232 : call outfld('FSDS ',fsds ,pcols,lchnk)
1242 749232 : call outfld('FSNIRTOA',rd%fsnirt ,pcols,lchnk)
1243 749232 : call outfld('FSNRTOAC',rd%fsnrtc ,pcols,lchnk)
1244 749232 : call outfld('FSNRTOAS',rd%fsnirtsq ,pcols,lchnk)
1245 749232 : call outfld('FSNT ',fsnt ,pcols,lchnk)
1246 749232 : call outfld('FSDTOA ',rd%fsdtoa ,pcols,lchnk)
1247 749232 : call outfld('FSNS ',fsns ,pcols,lchnk)
1248 749232 : call outfld('FSNTC ',rd%fsntc ,pcols,lchnk)
1249 749232 : call outfld('FSNSC ',rd%fsnsc ,pcols,lchnk)
1250 749232 : call outfld('FSDSC ',rd%fsdsc ,pcols,lchnk)
1251 749232 : call outfld('FSNTOA ',rd%fsntoa ,pcols,lchnk)
1252 749232 : call outfld('FSUTOA ',rd%fsutoa ,pcols,lchnk)
1253 749232 : call outfld('FSNTOAC ',rd%fsntoac ,pcols,lchnk)
1254 749232 : call outfld('SOLS ',cam_out%sols ,pcols,lchnk)
1255 749232 : call outfld('SOLL ',cam_out%soll ,pcols,lchnk)
1256 749232 : call outfld('SOLSD ',cam_out%solsd ,pcols,lchnk)
1257 749232 : call outfld('SOLLD ',cam_out%solld ,pcols,lchnk)
1258 749232 : call outfld('FSN200 ',rd%fsn200 ,pcols,lchnk)
1259 749232 : call outfld('FSN200C ',rd%fsn200c ,pcols,lchnk)
1260 749232 : call outfld('FSNR' ,rd%fsnr ,pcols,lchnk)
1261 749232 : call outfld('SWCF ',rd%swcf ,pcols,lchnk)
1262 :
1263 749232 : call outfld('TOT_CLD_VISTAU ',rd%tot_cld_vistau ,pcols,lchnk)
1264 749232 : call outfld('TOT_ICLD_VISTAU ',rd%tot_icld_vistau ,pcols,lchnk)
1265 749232 : call outfld('LIQ_ICLD_VISTAU ',rd%liq_icld_vistau ,pcols,lchnk)
1266 749232 : call outfld('ICE_ICLD_VISTAU ',rd%ice_icld_vistau ,pcols,lchnk)
1267 :
1268 1495368 : end subroutine radiation_output_sw
1269 :
1270 : !===============================================================================
1271 :
1272 749232 : subroutine radiation_output_lw(state, rd, cam_out, flns, flnt, qrl)
1273 :
1274 : ! Dump longwave radiation information to history tape buffer (diagnostics)
1275 :
1276 : type(physics_state), intent(in) :: state
1277 : type(rad_out_t), intent(in) :: rd
1278 : type(cam_out_t), intent(in) :: cam_out
1279 : real(r8), intent(in) :: flns(pcols) ! Srf longwave cooling (up-down) flux
1280 : real(r8), intent(in) :: flnt(pcols) ! Net outgoing lw flux at model top
1281 : real(r8), pointer :: qrl(:,:) ! longwave radiative heating rate
1282 :
1283 : ! Local variables
1284 : integer :: lchnk, ncol
1285 : real(r8) :: ftem(pcols,pver)
1286 : !----------------------------------------------------------------------
1287 :
1288 749232 : lchnk = state%lchnk
1289 749232 : ncol = state%ncol
1290 :
1291 326020464 : call outfld('QRL ',qrl(:ncol,:)/cpair,ncol,lchnk)
1292 326020464 : call outfld('QRLC ',rd%qrlc(:ncol,:)/cpair,ncol,lchnk)
1293 749232 : call outfld('FLNT ',flnt ,pcols,lchnk)
1294 749232 : call outfld('FLUT ',rd%flut ,pcols,lchnk)
1295 749232 : call outfld('FLUTC ',rd%flutc ,pcols,lchnk)
1296 749232 : call outfld('FLNTC ',rd%flntc ,pcols,lchnk)
1297 749232 : call outfld('FLNS ',flns ,pcols,lchnk)
1298 749232 : call outfld('FLDS ',cam_out%flwds ,pcols,lchnk)
1299 749232 : call outfld('FLNSC ',rd%flnsc ,pcols,lchnk)
1300 749232 : call outfld('FLDSC ',rd%fldsc ,pcols,lchnk)
1301 749232 : call outfld('LWCF ',rd%lwcf ,pcols,lchnk)
1302 749232 : call outfld('FLN200 ',rd%fln200,pcols,lchnk)
1303 749232 : call outfld('FLN200C ',rd%fln200c,pcols,lchnk)
1304 749232 : call outfld('FLNR ' ,rd%flnr,pcols,lchnk)
1305 :
1306 749232 : end subroutine radiation_output_lw
1307 :
1308 : !===============================================================================
1309 :
1310 749232 : subroutine radinp(ncol, pmid, pint, pmidrd, pintrd, eccf)
1311 :
1312 : use shr_orb_mod
1313 : use time_manager, only: get_curr_calday
1314 :
1315 : !------------------------------Arguments--------------------------------
1316 : integer, intent(in) :: ncol ! number of atmospheric columns
1317 :
1318 : real(r8), intent(in) :: pmid(pcols,pver) ! Pressure at model mid-levels (pascals)
1319 : real(r8), intent(in) :: pint(pcols,pverp) ! Pressure at model interfaces (pascals)
1320 :
1321 : real(r8), intent(out) :: pmidrd(pcols,pver) ! Pressure at mid-levels (dynes/cm*2)
1322 : real(r8), intent(out) :: pintrd(pcols,pverp) ! Pressure at interfaces (dynes/cm*2)
1323 : real(r8), intent(out) :: eccf ! Earth-sun distance factor
1324 :
1325 : !---------------------------Local variables-----------------------------
1326 : integer :: i, k
1327 : real(r8) :: calday ! current calendar day
1328 : real(r8) :: delta ! Solar declination angle
1329 : !-----------------------------------------------------------------------
1330 :
1331 749232 : calday = get_curr_calday()
1332 : call shr_orb_decl (calday ,eccen ,mvelpp ,lambm0 ,obliqr , &
1333 749232 : delta ,eccf)
1334 :
1335 : ! Convert pressure from pascals to dynes/cm2
1336 20229264 : do k=1,pver
1337 326020464 : do i=1,ncol
1338 305791200 : pmidrd(i,k) = pmid(i,k)*10.0_r8
1339 325271232 : pintrd(i,k) = pint(i,k)*10.0_r8
1340 : end do
1341 : end do
1342 12510432 : do i=1,ncol
1343 12510432 : pintrd(i,pverp) = pint(i,pverp)*10.0_r8
1344 : end do
1345 :
1346 749232 : end subroutine radinp
1347 :
1348 : !===============================================================================
1349 :
1350 1498464 : subroutine calc_col_mean(state, mmr_pointer, mean_value)
1351 :
1352 : ! Compute the column mean.
1353 :
1354 749232 : use cam_logfile, only: iulog
1355 :
1356 : type(physics_state), intent(in) :: state
1357 : real(r8), dimension(:,:), pointer :: mmr_pointer ! mass mixing ratio (lev)
1358 : real(r8), dimension(pcols), intent(out) :: mean_value ! column mean mmr
1359 :
1360 : integer :: i, k, ncol
1361 : real(r8) :: ptot(pcols)
1362 : !-----------------------------------------------------------------------
1363 :
1364 1498464 : ncol = state%ncol
1365 1498464 : mean_value = 0.0_r8
1366 1498464 : ptot = 0.0_r8
1367 :
1368 40458528 : do k=1,pver
1369 652040928 : do i=1,ncol
1370 611582400 : mean_value(i) = mean_value(i) + mmr_pointer(i,k)*state%pdeldry(i,k)
1371 650542464 : ptot(i) = ptot(i) + state%pdeldry(i,k)
1372 : end do
1373 : end do
1374 25020864 : do i=1,ncol
1375 25020864 : mean_value(i) = mean_value(i) / ptot(i)
1376 : end do
1377 :
1378 1498464 : end subroutine calc_col_mean
1379 :
1380 : !===============================================================================
1381 :
1382 0 : end module radiation
1383 :
|
