Line data Source code
1 : module metdata
2 : !-----------------------------------------------------------------------
3 : !
4 : ! BOP
5 : !
6 : ! !MODULE: metdata
7 : !
8 : ! !DESCRIPTION
9 : ! Handles reading and interpolating offline meteorological data which
10 : ! is used to drive the dynamics.
11 : !
12 : ! !USES
13 : use shr_kind_mod, only: r8 => shr_kind_r8, r4 => shr_kind_r4
14 : use shr_cal_mod, only: shr_cal_gregorian
15 : use time_manager, only: get_curr_date, get_step_size, timemgr_is_caltype
16 : use spmd_utils, only: masterproc
17 : use ppgrid, only: pcols, pver, begchunk, endchunk
18 : use time_manager, only: get_curr_calday, get_curr_date, get_step_size
19 : use cam_abortutils, only: endrun
20 : use dynamics_vars, only: T_FVDYCORE_GRID
21 :
22 : #if ( defined SPMD )
23 : use mpishorthand, only: mpicom, mpir8, mpiint,mpichar
24 : use mod_comm, only: mp_sendirr,mp_recvirr
25 : #endif
26 : use perf_mod
27 : use cam_logfile, only: iulog
28 : use pio, only: file_desc_t, pio_put_att, pio_global, pio_get_att, pio_inq_att, &
29 : pio_inq_dimid, pio_inq_dimlen, pio_closefile, pio_get_var, pio_inq_varid, &
30 : pio_offset_kind
31 : use cam_pio_utils, only: cam_pio_openfile
32 :
33 :
34 : implicit none
35 :
36 : private ! all unless made public
37 : save
38 :
39 : ! !PUBLIC MEMBERS
40 :
41 : public :: metdata_dyn_init ! subroutine to open files, allocate blocked arrays, etc
42 : public :: metdata_phys_init ! subroutine to allocate chunked arrays
43 : public :: advance_met ! subroutine to read more data and interpolate
44 : public :: get_met_fields ! interface to set meteorology fields
45 : public :: get_met_srf1
46 : public :: get_met_srf2
47 : public :: get_us_vs
48 : public :: metdata_readnl
49 : public :: met_winds_on_walls
50 : public :: write_met_restart
51 : public :: read_met_restart
52 : public :: met_rlx
53 : public :: met_fix_mass
54 : public :: met_srf_feedback
55 :
56 : interface write_met_restart
57 : Module procedure write_met_restart_bin
58 : Module procedure write_met_restart_pio
59 : end interface
60 :
61 : interface read_met_restart
62 : Module procedure read_met_restart_bin
63 : Module procedure read_met_restart_pio
64 : end interface
65 :
66 :
67 : !------------------------------------------------------------------
68 : ! Interface to access the meteorology fields. Possible invocations
69 : ! are as follows:
70 : ! call get_met_fields( physics_state, us, vs , tend, dt )
71 : ! call get_met_fields( u, v )
72 : ! call get_met_fields( cam_in_t )
73 : !------------------------------------------------------------------
74 : Interface get_met_fields ! overload accessors
75 : Module Procedure get_dyn_flds
76 : Module Procedure get_uv_centered
77 : Module Procedure get_ps
78 : Module Procedure get_ocn_ice_frcs
79 : End Interface
80 :
81 : real(r8), allocatable :: met_ps_next(:,:) ! PS interpolated to next timestep
82 : real(r8), allocatable :: met_ps_curr(:,:) ! PS interpolated to next timestep
83 :
84 : logical :: met_cell_wall_winds = .false. ! true => met data winds are defined on model grid cell walls
85 : logical :: met_remove_file = .false. ! delete metdata file when finished with it
86 :
87 : character(len=16) :: met_shflx_name = 'SHFLX'
88 : character(len=16) :: met_qflx_name = 'QFLX'
89 : real(r8) :: met_snowh_factor = 1._r8
90 : real(r8) :: met_shflx_factor = 1._r8
91 : real(r8) :: met_qflx_factor = 1._r8
92 : logical :: met_srf_feedback = .true.
93 : logical :: met_srf_nudge_flux = .true. ! wsx, wsy, shf, and cflx nudged rather than forced.
94 : ! This is done primarily to prevent unrealistic
95 : ! surface temperatures.
96 :
97 : logical :: met_srf_land = .true. ! nudge surface fields over land (if false ignores
98 : ! all surface nudging for gridboxes with LANDFRAC=1)
99 : logical :: met_srf_land_scale = .false. ! when met_srf_land is false, nudges proportional
100 : ! to the non land fraction, rather than just where
101 : ! LANDFRAC=1
102 : logical :: met_srf_rad = .false. ! nudge albedo and lwup?
103 : logical :: met_srf_refs = .false. ! nudge 2m Q and T and 10m wind
104 : logical :: met_srf_sst = .false. ! nudge sea surface temperature
105 : logical :: met_srf_tau = .true. ! nudge taux and tauy
106 : logical :: met_nudge_temp = .true. ! nudge atmospheric temperature
107 :
108 :
109 : ! radiation/albedo surface field fill value (where there is no sunlight) read in from input data file
110 : real(r8) :: srf_fill_value
111 :
112 : ! !REVISION HISTORY:
113 : ! 31 Oct 2003 Francis Vitt Creation
114 : ! 05 Feb 2004 F Vitt Removed reading/inperpolating PS for current timestep
115 : ! -- only met_ps_next is needed
116 : ! 10 Nov 2004 F Vitt Implemented ability to read from series of files
117 : ! 16 Dec 2004 F Vitt Added offline_met_defaultopts and offline_met_setopts
118 : ! 14 Jul 2005 W Sawyer Removed pmgrid, spmd_dyn dependencies
119 : ! 12 Apr 2006 W Sawyer Removed unneeded ghosting of met_us, met_vs
120 : ! 08 Apr 2010 J Edwards Replaced serial netcdf calls with pio interface
121 : !
122 : ! EOP
123 : !-----------------------------------------------------------------------
124 : ! $Id$
125 : ! $Author$
126 : !-----------------------------------------------------------------------
127 :
128 : type input2d
129 : real(r8), dimension(:,:), pointer :: data => null()
130 : endtype input2d
131 :
132 : type input3d
133 : real(r8), dimension(:,:,:), pointer :: data => null()
134 : endtype input3d
135 :
136 : real(r8), allocatable :: met_t(:,:,:) ! interpolated temperature
137 : real(r8), allocatable :: met_u(:,:,:) ! interpolated zonal wind
138 : real(r8), allocatable :: met_v(:,:,:) ! interpolated meridional wind
139 : real(r8), allocatable :: met_us(:,:,:) ! interpolated zonal wind -staggered
140 : real(r8), allocatable :: met_vs(:,:,:) ! interpolated meridional wind -staggered
141 : real(r8), allocatable :: met_q(:,:,:) ! interpolated water vapor
142 :
143 : real(r8), allocatable :: met_lhflx(:,:)! interpolated latent heat flux
144 : real(r8), allocatable :: met_shflx(:,:)! interpolated sensible heat flux
145 : real(r8), allocatable :: met_qflx(:,:) ! interpolated water vapor flux
146 : real(r8), allocatable :: met_taux(:,:) ! interpolated
147 : real(r8), allocatable :: met_tauy(:,:) ! interpolated
148 : real(r8), allocatable :: met_snowh(:,:) ! interpolated snow height
149 :
150 : real(r8), allocatable :: met_ts(:,:) ! interpolated
151 :
152 : real(r8), allocatable :: met_asdir(:,:) ! interpolated VIS direct albedo
153 : real(r8), allocatable :: met_asdif(:,:) ! interpolated VIS diffuse albedo
154 : real(r8), allocatable :: met_aldir(:,:) ! interpolated NIR direct albedo
155 : real(r8), allocatable :: met_aldif(:,:) ! interpolated NIR diffuse albedo
156 : real(r8), allocatable :: met_lwup(:,:) ! interpolated upwelling LW flux
157 : real(r8), allocatable :: met_sst(:,:) ! interpolated sea surface temperature
158 : real(r8), allocatable :: met_icefrac(:,:) ! interpolated ice fraction
159 : real(r8), allocatable :: met_qref(:,:) ! interpolated reference (2m) specific humidity
160 : real(r8), allocatable :: met_tref(:,:) ! interpolated reference (2m) temperature
161 : real(r8), allocatable :: met_u10(:,:) ! interpolated 10m wind speed
162 :
163 : type(input3d) :: met_ti(2)
164 : type(input3d) :: met_ui(2)
165 : type(input3d) :: met_vi(2)
166 : type(input3d) :: met_usi(2)
167 : type(input3d) :: met_vsi(2)
168 : type(input3d) :: met_qi(2)
169 :
170 : type(input2d) :: met_psi_next(2)
171 : type(input2d) :: met_psi_curr(2)
172 : type(input2d) :: met_lhflxi(2)
173 : type(input2d) :: met_shflxi(2)
174 : type(input2d) :: met_qflxi(2)
175 : type(input2d) :: met_tauxi(2)
176 : type(input2d) :: met_tauyi(2)
177 : type(input2d) :: met_tsi(2)
178 : type(input2d) :: met_snowhi(2)
179 : type(input2d) :: met_asdiri(2)
180 : type(input2d) :: met_asdifi(2)
181 : type(input2d) :: met_aldiri(2)
182 : type(input2d) :: met_aldifi(2)
183 : type(input2d) :: met_lwupi(2)
184 : type(input2d) :: met_ssti(2)
185 : type(input2d) :: met_icefraci(2)
186 : type(input2d) :: met_qrefi(2)
187 : type(input2d) :: met_trefi(2)
188 : type(input2d) :: met_u10i(2)
189 :
190 : integer :: dateid ! var id of the date in the netCDF
191 : integer :: secid ! var id of the sec data
192 : real(r8) :: datatimem = -1.e36_r8 ! time of prv. values read in
193 : real(r8) :: datatimep = -1.e36_r8 ! time of nxt. values read in
194 : real(r8) :: datatimemn = -1.e36_r8 ! time of prv. values read in for next timestep
195 : real(r8) :: datatimepn = -1.e36_r8 ! time of nxt. values read in for next timestep
196 :
197 : integer, parameter :: nm=1 ! array index for previous (minus) data
198 : integer, parameter :: np=2 ! array indes for next (plus) data
199 :
200 : real(r8) :: curr_mod_time ! model time - calendar day
201 : real(r8) :: next_mod_time ! model time - calendar day - next time step
202 :
203 : character(len=256) :: curr_filename, next_filename, met_data_file
204 : character(len=256) :: met_filenames_list = ''
205 : character(len=256) :: met_data_path = ''
206 : type(file_desc_t) :: curr_fileid, next_fileid ! the id of the NetCDF file
207 : real(r8), pointer, dimension(:) :: curr_data_times => null()
208 : real(r8), pointer, dimension(:) :: next_data_times => null()
209 :
210 : real(r8) :: alpha = 1.0_r8 ! don't read in water vapor
211 : ! real(r8), private :: alpha = 0.0 ! read in water vaper each time step
212 :
213 : real(r8), parameter :: D0_0 = 0.0_r8
214 : real(r8), parameter :: D0_5 = 0.5_r8
215 : real(r8), parameter :: D0_75 = 0.75_r8
216 : real(r8), parameter :: D1_0 = 1.0_r8
217 : real(r8), parameter :: days_per_month = 30.6_r8
218 : real(r8), parameter :: days_per_non_leapyear = 365.0_r8
219 : real(r8), parameter :: days_per_year = 365.25_r8
220 : real(r8), parameter :: seconds_per_day = 86400.0_r8
221 : real(r8), parameter :: fill_value = -9999.0_r8
222 :
223 : logical :: online_test = .false.
224 : logical :: debug = .false.
225 :
226 : real(r8) :: met_rlx(pver) = 0._r8
227 : integer :: met_levels
228 : integer :: num_met_levels
229 :
230 : real(r8) :: met_rlx_top = 60._r8
231 : real(r8) :: met_rlx_bot = 50._r8
232 :
233 : real(r8) :: met_rlx_bot_top = 0._r8
234 : real(r8) :: met_rlx_bot_bot = 0._r8
235 :
236 : real(r8) :: met_rlx_time = 0._r8
237 :
238 : #if ( defined OFFLINE_DYN )
239 : logical :: met_fix_mass = .true.
240 : #else
241 : logical :: met_fix_mass = .false.
242 : #endif
243 : logical :: has_ts = .false.
244 : logical :: has_lhflx = .false. ! Is LHFLX present in the met file?
245 :
246 : contains
247 :
248 : !-------------------------------------------------------------------------
249 : ! meteorology data options
250 : !-------------------------------------------------------------------------
251 0 : subroutine metdata_readnl(nlfile)
252 :
253 : use namelist_utils, only: find_group_name
254 : use units, only: getunit, freeunit
255 : use mpishorthand
256 :
257 : character(len=*), intent(in) :: nlfile ! filepath for file containing namelist input
258 :
259 : ! Local variables
260 : integer :: unitn, ierr
261 : character(len=*), parameter :: subname = 'metdata_readnl'
262 :
263 : namelist /metdata_nl/ &
264 : met_data_file, &
265 : met_data_path, &
266 : met_remove_file, &
267 : met_cell_wall_winds, &
268 : met_filenames_list, &
269 : met_rlx_top, &
270 : met_rlx_bot, &
271 : met_rlx_bot_top, &
272 : met_rlx_bot_bot, &
273 : met_rlx_time, &
274 : met_fix_mass, &
275 : met_shflx_name, &
276 : met_shflx_factor, &
277 : met_snowh_factor, &
278 : met_qflx_name, &
279 : met_qflx_factor, &
280 : met_srf_feedback, &
281 : met_srf_nudge_flux, &
282 : met_srf_land, &
283 : met_srf_land_scale, &
284 : met_srf_rad, &
285 : met_srf_refs, &
286 : met_srf_sst, &
287 : met_srf_tau, &
288 : met_nudge_temp
289 :
290 : ! Read namelist
291 0 : if (masterproc) then
292 0 : unitn = getunit()
293 0 : open( unitn, file=trim(nlfile), status='old' )
294 0 : call find_group_name(unitn, 'metdata_nl', status=ierr)
295 0 : if (ierr == 0) then
296 0 : read(unitn, metdata_nl, iostat=ierr)
297 0 : if (ierr /= 0) then
298 0 : call endrun(subname // ':: ERROR reading namelist')
299 : end if
300 : end if
301 0 : close(unitn)
302 0 : call freeunit(unitn)
303 : end if
304 :
305 : #if ( defined SPMD )
306 :
307 : ! Broadcast namelist variables
308 :
309 0 : call mpibcast (met_data_file ,len(met_data_file) ,mpichar,0,mpicom)
310 0 : call mpibcast (met_data_path ,len(met_data_path) ,mpichar,0,mpicom)
311 0 : call mpibcast (met_remove_file ,1 ,mpilog, 0, mpicom )
312 0 : call mpibcast (met_cell_wall_winds,1 ,mpilog, 0, mpicom )
313 0 : call mpibcast (met_filenames_list ,len(met_filenames_list),mpichar,0,mpicom)
314 0 : call mpibcast (met_rlx_top, 1 ,mpir8, 0, mpicom )
315 0 : call mpibcast (met_rlx_bot, 1 ,mpir8, 0, mpicom )
316 0 : call mpibcast (met_rlx_bot_top, 1 ,mpir8, 0, mpicom )
317 0 : call mpibcast (met_rlx_bot_bot, 1 ,mpir8, 0, mpicom )
318 0 : call mpibcast (met_rlx_time, 1 ,mpir8, 0, mpicom )
319 0 : call mpibcast (met_fix_mass, 1 ,mpilog, 0, mpicom )
320 0 : call mpibcast (met_qflx_name ,len(met_qflx_name), mpichar,0,mpicom)
321 0 : call mpibcast (met_shflx_name ,len(met_shflx_name), mpichar,0,mpicom)
322 0 : call mpibcast (met_qflx_factor ,1, mpir8, 0, mpicom )
323 0 : call mpibcast (met_shflx_factor ,1, mpir8, 0, mpicom )
324 0 : call mpibcast (met_snowh_factor ,1, mpir8, 0, mpicom )
325 0 : call mpibcast (met_srf_feedback ,1 ,mpilog, 0, mpicom )
326 0 : call mpibcast (met_srf_nudge_flux ,1 ,mpilog, 0, mpicom )
327 0 : call mpibcast (met_srf_land ,1 ,mpilog, 0, mpicom )
328 0 : call mpibcast (met_srf_land_scale ,1 ,mpilog, 0, mpicom )
329 0 : call mpibcast (met_srf_rad ,1 ,mpilog, 0, mpicom )
330 0 : call mpibcast (met_srf_refs ,1 ,mpilog, 0, mpicom )
331 0 : call mpibcast (met_srf_sst ,1 ,mpilog, 0, mpicom )
332 0 : call mpibcast (met_srf_tau ,1 ,mpilog, 0, mpicom )
333 0 : call mpibcast (met_nudge_temp ,1 ,mpilog, 0, mpicom )
334 : #endif
335 :
336 0 : if (masterproc) then
337 0 : write(iulog,*)'Time-variant meteorological dataset (met_data_file) is: ', trim(met_data_file)
338 0 : write(iulog,*)'Meteorological data file will be removed (met_remove_file): ', met_remove_file
339 0 : write(iulog,*)'Meteorological winds are on cell walls (met_cell_wall_winds): ', met_cell_wall_winds
340 0 : write(iulog,*)'Meteorological file names list file: ', trim(met_filenames_list)
341 0 : write(iulog,*)'Meteorological relax ramp region top at top is (km): ', met_rlx_top
342 0 : write(iulog,*)'Meteorological relax ramp region bottom at top is (km): ', met_rlx_bot
343 0 : write(iulog,*)'Meteorological relax ramp region top at bottom is (km): ', met_rlx_bot_top
344 0 : write(iulog,*)'Meteorological relax ramp region bottom at bottom is (km): ', met_rlx_bot_bot
345 0 : write(iulog,*)'Meteorological relaxation time (hours): ',met_rlx_time
346 0 : write(iulog,*)'Offline driver mass fixer is trurned on (met_fix_mass): ',met_fix_mass
347 0 : write(iulog,*)'Meteorological shflx field name : ', trim(met_shflx_name)
348 0 : write(iulog,*)'Meteorological shflx multiplication factor : ', met_shflx_factor
349 0 : write(iulog,*)'Meteorological qflx field name : ', trim(met_qflx_name)
350 0 : write(iulog,*)'Meteorological qflx multiplication factor : ', met_qflx_factor
351 0 : write(iulog,*)'Meteorological snowh multiplication factor : ', met_snowh_factor
352 0 : write(iulog,*)'Meteorological allow srf models feedbacks : ', met_srf_feedback
353 0 : write(iulog,*)'Meteorological allow srf land nudging : ', met_srf_land
354 0 : write(iulog,*)'Meteorological scale srf land nudging : ', met_srf_land_scale
355 0 : write(iulog,*)'Meteorological allow srf radiation nudging : ', met_srf_rad
356 0 : write(iulog,*)'Meteorological allow srf reference field nudging : ', met_srf_refs
357 0 : write(iulog,*)'Meteorological allow srf sst nudging : ', met_srf_sst
358 0 : write(iulog,*)'Meteorological allow srf tau nudging : ', met_srf_tau
359 0 : write(iulog,*)'Meteorological allow atm tempature nudging : ',met_nudge_temp
360 : endif
361 :
362 0 : end subroutine metdata_readnl
363 :
364 : !--------------------------------------------------------------------------
365 : ! Opens file, allocates arrays
366 : !--------------------------------------------------------------------------
367 0 : subroutine metdata_dyn_init(grid)
368 : use infnan, only : nan, assignment(=)
369 : use cam_control_mod, only : restart_run
370 : implicit none
371 :
372 : ! !INPUT PARAMETERS:
373 : type (T_FVDYCORE_GRID), intent(in) :: grid
374 :
375 :
376 : integer :: im, km, jfirst, jlast, kfirst, klast
377 : integer :: ng_d, ng_s
378 :
379 0 : im = grid%im
380 0 : km = grid%km
381 0 : jfirst = grid%jfirst
382 0 : jlast = grid%jlast
383 0 : kfirst = grid%kfirst
384 0 : klast = grid%klast
385 0 : ng_d = grid%ng_d
386 0 : ng_s = grid%ng_s
387 :
388 :
389 0 : if (.not. restart_run) then ! initial run or branch run
390 0 : curr_filename = met_data_file
391 0 : next_filename = ''
392 : else
393 : ! restart run
394 : ! curr_filename & next_filename already set by restart_dynamics
395 : endif
396 :
397 0 : call open_met_datafile( curr_filename, curr_fileid, curr_data_times, met_data_path, check_dims=.true., grid=grid)
398 :
399 0 : if ( len_trim(next_filename) > 0 ) &
400 0 : call open_met_datafile( next_filename, next_fileid, next_data_times, met_data_path )
401 :
402 : !
403 : ! allocate space for data arrays ...
404 : !
405 : ! dynamics grid
406 :
407 0 : allocate( met_psi_next(nm)%data(im, jfirst:jlast) )
408 0 : allocate( met_psi_next(np)%data(im, jfirst:jlast) )
409 0 : allocate( met_psi_curr(nm)%data(im, jfirst:jlast) )
410 0 : allocate( met_psi_curr(np)%data(im, jfirst:jlast) )
411 0 : allocate( met_ps_next(im, jfirst:jlast) )
412 0 : allocate( met_ps_curr(im, jfirst:jlast) )
413 :
414 0 : allocate( met_us(im, jfirst-ng_d:jlast+ng_s, kfirst:klast) )
415 0 : allocate( met_vs(im, jfirst-ng_s:jlast+ng_d, kfirst:klast) )
416 :
417 0 : met_us = nan
418 0 : met_vs = nan
419 :
420 0 : if (met_cell_wall_winds) then
421 0 : allocate( met_usi(nm)%data(im, jfirst:jlast, kfirst:klast) )
422 0 : allocate( met_usi(np)%data(im, jfirst:jlast, kfirst:klast) )
423 0 : allocate( met_vsi(nm)%data(im, jfirst:jlast, kfirst:klast) )
424 0 : allocate( met_vsi(np)%data(im, jfirst:jlast, kfirst:klast) )
425 : endif
426 :
427 0 : if (.not. met_cell_wall_winds) then
428 :
429 0 : allocate( met_u(im, jfirst-ng_d:jlast+ng_d, kfirst:klast) )
430 0 : allocate( met_ui(nm)%data(im, jfirst:jlast, kfirst:klast) )
431 0 : allocate( met_ui(np)%data(im, jfirst:jlast, kfirst:klast) )
432 :
433 0 : allocate( met_v(im, jfirst-ng_s:jlast+ng_d, kfirst:klast) )
434 0 : allocate( met_vi(nm)%data(im, jfirst:jlast, kfirst:klast) )
435 0 : allocate( met_vi(np)%data(im, jfirst:jlast, kfirst:klast) )
436 :
437 : endif
438 :
439 0 : end subroutine metdata_dyn_init
440 :
441 : !=================================================================================
442 :
443 0 : subroutine metdata_phys_init
444 : use infnan, only : nan, assignment(=)
445 : use cam_history, only : addfld, horiz_only
446 :
447 0 : call addfld ('MET_RLX', (/ 'lev' /), 'A', ' ', 'Meteorology relax function', gridname='fv_centers')
448 0 : call addfld ('MET_TAUX', horiz_only, 'A', ' ', 'Meteorology taux', gridname='physgrid')
449 0 : call addfld ('MET_TAUY', horiz_only, 'A', ' ', 'Meteorology tauy', gridname='physgrid')
450 0 : call addfld ('MET_LHFX', horiz_only, 'A', ' ', 'Meteorology lhflx', gridname='physgrid')
451 0 : call addfld ('MET_SHFX', horiz_only, 'A', ' ', 'Meteorology shflx', gridname='physgrid')
452 0 : call addfld ('MET_QFLX', horiz_only, 'A', ' ', 'Meteorology qflx', gridname='physgrid')
453 0 : call addfld ('MET_PS', horiz_only, 'A', ' ', 'Meteorology PS', gridname='fv_centers')
454 0 : call addfld ('MET_T', (/ 'lev' /), 'A', ' ', 'Meteorology T', gridname='physgrid')
455 0 : call addfld ('MET_U', (/ 'lev' /), 'A', ' ', 'Meteorology U', gridname='fv_centers')
456 0 : call addfld ('MET_V', (/ 'lev' /), 'A', ' ', 'Meteorology V', gridname='fv_centers')
457 0 : call addfld ('MET_SNOWH', horiz_only, 'A', ' ', 'Meteorology snow height', gridname='physgrid')
458 :
459 0 : call addfld ('MET_TS', horiz_only, 'A', 'K', 'Meteorology TS', gridname='physgrid')
460 0 : call addfld ('MET_OCNFRC', horiz_only, 'A', 'fraction', 'Ocean frac derived from met TS', gridname='physgrid')
461 0 : call addfld ('MET_ICEFRC', horiz_only, 'A', 'fraction', 'Sea ice frac derived from met TS', gridname='physgrid')
462 :
463 0 : call addfld ('MET_ASDIR', horiz_only, 'A', '1', 'Meteorology ASDIR', gridname='physgrid')
464 0 : call addfld ('MET_ASDIF', horiz_only, 'A', '1', 'Meteorology ASDIF', gridname='physgrid')
465 0 : call addfld ('MET_ALDIR', horiz_only, 'A', '1', 'Meteorology ALDIR', gridname='physgrid')
466 0 : call addfld ('MET_ALDIF', horiz_only, 'A', '1', 'Meteorology ALDIF', gridname='physgrid')
467 0 : call addfld ('MET_LWUP', horiz_only, 'A', 'Wm-2', 'Meteorology LWUP', gridname='physgrid')
468 0 : call addfld ('MET_SST', horiz_only, 'A', 'K', 'Meteorology SST', gridname='physgrid')
469 0 : call addfld ('MET_ICEFRAC', horiz_only, 'A', 'K', 'Meteorology ICEFRAC', gridname='physgrid')
470 0 : call addfld ('MET_QREF', horiz_only, 'A', 'kg/kg', 'Meteorology QREF', gridname='physgrid')
471 0 : call addfld ('MET_TREF', horiz_only, 'A', 'K', 'Meteorology TREF', gridname='physgrid')
472 0 : call addfld ('MET_U10', horiz_only, 'A', 'ms-1', 'Meteorology U10', gridname='physgrid')
473 :
474 : ! allocate chunked arrays
475 :
476 0 : allocate( met_ti(nm)%data(pcols,pver,begchunk:endchunk) )
477 0 : allocate( met_ti(np)%data(pcols,pver,begchunk:endchunk) )
478 0 : allocate( met_t(pcols,pver,begchunk:endchunk) )
479 :
480 0 : allocate( met_qi(nm)%data(pcols,pver,begchunk:endchunk) )
481 0 : allocate( met_qi(np)%data(pcols,pver,begchunk:endchunk) )
482 0 : allocate( met_q(pcols,pver,begchunk:endchunk) )
483 :
484 0 : allocate( met_lhflxi(nm)%data(pcols,begchunk:endchunk) )
485 0 : allocate( met_lhflxi(np)%data(pcols,begchunk:endchunk) )
486 0 : allocate( met_lhflx(pcols,begchunk:endchunk) )
487 :
488 0 : allocate( met_shflxi(nm)%data(pcols,begchunk:endchunk) )
489 0 : allocate( met_shflxi(np)%data(pcols,begchunk:endchunk) )
490 0 : allocate( met_shflx(pcols,begchunk:endchunk) )
491 :
492 0 : allocate( met_qflxi(nm)%data(pcols,begchunk:endchunk) )
493 0 : allocate( met_qflxi(np)%data(pcols,begchunk:endchunk) )
494 0 : allocate( met_qflx(pcols,begchunk:endchunk) )
495 :
496 0 : allocate( met_tauxi(nm)%data(pcols,begchunk:endchunk) )
497 0 : allocate( met_tauxi(np)%data(pcols,begchunk:endchunk) )
498 0 : allocate( met_taux(pcols,begchunk:endchunk) )
499 :
500 0 : allocate( met_tauyi(nm)%data(pcols,begchunk:endchunk) )
501 0 : allocate( met_tauyi(np)%data(pcols,begchunk:endchunk) )
502 0 : allocate( met_tauy(pcols,begchunk:endchunk) )
503 :
504 0 : allocate( met_tsi(nm)%data(pcols,begchunk:endchunk) )
505 0 : allocate( met_tsi(np)%data(pcols,begchunk:endchunk) )
506 0 : allocate( met_ts(pcols,begchunk:endchunk) )
507 0 : met_ts(:,:) = nan
508 :
509 0 : if(.not.met_srf_feedback) then
510 0 : allocate( met_snowhi(nm)%data(pcols,begchunk:endchunk) )
511 0 : allocate( met_snowhi(np)%data(pcols,begchunk:endchunk) )
512 0 : allocate( met_snowh(pcols,begchunk:endchunk) )
513 0 : met_snowh(:,:) = nan
514 : endif
515 :
516 0 : if(met_srf_rad) then
517 0 : allocate( met_asdiri(nm)%data(pcols,begchunk:endchunk) )
518 0 : allocate( met_asdiri(np)%data(pcols,begchunk:endchunk) )
519 0 : allocate( met_asdir(pcols,begchunk:endchunk) )
520 :
521 0 : allocate( met_asdifi(nm)%data(pcols,begchunk:endchunk) )
522 0 : allocate( met_asdifi(np)%data(pcols,begchunk:endchunk) )
523 0 : allocate( met_asdif(pcols,begchunk:endchunk) )
524 :
525 0 : allocate( met_aldiri(nm)%data(pcols,begchunk:endchunk) )
526 0 : allocate( met_aldiri(np)%data(pcols,begchunk:endchunk) )
527 0 : allocate( met_aldir(pcols,begchunk:endchunk) )
528 :
529 0 : allocate( met_aldifi(nm)%data(pcols,begchunk:endchunk) )
530 0 : allocate( met_aldifi(np)%data(pcols,begchunk:endchunk) )
531 0 : allocate( met_aldif(pcols,begchunk:endchunk) )
532 :
533 0 : allocate( met_lwupi(nm)%data(pcols,begchunk:endchunk) )
534 0 : allocate( met_lwupi(np)%data(pcols,begchunk:endchunk) )
535 0 : allocate( met_lwup(pcols,begchunk:endchunk) )
536 : end if
537 :
538 0 : if(met_srf_refs) then
539 0 : allocate( met_qrefi(nm)%data(pcols,begchunk:endchunk) )
540 0 : allocate( met_qrefi(np)%data(pcols,begchunk:endchunk) )
541 0 : allocate( met_qref(pcols,begchunk:endchunk) )
542 :
543 0 : allocate( met_trefi(nm)%data(pcols,begchunk:endchunk) )
544 0 : allocate( met_trefi(np)%data(pcols,begchunk:endchunk) )
545 0 : allocate( met_tref(pcols,begchunk:endchunk) )
546 :
547 0 : allocate( met_u10i(nm)%data(pcols,begchunk:endchunk) )
548 0 : allocate( met_u10i(np)%data(pcols,begchunk:endchunk) )
549 0 : allocate( met_u10(pcols,begchunk:endchunk) )
550 : end if
551 :
552 0 : if(met_srf_sst) then
553 0 : allocate( met_ssti(nm)%data(pcols,begchunk:endchunk) )
554 0 : allocate( met_ssti(np)%data(pcols,begchunk:endchunk) )
555 0 : allocate( met_sst(pcols,begchunk:endchunk) )
556 :
557 0 : allocate( met_icefraci(nm)%data(pcols,begchunk:endchunk) )
558 0 : allocate( met_icefraci(np)%data(pcols,begchunk:endchunk) )
559 0 : allocate( met_icefrac(pcols,begchunk:endchunk) )
560 : end if
561 :
562 0 : call set_met_rlx()
563 :
564 : ! initialize phys surface fields...
565 0 : call get_model_time()
566 0 : call check_files()
567 0 : call read_phys_srf_flds()
568 0 : call interp_phys_srf_flds()
569 0 : datatimem = -1.e36_r8
570 0 : datatimep = -1.e36_r8
571 0 : end subroutine metdata_phys_init
572 :
573 :
574 : !-----------------------------------------------------------------------
575 : ! Reads more data if needed and interpolates data to current model time
576 : !-----------------------------------------------------------------------
577 0 : subroutine advance_met(grid)
578 0 : use cam_history, only : outfld
579 : implicit none
580 :
581 : type (T_FVDYCORE_GRID), intent(in) :: grid
582 :
583 0 : real(r8) :: met_rlx_2d(grid%ifirstxy:grid%ilastxy,grid%km)
584 : integer :: i,j,k, idim
585 :
586 0 : call t_startf('MET__advance')
587 :
588 : !
589 : !
590 0 : call get_model_time()
591 :
592 0 : if ( ( curr_mod_time > datatimep ) .or. &
593 : ( next_mod_time > datatimepn ) ) then
594 0 : call check_files()
595 : endif
596 :
597 0 : if ( curr_mod_time > datatimep ) then
598 0 : call read_next_metdata(grid)
599 : end if
600 :
601 0 : if ( next_mod_time > datatimepn ) then
602 0 : call read_next_ps(grid)
603 : end if
604 :
605 : ! need to inperpolate the data, regardless !
606 : ! each mpi tasks needs to interpolate
607 0 : call interpolate_metdata(grid)
608 :
609 0 : call t_stopf('MET__advance')
610 :
611 0 : idim = grid%ilastxy - grid%ifirstxy + 1
612 0 : do j = grid%jfirstxy, grid%jlastxy
613 0 : do k = 1, grid%km
614 0 : do i = grid%ifirstxy, grid%ilastxy
615 0 : met_rlx_2d(i,k) = met_rlx(k)
616 : enddo
617 : enddo
618 0 : call outfld('MET_RLX',met_rlx_2d, idim, j)
619 : enddo
620 0 : end subroutine advance_met
621 :
622 : !-------------------------------------------------------------------
623 : ! Method to get some the meteorology data.
624 : ! Sets the following cam_in_t member fields to the
625 : ! meteorology data :
626 : ! qflx
627 : ! lhflx
628 : ! shflx
629 : ! taux
630 : ! tauy
631 : ! snowh
632 : !-------------------------------------------------------------------
633 0 : subroutine get_met_srf2( cam_in )
634 0 : use camsrfexch, only: cam_in_t
635 : use phys_grid, only: get_ncols_p
636 : use cam_history, only: outfld
637 : use shr_const_mod, only: shr_const_stebol
638 : use physconst, only: latvap, latice
639 :
640 : implicit none
641 :
642 : type(cam_in_t), intent(inout), dimension(begchunk:endchunk) :: cam_in
643 :
644 : integer :: c,ncol,i
645 : real(r8) :: met_rlx_sfc(pcols)
646 : real(r8) :: lcl_rlx(pcols)
647 :
648 0 : do c=begchunk,endchunk
649 0 : ncol = get_ncols_p(c)
650 :
651 : ! Nudge or force the surface fields?
652 0 : if (met_srf_nudge_flux) then
653 0 : met_rlx_sfc(:ncol) = met_rlx(pver)
654 : else
655 0 : met_rlx_sfc(:ncol) = 1._r8
656 : end if
657 :
658 : ! Don't nudge the surface for locations that are entirely land?
659 0 : if (.not. met_srf_land) then
660 :
661 : ! Nudging land and forcing ocean.
662 0 : if (met_srf_land_scale) then
663 0 : met_rlx_sfc(:ncol) = (1._r8 - cam_in(c)%landfrac(:ncol)) * &
664 : met_rlx_sfc(:ncol) + &
665 0 : cam_in(c)%landfrac(:ncol) * met_rlx(pver)
666 : else
667 0 : where(cam_in(c)%landfrac(:ncol) == 1._r8) met_rlx_sfc(:ncol) = 0._r8
668 : end if
669 : end if
670 :
671 0 : if (met_srf_tau) then
672 0 : cam_in(c)%wsx(:ncol) = (1._r8-met_rlx_sfc(:ncol)) * cam_in(c)%wsx(:ncol) + met_rlx_sfc(:ncol) * met_taux(:ncol,c)
673 0 : cam_in(c)%wsy(:ncol) = (1._r8-met_rlx_sfc(:ncol)) * cam_in(c)%wsy(:ncol) + met_rlx_sfc(:ncol) * met_tauy(:ncol,c)
674 : end if
675 :
676 0 : cam_in(c)%shf(:ncol) = (1._r8-met_rlx_sfc(:ncol)) * cam_in(c)%shf(:ncol) + &
677 0 : met_rlx_sfc(:ncol) * (met_shflx(:ncol,c) * met_shflx_factor)
678 : cam_in(c)%cflx(:ncol,1) = (1._r8-met_rlx_sfc(:ncol)) * cam_in(c)%cflx(:ncol,1) + &
679 0 : met_rlx_sfc(:ncol) * (met_qflx(:ncol,c) * met_qflx_factor)
680 :
681 : ! If present, nudge the latent heat; otherwise, make an approximation by scaling the
682 : ! water vapor flux.
683 0 : if (has_lhflx) then
684 : cam_in(c)%lhf(:ncol) = (1._r8-met_rlx_sfc(:ncol)) * cam_in(c)%lhf(:ncol) + &
685 0 : met_rlx_sfc(:ncol) * (met_lhflx(:ncol,c) * met_qflx_factor)
686 : else
687 : cam_in(c)%lhf(:ncol) = (1._r8-met_rlx_sfc(:ncol)) * cam_in(c)%lhf(:ncol) + &
688 0 : met_rlx_sfc(:ncol) * (met_qflx(:ncol,c) * met_qflx_factor * latvap)
689 : end if
690 :
691 0 : if (met_srf_rad) then
692 :
693 : ! There can be fill values in the albedos from the met fields, so use the cam_in value
694 : ! if fill is detected. This could be jumpy if that value gets used, but should be in
695 : ! an area with no downwelling solar. Time interpolate around the terminator could cause
696 : ! problems, but the interpolation provides a non-fill value if either endpoint of the
697 : ! interpolation is not fill.
698 0 : where(met_asdir(:ncol,c) == srf_fill_value)
699 : lcl_rlx(:ncol) = 0._r8
700 : elsewhere
701 : lcl_rlx(:ncol) = met_rlx_sfc(:ncol)
702 : end where
703 0 : cam_in(c)%asdir(:ncol) = (1._r8-lcl_rlx(:ncol)) * cam_in(c)%asdir(:ncol) + lcl_rlx(:ncol) * met_asdir(:ncol,c)
704 :
705 0 : where(met_asdif(:ncol,c) == srf_fill_value)
706 : lcl_rlx(:ncol) = 0._r8
707 : elsewhere
708 : lcl_rlx(:ncol) = met_rlx_sfc(:ncol)
709 : end where
710 0 : cam_in(c)%asdif(:ncol) = (1._r8-lcl_rlx(:ncol)) * cam_in(c)%asdif(:ncol) + lcl_rlx(:ncol) * met_asdif(:ncol,c)
711 :
712 0 : where(met_aldir(:ncol,c) == srf_fill_value)
713 : lcl_rlx(:ncol) = 0._r8
714 : elsewhere
715 : lcl_rlx(:ncol) = met_rlx_sfc(:ncol)
716 : end where
717 0 : cam_in(c)%aldir(:ncol) = (1._r8-lcl_rlx(:ncol)) * cam_in(c)%aldir(:ncol) + lcl_rlx(:ncol) * met_aldir(:ncol,c)
718 :
719 0 : where(met_aldif(:ncol,c) == srf_fill_value)
720 : lcl_rlx(:ncol) = 0._r8
721 : elsewhere
722 : lcl_rlx(:ncol) = met_rlx_sfc(:ncol)
723 : end where
724 0 : cam_in(c)%aldif(:ncol) = (1._r8-lcl_rlx(:ncol)) * cam_in(c)%aldif(:ncol) + lcl_rlx(:ncol) * met_aldif(:ncol,c)
725 :
726 0 : cam_in(c)%lwup(:ncol) = (1._r8-met_rlx_sfc(:ncol)) * cam_in(c)%lwup(:ncol) + met_rlx_sfc(:ncol) * met_lwup(:ncol,c)
727 : end if
728 :
729 0 : if (met_srf_refs) then
730 0 : cam_in(c)%qref(:ncol) = (1._r8-met_rlx_sfc(:ncol)) * cam_in(c)%qref(:ncol) + met_rlx_sfc(:ncol) * met_qref(:ncol,c)
731 0 : cam_in(c)%tref(:ncol) = (1._r8-met_rlx_sfc(:ncol)) * cam_in(c)%tref(:ncol) + met_rlx_sfc(:ncol) * met_tref(:ncol,c)
732 0 : cam_in(c)%u10(:ncol) = (1._r8-met_rlx_sfc(:ncol)) * cam_in(c)%u10(:ncol) + met_rlx_sfc(:ncol) * met_u10(:ncol,c)
733 : end if
734 :
735 0 : if (met_srf_sst) then
736 :
737 : ! Meteorological sst is 0 over 100% land, so use the cam_in value if the meteorology thinks
738 : ! it is land.
739 0 : where(met_sst(:ncol,c) == srf_fill_value)
740 : lcl_rlx(:ncol) = 0._r8
741 : elsewhere
742 : lcl_rlx(:ncol) = met_rlx_sfc(:ncol)
743 : end where
744 0 : cam_in(c)%sst(:ncol) = (1._r8-lcl_rlx(:ncol)) * cam_in(c)%sst(:ncol) + lcl_rlx(:ncol) * met_sst(:ncol,c)
745 :
746 0 : where(met_icefrac(:ncol,c) == srf_fill_value)
747 : lcl_rlx(:ncol) = 0._r8
748 : elsewhere
749 : lcl_rlx(:ncol) = met_rlx_sfc(:ncol)
750 : end where
751 0 : cam_in(c)%icefrac(:ncol) = (1._r8-lcl_rlx(:ncol)) * cam_in(c)%icefrac(:ncol) + lcl_rlx(:ncol) * met_icefrac(:ncol,c)
752 :
753 : end if
754 : end do ! Chunk loop
755 :
756 0 : if (debug) then
757 0 : if (masterproc) then
758 0 : write(iulog,*)'METDATA maxval(met_taux),minval(met_taux): ',maxval(met_taux),minval(met_taux)
759 0 : write(iulog,*)'METDATA maxval(met_tauy),minval(met_tauy): ',maxval(met_tauy),minval(met_tauy)
760 0 : write(iulog,*)'METDATA maxval(met_lhflx),minval(met_lhflx): ',maxval(met_lhflx),minval(met_lhflx)
761 0 : write(iulog,*)'METDATA maxval(met_shflx),minval(met_shflx): ',maxval(met_shflx),minval(met_shflx)
762 0 : write(iulog,*)'METDATA maxval(met_qflx),minval(met_qflx): ',maxval(met_qflx),minval(met_qflx)
763 0 : write(iulog,*)'METDATA maxval(met_asdir),minval(met_asdir): ',maxval(met_asdir),minval(met_asdir)
764 0 : write(iulog,*)'METDATA maxval(met_asdif),minval(met_asdif): ',maxval(met_asdif),minval(met_asdif)
765 0 : write(iulog,*)'METDATA maxval(met_aldir),minval(met_aldir): ',maxval(met_aldir),minval(met_aldir)
766 0 : write(iulog,*)'METDATA maxval(met_aldif),minval(met_aldif): ',maxval(met_aldif),minval(met_aldif)
767 0 : write(iulog,*)'METDATA maxval(met_lwup),minval(met_lwup): ',maxval(met_lwup),minval(met_lwup)
768 0 : write(iulog,*)'METDATA maxval(met_qref),minval(met_qref): ',maxval(met_qref),minval(met_qref)
769 0 : write(iulog,*)'METDATA maxval(met_tref),minval(met_tref): ',maxval(met_tref),minval(met_tref)
770 0 : write(iulog,*)'METDATA maxval(met_u10),minval(met_u10): ',maxval(met_u10),minval(met_u10)
771 0 : write(iulog,*)'METDATA maxval(met_sst),minval(met_sst): ',maxval(met_sst),minval(met_sst)
772 0 : write(iulog,*)'METDATA maxval(met_icefrac),minval(met_icefrac): ',maxval(met_icefrac),minval(met_icefrac)
773 : endif
774 : endif
775 :
776 0 : do c = begchunk, endchunk
777 0 : call outfld('MET_TAUX',cam_in(c)%wsx , pcols ,c )
778 0 : call outfld('MET_TAUY',cam_in(c)%wsy , pcols ,c )
779 0 : call outfld('MET_LHFX',cam_in(c)%lhf , pcols ,c )
780 0 : call outfld('MET_SHFX',cam_in(c)%shf , pcols ,c )
781 0 : call outfld('MET_QFLX',cam_in(c)%cflx(:,1) , pcols ,c )
782 0 : call outfld('MET_ASDIR',cam_in(c)%asdir , pcols ,c )
783 0 : call outfld('MET_ASDIF',cam_in(c)%asdif , pcols ,c )
784 0 : call outfld('MET_ALDIR',cam_in(c)%aldir , pcols ,c )
785 0 : call outfld('MET_ALDIF',cam_in(c)%aldif , pcols ,c )
786 0 : call outfld('MET_LWUP',cam_in(c)%lwup , pcols ,c )
787 0 : call outfld('MET_QREF',cam_in(c)%qref , pcols ,c )
788 0 : call outfld('MET_TREF',cam_in(c)%tref , pcols ,c )
789 0 : call outfld('MET_U10',cam_in(c)%u10 , pcols ,c )
790 0 : call outfld('MET_SST',cam_in(c)%sst , pcols ,c )
791 0 : call outfld('MET_ICEFRAC',cam_in(c)%icefrac , pcols ,c )
792 : enddo
793 :
794 0 : end subroutine get_met_srf2
795 :
796 : !-------------------------------------------------------------------
797 : !-------------------------------------------------------------------
798 0 : subroutine get_met_srf1( cam_in )
799 0 : use camsrfexch, only: cam_in_t
800 : use phys_grid, only: get_ncols_p
801 : use cam_history, only: outfld
802 : use shr_const_mod, only: shr_const_stebol
803 :
804 : implicit none
805 :
806 : type(cam_in_t), intent(inout), dimension(begchunk:endchunk) :: cam_in
807 :
808 : integer :: c,ncol,i
809 :
810 0 : if (met_srf_feedback) return
811 0 : if (.not.has_ts) then
812 0 : call endrun('The meteorolgy input must have TS to run with met_srf_feedback set to FALSE')
813 : endif
814 :
815 0 : do c=begchunk,endchunk
816 0 : ncol = get_ncols_p(c)
817 0 : cam_in(c)%ts(:ncol) = met_ts(:ncol,c)
818 0 : do i = 1,ncol
819 0 : cam_in(c)%snowhland(i) = met_snowh(i,c)*cam_in(c)%landfrac(i) * met_snowh_factor
820 : enddo
821 : end do ! Chunk loop
822 :
823 0 : if (debug) then
824 0 : if (masterproc) then
825 0 : write(iulog,*)'METDATA maxval(met_ts),minval(met_ts): ',maxval(met_ts),minval(met_ts)
826 0 : write(iulog,*)'METDATA maxval(met_snowh),minval(met_snowh): ',maxval(met_snowh),minval(met_snowh)
827 : endif
828 : endif
829 :
830 0 : do c = begchunk, endchunk
831 0 : call outfld('MET_SNOWH',cam_in(c)%snowhland, pcols ,c )
832 : enddo
833 :
834 0 : end subroutine get_met_srf1
835 :
836 : !-------------------------------------------------------------------
837 : !-------------------------------------------------------------------
838 0 : subroutine get_ocn_ice_frcs( lndfrc, ocnfrc, icefrc, lchnk, ncol )
839 :
840 0 : use shr_const_mod, only: SHR_CONST_TKFRZSW
841 : use shr_const_mod, only: SHR_CONST_TKFRZ
842 : use cam_history, only: outfld
843 :
844 : ! args
845 : real(r8), intent( in) :: lndfrc (pcols)
846 : real(r8), intent(out) :: ocnfrc (pcols)
847 : real(r8), intent(out) :: icefrc (pcols)
848 :
849 : integer, intent(in) :: lchnk
850 : integer, intent(in) :: ncol
851 :
852 : ! local vars
853 : integer :: i
854 :
855 0 : if (met_srf_sst) then
856 0 : do i = 1,ncol
857 :
858 : ! If configured for using SST, and ICEFRAC, then get icefrc
859 : ! directly from the meteorological data.
860 0 : if (met_icefrac(i,lchnk) == srf_fill_value) then
861 0 : icefrc(i) = 0._r8
862 : else
863 0 : icefrc(i) = min(met_icefrac(i,lchnk), 1._r8 - lndfrc(i))
864 : end if
865 0 : ocnfrc(i) = 1._r8 - lndfrc(i) - icefrc(i)
866 : enddo
867 : else
868 :
869 0 : if (.not.has_ts) then
870 0 : if (masterproc) then
871 0 : write(iulog,*) 'get_ocn_ice_frcs: TS is not in the met dataset and cannot set ocnfrc and icefrc'
872 0 : call endrun('get_ocn_ice_frcs: TS is not in the met dataset')
873 : endif
874 : endif
875 :
876 0 : do i = 1,ncol
877 :
878 0 : if ( met_ts(i,lchnk) < SHR_CONST_TKFRZ-2._r8 ) then
879 0 : ocnfrc(i) = 0._r8
880 0 : icefrc(i) = 1._r8 - lndfrc(i)
881 : else
882 0 : icefrc(i) = 0._r8
883 0 : ocnfrc(i) = 1._r8 - lndfrc(i)
884 : endif
885 :
886 : enddo
887 : end if
888 :
889 0 : call outfld('MET_TS', met_ts(:ncol,lchnk) , ncol ,lchnk )
890 0 : call outfld('MET_OCNFRC', ocnfrc(:ncol) , ncol ,lchnk )
891 0 : call outfld('MET_ICEFRC', icefrc(:ncol) , ncol ,lchnk )
892 :
893 0 : endsubroutine get_ocn_ice_frcs
894 :
895 : !-------------------------------------------------------------------
896 : ! allows access to physics state fields
897 : ! q : water vapor
898 : ! ps : surface pressure
899 : ! t : temperature
900 : !-------------------------------------------------------------------
901 0 : subroutine get_dyn_flds( state, tend, dt )
902 :
903 0 : use physics_types, only: physics_state, physics_tend, physics_dme_adjust
904 : use ppgrid, only: pcols, pver, begchunk, endchunk
905 : use phys_grid, only: get_ncols_p
906 : use cam_history, only: outfld
907 : use air_composition,only: thermodynamic_active_species_liq_num, thermodynamic_active_species_ice_num
908 : use air_composition,only: thermodynamic_active_species_liq_idx,thermodynamic_active_species_ice_idx
909 :
910 : implicit none
911 :
912 : type(physics_state), intent(inout), dimension(begchunk:endchunk) :: state
913 : type(physics_tend ), intent(inout), dimension(begchunk:endchunk) :: tend
914 : real(r8), intent(in ) :: dt ! model physics timestep
915 :
916 : integer :: lats(pcols) ! array of latitude indices
917 : integer :: lons(pcols) ! array of longitude indices
918 : integer :: c, ncol, i,j,k
919 : integer :: m_cnst,m
920 : real(r8):: qini(pcols,pver) ! initial specific humidity
921 : real(r8):: totliqini(pcols,pver) ! initial total liquid
922 : real(r8):: toticeini(pcols,pver) ! initial total ice
923 :
924 : real(r8) :: tmp(pcols,pver)
925 :
926 0 : call t_startf('MET__GET_DYN2')
927 :
928 0 : do c = begchunk, endchunk
929 0 : ncol = get_ncols_p(c)
930 : !
931 : ! update water variables
932 : !
933 0 : qini(:ncol,:pver) = state(c)%q(:ncol,:pver,1)
934 0 : totliqini = 0.0_r8
935 0 : do m_cnst=1,thermodynamic_active_species_liq_num
936 0 : m = thermodynamic_active_species_liq_idx(m_cnst)
937 0 : totliqini(:ncol,:pver) = totliqini(:ncol,:pver)+state(c)%q(:ncol,:pver,m)
938 : end do
939 0 : toticeini = 0.0_r8
940 0 : do m_cnst=1,thermodynamic_active_species_ice_num
941 0 : m = thermodynamic_active_species_ice_idx(m_cnst)
942 0 : toticeini(:ncol,:pver) = toticeini(:ncol,:pver)+state(c)%q(:ncol,:pver,m)
943 : end do
944 :
945 0 : do k=1,pver
946 0 : do i=1,ncol
947 0 : if (met_nudge_temp) then
948 0 : state(c)%t(i,k) = (1._r8-met_rlx(k))*state(c)%t(i,k) + met_rlx(k)*met_t(i,k,c)
949 : end if
950 : ! at this point tracer mixing ratios have already been
951 : ! converted from dry to moist
952 0 : state(c)%q(i,k,1) = alpha*state(c)%q(i,k,1) + (D1_0-alpha)*met_q(i,k,c)
953 :
954 0 : if ((state(c)%q(i,k,1) < D0_0).and. (alpha .ne. D1_0 )) state(c)%q(i,k,1) = D0_0
955 :
956 : end do
957 :
958 : end do
959 :
960 : ! now adjust mass of each layer now that water vapor has changed
961 0 : if (( .not. online_test ) .and. (alpha .ne. D1_0 )) then
962 : call physics_dme_adjust(state(c), tend(c), qini, totliqini, toticeini, dt)
963 : endif
964 :
965 : end do
966 :
967 0 : if (debug) then
968 0 : if (masterproc) then
969 0 : write(iulog,*)'METDATA maxval(met_t),minval(met_t): ', maxval(met_t),minval(met_t)
970 0 : write(iulog,*)'METDATA maxval(met_ps_next),minval(met_ps_next): ', maxval(met_ps_next),minval(met_ps_next)
971 : endif
972 : endif
973 :
974 0 : do c = begchunk, endchunk
975 0 : call outfld('MET_T ',state(c)%t , pcols ,c )
976 : enddo
977 0 : call t_stopf('MET__GET_DYN2')
978 :
979 0 : end subroutine get_dyn_flds
980 :
981 : !------------------------------------------------------------------------
982 : ! get the meteorological winds on the grid cell centers (A-grid)
983 : !------------------------------------------------------------------------
984 0 : subroutine get_uv_centered( grid, u, v )
985 :
986 0 : use cam_history, only: outfld
987 :
988 : implicit none
989 :
990 : type (T_FVDYCORE_GRID), intent(in) :: grid
991 : real(r8), intent(out) :: u(grid%im, grid%jfirst-grid%ng_d:grid%jlast+grid%ng_d, &
992 : grid%kfirst:grid%klast) ! u-wind on A-grid
993 : real(r8), intent(out) :: v(grid%im, grid%jfirst-grid%ng_s:grid%jlast+grid%ng_d, &
994 : grid%kfirst:grid%klast) ! v-wind on A-grid
995 :
996 : integer :: i,j,k
997 :
998 : integer :: jm, jfirst, jlast, jfirstxy, jlastxy, kfirst, klast, ng_d, ng_s, ifirstxy, ilastxy
999 :
1000 0 : real(r8) :: u3s_tmp(grid%ifirstxy:grid%ilastxy,grid%km)
1001 0 : real(r8) :: v3s_tmp(grid%ifirstxy:grid%ilastxy,grid%km)
1002 :
1003 0 : jm = grid%jm
1004 0 : jfirstxy= grid%jfirstxy
1005 0 : jlastxy = grid%jlastxy
1006 0 : jfirst = grid%jfirst
1007 0 : jlast = grid%jlast
1008 0 : kfirst = grid%kfirst
1009 0 : klast = grid%klast
1010 0 : ifirstxy= grid%ifirstxy
1011 0 : ilastxy = grid%ilastxy
1012 :
1013 0 : ng_d = grid%ng_d
1014 0 : ng_s = grid%ng_s
1015 :
1016 0 : u(:,:,:) = D0_0
1017 0 : v(:,:,:) = D0_0
1018 :
1019 0 : u( :, max(1,jfirst-ng_d):min(jm,jlast+ng_d), kfirst:klast ) = &
1020 0 : met_u(:, max(1,jfirst-ng_d):min(jm,jlast+ng_d), kfirst:klast )
1021 :
1022 0 : v( :, max(1,jfirst-ng_s):min(jm,jlast+ng_d), kfirst:klast ) = &
1023 0 : met_v(:, max(1,jfirst-ng_s):min(jm,jlast+ng_d), kfirst:klast )
1024 :
1025 0 : if (masterproc) then
1026 0 : if (debug) write(iulog,*)'METDATA maxval(u),minval(u),maxval(v),minval(v) : ',&
1027 0 : maxval(u(:, max(1,jfirst-ng_d):min(jm,jlast+ng_d), kfirst:klast )),&
1028 0 : minval(u(:, max(1,jfirst-ng_d):min(jm,jlast+ng_d), kfirst:klast )),&
1029 0 : maxval(v(:, max(1,jfirst-ng_s):min(jm,jlast+ng_d), kfirst:klast )),&
1030 0 : minval(v(:, max(1,jfirst-ng_s):min(jm,jlast+ng_d), kfirst:klast ))
1031 : endif
1032 :
1033 0 : if ( grid%twod_decomp == 0 ) then
1034 0 : do j = jfirst, jlast
1035 0 : do k = kfirst, klast
1036 0 : do i = 1, grid%im
1037 0 : u3s_tmp(i,k) = u(i,j,k)
1038 0 : v3s_tmp(i,k) = v(i,j,k)
1039 : enddo
1040 : enddo
1041 0 : call outfld ('MET_U ', u3s_tmp, grid%im, j )
1042 0 : call outfld ('MET_V ', v3s_tmp, grid%im, j )
1043 : enddo
1044 : endif
1045 :
1046 0 : end subroutine get_uv_centered
1047 :
1048 : !------------------------------------------------------------------------
1049 : ! get the meteorological surface pressure interp to dyn substep
1050 : !------------------------------------------------------------------------
1051 0 : subroutine get_ps( grid, ps, nsubsteps, n )
1052 :
1053 0 : use cam_history, only: outfld
1054 :
1055 : implicit none
1056 :
1057 : type (T_FVDYCORE_GRID), intent(in) :: grid
1058 : real(r8), intent(out) :: ps(grid%im, grid%jfirst:grid%jlast)
1059 : integer, intent(in) :: nsubsteps
1060 : integer, intent(in) :: n
1061 :
1062 : real(r8) :: num1, num2
1063 : integer :: j
1064 :
1065 0 : num1 = n
1066 0 : num2 = nsubsteps
1067 :
1068 0 : ps(:,:) = met_ps_curr(:,:) + num1*(met_ps_next(:,:)-met_ps_curr(:,:))/num2
1069 :
1070 0 : if ( grid%twod_decomp == 0 ) then
1071 0 : do j = grid%jfirst, grid%jlast
1072 0 : call outfld('MET_PS',ps(:,j), grid%im ,j )
1073 : enddo
1074 : endif
1075 0 : end subroutine get_ps
1076 :
1077 : !------------------------------------------------------------------------
1078 : ! get the meteorological winds on the grid cell walls (vorticity winds)
1079 : ! us : staggered zonal wind
1080 : ! vs : staggered meridional wind
1081 : !------------------------------------------------------------------------
1082 0 : subroutine get_us_vs( grid, us, vs )
1083 :
1084 : implicit none
1085 :
1086 : type (T_FVDYCORE_GRID), intent(in) :: grid
1087 : real(r8), intent(inout) :: us(grid%im, grid%jfirst-grid%ng_d:grid%jlast+grid%ng_s, &
1088 : grid%kfirst:grid%klast) ! u-wind on d-grid
1089 : real(r8), intent(inout) :: vs(grid%im, grid%jfirst-grid%ng_s:grid%jlast+grid%ng_d, &
1090 : grid%kfirst:grid%klast) ! v-wind on d-grid
1091 :
1092 : integer :: i,j,k
1093 :
1094 : integer :: jm, jfirst, jlast, kfirst, klast, ng_d, ng_s
1095 :
1096 0 : jm = grid%jm
1097 0 : jfirst = grid%jfirst
1098 0 : jlast = grid%jlast
1099 0 : kfirst = grid%kfirst
1100 0 : klast = grid%klast
1101 0 : ng_d = grid%ng_d
1102 0 : ng_s = grid%ng_s
1103 :
1104 0 : call t_startf('MET__get_us_vs')
1105 :
1106 : ! vertical relaxation (blending) occurs in dyn_run (dyn_comp.F90)
1107 :
1108 0 : us(:,:,:) = 1.e36_r8
1109 0 : vs(:,:,:) = 1.e36_r8
1110 0 : us( :, max(2,jfirst): min(jm,jlast), kfirst:klast) = &
1111 0 : met_us( :, max(2,jfirst): min(jm,jlast), kfirst:klast)
1112 0 : vs( :, max(1,jfirst): min(jm,jlast), kfirst:klast) = &
1113 0 : met_vs( :, max(1,jfirst): min(jm,jlast), kfirst:klast)
1114 0 : if (masterproc) then
1115 0 : if (debug) write(iulog,*)grid%iam,': METDATA maxval(us),minval(us),maxval(vs),minval(vs) : ',&
1116 0 : maxval(us( :, max(2,jfirst): min(jm,jlast), kfirst:klast)),&
1117 0 : minval(us( :, max(2,jfirst): min(jm,jlast), kfirst:klast)),&
1118 0 : maxval(vs( :, max(1,jfirst): min(jm,jlast), kfirst:klast)),&
1119 0 : minval(vs( :, max(1,jfirst): min(jm,jlast), kfirst:klast))
1120 : endif
1121 :
1122 : !!$ if (debug) then
1123 : !!$ u3s_tmp = 1.e36
1124 : !!$ do j = jfirst, jlast
1125 : !!$ do k = kfirst, klast
1126 : !!$ do i = 1, im
1127 : !!$ if (j >= 2) u3s_tmp(i,k) = us(i,j,k)
1128 : !!$ v3s_tmp(i,k) = vs(i,j,k)
1129 : !!$ enddo
1130 : !!$ enddo
1131 : !!$ call outfld ('MET_US ', u3s_tmp, im, j )
1132 : !!$ call outfld ('MET_VS ', v3s_tmp, im, j )
1133 : !!$ enddo
1134 : !!$ endif
1135 : !!$
1136 0 : call t_stopf('MET__get_us_vs')
1137 :
1138 0 : end subroutine get_us_vs
1139 :
1140 : !-------------------------------------------------------------------------
1141 : ! writes file names to restart file
1142 : !-------------------------------------------------------------------------
1143 :
1144 0 : subroutine write_met_restart_pio(File)
1145 : type(file_desc_t), intent(inout) :: File
1146 : integer :: ierr
1147 0 : ierr = pio_put_att(File, PIO_GLOBAL, 'current_metdata_filename', curr_filename)
1148 0 : ierr = pio_put_att(File, PIO_GLOBAL, 'next_metdata_filename', next_filename)
1149 :
1150 0 : end subroutine write_met_restart_pio
1151 0 : subroutine read_met_restart_pio(File)
1152 : type(file_desc_t), intent(inout) :: File
1153 :
1154 : integer :: ierr, xtype
1155 : integer(pio_offset_kind) :: slen
1156 :
1157 0 : ierr = pio_inq_att(File, PIO_GLOBAL, 'current_metdata_filename',xtype, slen)
1158 0 : ierr = pio_get_att(File, PIO_GLOBAL, 'current_metdata_filename', curr_filename)
1159 0 : curr_filename(slen+1:256) = ''
1160 :
1161 0 : ierr = pio_inq_att(File, PIO_GLOBAL, 'next_metdata_filename',xtype, slen)
1162 0 : ierr = pio_get_att(File, PIO_GLOBAL, 'next_metdata_filename', next_filename)
1163 0 : next_filename(slen+1:256) = ''
1164 :
1165 0 : end subroutine read_met_restart_pio
1166 :
1167 0 : subroutine write_met_restart_bin( nrg )
1168 : implicit none
1169 : integer,intent(in) :: nrg ! Unit number
1170 : integer :: ioerr ! error status
1171 :
1172 0 : if (masterproc) then
1173 0 : write(nrg, iostat=ioerr) curr_filename
1174 0 : if (ioerr /= 0 ) then
1175 0 : write(iulog,*) 'WRITE ioerror ',ioerr,' on i/o unit = ',nrg
1176 0 : call endrun ('WRITE_RESTART_DYNAMICS')
1177 : end if
1178 0 : write(nrg, iostat=ioerr) next_filename
1179 0 : if (ioerr /= 0 ) then
1180 0 : write(iulog,*) 'WRITE ioerror ',ioerr,' on i/o unit = ',nrg
1181 0 : call endrun ('WRITE_RESTART_DYNAMICS')
1182 : end if
1183 : end if
1184 0 : end subroutine write_met_restart_bin
1185 :
1186 : !-------------------------------------------------------------------------
1187 : ! reads file names from restart file
1188 : !-------------------------------------------------------------------------
1189 0 : subroutine read_met_restart_bin( nrg )
1190 : implicit none
1191 : integer,intent(in) :: nrg ! Unit number
1192 : integer :: ioerr ! error status
1193 :
1194 0 : if (masterproc) then
1195 0 : read(nrg, iostat=ioerr) curr_filename
1196 0 : if (ioerr /= 0 ) then
1197 0 : write(iulog,*) 'READ ioerror ',ioerr,' on i/o unit = ',nrg
1198 0 : call endrun ('READ_RESTART_DYNAMICS')
1199 : end if
1200 0 : read(nrg, iostat=ioerr) next_filename
1201 0 : if (ioerr /= 0 ) then
1202 0 : write(iulog,*) 'READ ioerror ',ioerr,' on i/o unit = ',nrg
1203 0 : call endrun ('READ_RESTART_DYNAMICS')
1204 : end if
1205 : end if
1206 :
1207 : #if ( defined SPMD )
1208 0 : call mpibcast ( curr_filename ,len(curr_filename) ,mpichar,0,mpicom)
1209 0 : call mpibcast ( next_filename ,len(next_filename) ,mpichar,0,mpicom)
1210 : #endif
1211 0 : end subroutine read_met_restart_bin
1212 :
1213 : !-------------------------------------------------------------------------
1214 : ! returns true if the met winds are defined on cell walls
1215 : !-------------------------------------------------------------------------
1216 0 : function met_winds_on_walls()
1217 : logical :: met_winds_on_walls
1218 :
1219 0 : met_winds_on_walls = met_cell_wall_winds
1220 0 : end function met_winds_on_walls
1221 :
1222 : ! internal methods :
1223 :
1224 : !-------------------------------------------------------------------------
1225 : ! transfers cell-centered winds to cell walls
1226 : !-------------------------------------------------------------------------
1227 0 : subroutine transfer_windsToWalls(grid)
1228 :
1229 : implicit none
1230 :
1231 : type (T_FVDYCORE_GRID), intent(in) :: grid
1232 : integer :: i,j,k
1233 : integer :: im, jfirst, jlast, kfirst, klast
1234 :
1235 0 : im = grid%im
1236 0 : jfirst = grid%jfirst
1237 0 : jlast = grid%jlast
1238 0 : kfirst = grid%kfirst
1239 0 : klast = grid%klast
1240 :
1241 0 : call t_startf('MET__transfer_windsToWalls')
1242 :
1243 : !$omp parallel do private (i, j, k)
1244 0 : do k = kfirst, klast
1245 :
1246 0 : do j = jfirst+1,jlast
1247 0 : do i = 1,im
1248 0 : met_us(i,j,k) = ( met_u(i,j,k) + met_u(i,j-1,k) )*D0_5
1249 : end do
1250 : end do
1251 :
1252 : #if defined( SPMD )
1253 0 : if ( jfirst .gt. 1 ) then
1254 0 : do i = 1, im
1255 : ! met_u is alread ghosted at this point
1256 0 : met_us(i,jfirst,k) = ( met_u(i,jfirst,k) + met_u(i,jfirst-1,k) )*D0_5
1257 : enddo
1258 : endif
1259 : #endif
1260 :
1261 0 : do j = jfirst,jlast
1262 0 : met_vs(1,j,k) = ( met_v(1,j,k) + met_v(im,j,k) )*D0_5
1263 0 : do i = 2,im
1264 0 : met_vs(i,j,k) = ( met_v(i,j,k) + met_v(i-1,j,k) )*D0_5
1265 : end do
1266 : end do
1267 : end do
1268 :
1269 0 : call t_stopf('MET__transfer_windsToWalls')
1270 :
1271 0 : end subroutine transfer_windsToWalls
1272 :
1273 0 : subroutine get_model_time()
1274 : implicit none
1275 : integer yr, mon, day, ncsec ! components of a date
1276 :
1277 0 : call t_startf('MET__get_model_time')
1278 :
1279 0 : call get_curr_date(yr, mon, day, ncsec)
1280 :
1281 0 : curr_mod_time = get_time_float( yr, mon, day, ncsec )
1282 0 : next_mod_time = curr_mod_time + get_step_size()/seconds_per_day
1283 :
1284 0 : call t_stopf('MET__get_model_time')
1285 :
1286 0 : end subroutine get_model_time
1287 :
1288 : !------------------------------------------------------------------------------
1289 : !------------------------------------------------------------------------------
1290 0 : subroutine check_files()
1291 :
1292 : use shr_sys_mod, only: shr_sys_system
1293 : use ioFileMod, only: getfil
1294 :
1295 : implicit none
1296 :
1297 : !-----------------------------------------------------------------------
1298 : ! ... local variables
1299 : !-----------------------------------------------------------------------
1300 : character(len=256) :: ctmp
1301 : character(len=256) :: loc_fname
1302 : integer :: istat
1303 :
1304 :
1305 0 : if (next_mod_time > curr_data_times(size(curr_data_times))) then
1306 0 : if ( .not. associated(next_data_times) ) then
1307 : ! open next file...
1308 0 : next_filename = incr_filename( curr_filename )
1309 0 : call open_met_datafile( next_filename, next_fileid, next_data_times, met_data_path )
1310 : endif
1311 : endif
1312 :
1313 0 : if ( associated(next_data_times) ) then
1314 0 : if (curr_mod_time >= next_data_times(1)) then
1315 :
1316 : ! close current file ...
1317 0 : call pio_closefile( curr_fileid )
1318 0 : if (masterproc) then
1319 : ! remove if requested
1320 0 : if( met_remove_file ) then
1321 0 : call getfil( curr_filename, loc_fname, 0 )
1322 0 : write(iulog,*) 'check_files: removing file = ',trim(loc_fname)
1323 0 : ctmp = 'rm -f ' // trim(loc_fname)
1324 0 : write(iulog,*) 'check_files: fsystem issuing command - '
1325 0 : write(iulog,*) trim(ctmp)
1326 0 : call shr_sys_system( ctmp, istat )
1327 : end if
1328 : endif
1329 :
1330 0 : curr_filename = next_filename
1331 0 : curr_fileid = next_fileid
1332 :
1333 0 : deallocate( curr_data_times )
1334 0 : allocate( curr_data_times( size( next_data_times ) ) )
1335 0 : curr_data_times(:) = next_data_times(:)
1336 :
1337 0 : next_filename = ''
1338 :
1339 0 : deallocate( next_data_times )
1340 : nullify( next_data_times )
1341 :
1342 : endif
1343 : endif
1344 :
1345 0 : end subroutine check_files
1346 :
1347 : !-----------------------------------------------------------------------
1348 : !-----------------------------------------------------------------------
1349 0 : function incr_filename( filename )
1350 :
1351 : !-----------------------------------------------------------------------
1352 : ! ... Increment or decrement a date string withing a filename
1353 : ! the filename date section is assumed to be of the form
1354 : ! yyyy-dd-mm
1355 : !-----------------------------------------------------------------------
1356 :
1357 : use string_utils, only : incstr
1358 : use shr_file_mod, only : shr_file_getunit, shr_file_freeunit
1359 :
1360 : implicit none
1361 :
1362 :
1363 : character(len=*), intent(in) :: filename ! present dynamical dataset filename
1364 : character(len=256) :: incr_filename ! next filename in the sequence
1365 : character(len=*), parameter :: subname = 'incr_filename'
1366 :
1367 : ! set new next_filename ...
1368 :
1369 : !-----------------------------------------------------------------------
1370 : ! ... local variables
1371 : !-----------------------------------------------------------------------
1372 : integer :: pos, pos1, istat
1373 : character(len=256) :: fn_new, line
1374 : character(len=6) :: seconds
1375 : character(len=5) :: num
1376 : integer :: ios,unitnumber
1377 :
1378 0 : if ( len_trim(met_filenames_list) == 0) then
1379 : !-----------------------------------------------------------------------
1380 : ! ... ccm type filename
1381 : !-----------------------------------------------------------------------
1382 0 : pos = len_trim( filename )
1383 0 : fn_new = filename(:pos)
1384 0 : if (masterproc) write(iulog,*) 'incr_flnm: old filename = ',trim(fn_new)
1385 0 : if( fn_new(pos-2:) == '.nc' ) then
1386 0 : pos = pos - 3
1387 : end if
1388 0 : istat = incstr( fn_new(:pos), 1 )
1389 0 : if( istat /= 0 ) then
1390 0 : write(iulog,*) 'incr_flnm: incstr returned ', istat
1391 0 : write(iulog,*) ' while trying to decrement ',trim( fn_new )
1392 0 : call endrun (subname // ':: ERRROR - check atm.log for error message')
1393 : end if
1394 :
1395 : else
1396 :
1397 : ! open met_filenames_list
1398 0 : if (masterproc) write(iulog,*) 'incr_flnm: old filename = ',trim(filename)
1399 0 : if (masterproc) write(iulog,*) 'incr_flnm: open met_filenames_list : ',met_filenames_list
1400 0 : unitnumber = shr_file_getUnit()
1401 0 : open( unit=unitnumber, file=met_filenames_list, iostat=ios, status="OLD")
1402 0 : if (ios /= 0) then
1403 0 : call endrun('not able to open met_filenames_list file: '//met_filenames_list)
1404 : endif
1405 :
1406 : ! read file names
1407 0 : read( unit=unitnumber, fmt='(A)', iostat=ios ) line
1408 0 : if (ios /= 0) then
1409 0 : call endrun('not able to increment file name from met_filenames_list file: '//met_filenames_list)
1410 : endif
1411 0 : do while( trim(line) /= trim(filename) )
1412 0 : read( unit=unitnumber, fmt='(A)', iostat=ios ) line
1413 0 : if (ios /= 0) then
1414 0 : call endrun('not able to increment file name from met_filenames_list file: '//met_filenames_list)
1415 : endif
1416 : enddo
1417 :
1418 0 : read( unit=unitnumber, fmt='(A)', iostat=ios ) line
1419 0 : if (ios /= 0) then
1420 0 : call endrun('not able to increment file name from met_filenames_list file: '//met_filenames_list)
1421 : endif
1422 0 : fn_new = trim(line)
1423 :
1424 0 : close(unit=unitnumber)
1425 0 : call shr_file_freeUnit(unitnumber)
1426 : endif
1427 0 : incr_filename = trim(fn_new)
1428 0 : if (masterproc) write(iulog,*) 'incr_flnm: new filename = ',incr_filename
1429 :
1430 0 : end function incr_filename
1431 :
1432 : !------------------------------------------------------------------------------
1433 : !------------------------------------------------------------------------------
1434 0 : subroutine find_times( itms, fids, datatm, datatp, time )
1435 :
1436 : implicit none
1437 :
1438 : integer, intent(out) :: itms(2) ! record numbers that bracket time
1439 : type(file_desc_t), intent(out) :: fids(2) ! ids of files that contains these recs
1440 : real(r8), intent(in) :: time ! time of interest
1441 : real(r8), intent(out):: datatm, datatp
1442 : character(len=*), parameter :: subname = 'find_times'
1443 :
1444 : integer np1 ! current forward time index of dataset
1445 : integer n,i !
1446 : integer :: curr_tsize, next_tsize, all_tsize
1447 :
1448 0 : real(r8), allocatable, dimension(:):: all_data_times
1449 :
1450 0 : curr_tsize = size(curr_data_times)
1451 0 : next_tsize = 0
1452 0 : if ( associated(next_data_times)) next_tsize = size(next_data_times)
1453 :
1454 0 : all_tsize = curr_tsize + next_tsize
1455 :
1456 0 : allocate( all_data_times( all_tsize ) )
1457 :
1458 0 : all_data_times(:curr_tsize) = curr_data_times(:)
1459 0 : if (next_tsize > 0) all_data_times(curr_tsize+1:all_tsize) = next_data_times(:)
1460 :
1461 : ! find bracketing times
1462 0 : do n=1, all_tsize-1
1463 0 : np1 = n + 1
1464 0 : datatm = all_data_times(n)
1465 0 : datatp = all_data_times(np1)
1466 0 : if ( (time .ge. datatm) .and. (time .le. datatp) ) then
1467 : goto 20
1468 : endif
1469 : enddo
1470 :
1471 0 : write(iulog,*)'FIND_TIMES: Failed to find dates bracketing desired time =', time
1472 0 : write(iulog,*)' datatm = ',datatm
1473 0 : write(iulog,*)' datatp = ',datatp
1474 0 : write(iulog,*)' all_data_times = ',all_data_times
1475 :
1476 0 : call endrun (subname // ':: ERRROR - check atm.log for error message')
1477 :
1478 : 20 continue
1479 :
1480 0 : deallocate( all_data_times )
1481 :
1482 0 : itms(1) = n
1483 0 : itms(2) = np1
1484 0 : fids(:) = curr_fileid
1485 :
1486 0 : do i=1,2
1487 0 : if ( itms(i) > curr_tsize ) then
1488 0 : itms(i) = itms(i) - curr_tsize
1489 0 : fids(i) = next_fileid
1490 : endif
1491 : enddo
1492 :
1493 0 : end subroutine find_times
1494 :
1495 : !------------------------------------------------------------------------------
1496 : !------------------------------------------------------------------------------
1497 0 : subroutine read_next_ps(grid)
1498 : use ncdio_atm, only: infld
1499 :
1500 : implicit none
1501 :
1502 : type (T_FVDYCORE_GRID), intent(in) :: grid
1503 :
1504 : integer :: recnos(2)
1505 0 : type(file_desc_t) :: fids(2)
1506 : character(len=8) :: varname
1507 : integer :: ifirstxy, ilastxy, jfirstxy, jlastxy
1508 :
1509 0 : real(r8) :: wrk_xy(grid%ifirstxy:grid%ilastxy, grid%jfirstxy:grid%jlastxy )
1510 :
1511 : logical :: readvar
1512 :
1513 0 : if(online_test) then
1514 0 : varname='arch_PS'
1515 : else
1516 0 : varname='PS'
1517 : end if
1518 :
1519 0 : jfirstxy= grid%jfirstxy
1520 0 : jlastxy = grid%jlastxy
1521 0 : ifirstxy= grid%ifirstxy
1522 0 : ilastxy = grid%ilastxy
1523 :
1524 0 : call find_times( recnos, fids, datatimemn, datatimepn, next_mod_time )
1525 :
1526 : call infld(varname, fids(1), 'lon', 'lat', ifirstxy, ilastxy, jfirstxy, jlastxy, &
1527 0 : wrk_xy, readvar, gridname='fv_centers', timelevel=recnos(1))
1528 :
1529 : ! transpose xy -> yz decomposition
1530 0 : call transpose_xy2yz_2d( wrk_xy, met_psi_next(nm)%data, grid )
1531 :
1532 : call infld(varname, fids(2), 'lon', 'lat', ifirstxy, ilastxy, jfirstxy, jlastxy, &
1533 0 : wrk_xy, readvar, gridname='fv_centers', timelevel=recnos(2))
1534 :
1535 : ! transpose xy -> yz decomposition
1536 0 : call transpose_xy2yz_2d( wrk_xy, met_psi_next(np)%data, grid )
1537 :
1538 0 : if(masterproc) write(iulog,*)'READ_NEXT_PS: Read meteorological data '
1539 :
1540 0 : end subroutine read_next_ps
1541 :
1542 : !------------------------------------------------------------------------
1543 : !------------------------------------------------------------------------
1544 0 : subroutine transpose_xy2yz_2d( xy_2d, yz_2d, grid )
1545 :
1546 : implicit none
1547 : type (T_FVDYCORE_GRID), intent(in) :: grid
1548 : real(r8), intent(in) :: xy_2d(grid%ifirstxy:grid%ilastxy, grid%jfirstxy:grid%jlastxy)
1549 : real(r8), intent(out) :: yz_2d(1:grid%im, grid%jfirst:grid%jlast)
1550 :
1551 0 : real(r8) :: xy3(grid%ifirstxy:grid%ilastxy, grid%jfirstxy:grid%jlastxy, 1:grid%npr_z )
1552 : integer :: i,j,k
1553 :
1554 0 : if (grid%iam .lt. grid%npes_xy) then
1555 0 : if ( grid%twod_decomp == 1 ) then
1556 :
1557 : #if defined( SPMD )
1558 : !$omp parallel do private(i,j,k)
1559 0 : do k=1,grid%npr_z
1560 0 : do j=grid%jfirstxy,grid%jlastxy
1561 0 : do i=grid%ifirstxy,grid%ilastxy
1562 0 : xy3(i,j,k) = xy_2d(i,j)
1563 : enddo
1564 : enddo
1565 : enddo
1566 :
1567 : call mp_sendirr(grid%commxy, grid%xy2d_to_yz2d%SendDesc, &
1568 : grid%xy2d_to_yz2d%RecvDesc, xy3, yz_2d, &
1569 0 : modc=grid%modc_dynrun )
1570 : call mp_recvirr(grid%commxy, grid%xy2d_to_yz2d%SendDesc, &
1571 : grid%xy2d_to_yz2d%RecvDesc, xy3, yz_2d, &
1572 0 : modc=grid%modc_dynrun )
1573 : #endif
1574 :
1575 : else
1576 0 : yz_2d(:,:) = xy_2d(:,:)
1577 : endif
1578 : endif ! (grid%iam .lt. grid%npes_xy)
1579 :
1580 0 : end subroutine transpose_xy2yz_2d
1581 :
1582 : !------------------------------------------------------------------------
1583 : !------------------------------------------------------------------------
1584 0 : subroutine transpose_xy2yz_3d( xy_3d, yz_3d, grid )
1585 :
1586 : implicit none
1587 : type (T_FVDYCORE_GRID), intent(in) :: grid
1588 : real(r8), intent(in) :: xy_3d(grid%ifirstxy:grid%ilastxy, grid%jfirstxy:grid%jlastxy, 1:grid%km )
1589 : real(r8), intent(out) :: yz_3d(1:grid%im, grid%jfirst:grid%jlast, grid%kfirst:grid%klast)
1590 :
1591 0 : if (grid%iam .lt. grid%npes_xy) then
1592 0 : if ( grid%twod_decomp == 1 ) then
1593 : #if defined( SPMD )
1594 : call mp_sendirr( grid%commxy, grid%ijk_xy_to_yz%SendDesc, &
1595 : grid%ijk_xy_to_yz%RecvDesc, xy_3d, yz_3d, &
1596 0 : modc=grid%modc_dynrun )
1597 : call mp_recvirr( grid%commxy, grid%ijk_xy_to_yz%SendDesc, &
1598 : grid%ijk_xy_to_yz%RecvDesc, xy_3d, yz_3d, &
1599 0 : modc=grid%modc_dynrun )
1600 : #endif
1601 : else
1602 0 : yz_3d(:,:,:) = xy_3d(:,:,:)
1603 : endif
1604 : endif ! (grid%iam .lt. grid%npes_xy)
1605 :
1606 0 : end subroutine transpose_xy2yz_3d
1607 :
1608 :
1609 :
1610 : !------------------------------------------------------------------------
1611 : !------------------------------------------------------------------------
1612 0 : subroutine read_next_metdata(grid)
1613 : use ncdio_atm, only: infld
1614 : use cam_grid_support, only: cam_grid_check, cam_grid_id
1615 : use cam_grid_support, only: cam_grid_get_dim_names
1616 :
1617 : implicit none
1618 :
1619 : type (T_FVDYCORE_GRID), intent(in) :: grid
1620 : integer recnos(2), i !
1621 0 : type(file_desc_t) :: fids(2)
1622 :
1623 : character(len=8) :: Uname, Vname, Tname, Qname, psname
1624 : character(len=8) :: dim1name, dim2name
1625 : integer :: im, jm, km
1626 : logical :: readvar
1627 : integer :: ifirstxy, ilastxy, jfirstxy, jlastxy
1628 0 : real(r8), allocatable :: wrk2_xy(:,:)
1629 0 : real(r8), allocatable :: wrk3_xy(:,:,:)
1630 0 : real(r8), allocatable :: tmp_data(:,:,:)
1631 : integer :: elev1,blev1, elev2,blev2
1632 : integer :: elev3,blev3, elev4,blev4
1633 : integer :: grid_id ! grid ID for data mapping
1634 :
1635 0 : call t_startf('MET__read_next_metdata')
1636 :
1637 0 : jfirstxy= grid%jfirstxy
1638 0 : jlastxy = grid%jlastxy
1639 0 : ifirstxy= grid%ifirstxy
1640 0 : ilastxy = grid%ilastxy
1641 :
1642 0 : im = grid%im
1643 0 : jm = grid%jm
1644 0 : km = grid%km
1645 :
1646 0 : call find_times( recnos, fids, datatimem, datatimep, curr_mod_time )
1647 : !
1648 : ! Set up hyperslab corners
1649 : !
1650 :
1651 0 : if(online_test) then
1652 0 : Tname='arch_T'
1653 0 : Qname='arch_Q'
1654 0 : PSname='arch_PS'
1655 0 : if(met_cell_wall_winds) then
1656 0 : Uname='arch_US'
1657 0 : Vname='arch_VS'
1658 : else
1659 0 : Uname='arch_U'
1660 0 : Vname='arch_V'
1661 : end if
1662 : else
1663 0 : Tname='T'
1664 0 : Qname='Q'
1665 0 : PSname='PS'
1666 0 : if(met_cell_wall_winds) then
1667 0 : Uname='US'
1668 0 : Vname='VS'
1669 : else
1670 0 : Uname='U'
1671 0 : Vname='V'
1672 : end if
1673 :
1674 : end if
1675 :
1676 :
1677 0 : if ( num_met_levels>km ) then
1678 :
1679 0 : blev1 = 1
1680 0 : elev1 = km
1681 :
1682 0 : blev2 = num_met_levels-km+1
1683 0 : elev2 = num_met_levels
1684 :
1685 0 : blev3 = num_met_levels-km+1
1686 0 : elev3 = num_met_levels
1687 :
1688 0 : blev4 = 1
1689 0 : elev4 = num_met_levels
1690 :
1691 : else
1692 :
1693 0 : blev1 = km-num_met_levels+1
1694 0 : elev1 = km
1695 :
1696 0 : blev2 = 1
1697 0 : elev2 = num_met_levels
1698 :
1699 0 : blev3 = 1
1700 0 : elev3 = km
1701 :
1702 0 : blev4 = km-num_met_levels+1
1703 0 : elev4 = km
1704 :
1705 : endif
1706 :
1707 0 : allocate(tmp_data(pcols, 1:num_met_levels, begchunk:endchunk))
1708 0 : allocate(wrk2_xy(ifirstxy:ilastxy, jfirstxy:jlastxy))
1709 0 : allocate(wrk3_xy(ifirstxy:ilastxy, jfirstxy:jlastxy, 1:max(km,num_met_levels)))
1710 :
1711 : ! physgrid intput for FV is probably always lon/lat but let's be pedantic
1712 0 : grid_id = cam_grid_id('physgrid')
1713 0 : if (.not. cam_grid_check(grid_id)) then
1714 0 : call endrun('read_next_metdata: Internal error, no "physgrid" grid')
1715 : end if
1716 0 : call cam_grid_get_dim_names(grid_id, dim1name, dim2name)
1717 :
1718 0 : do i=1,2
1719 :
1720 0 : met_ti(i)%data = 0._r8
1721 :
1722 : call infld(Tname, fids(i), dim1name, 'lev', dim2name, 1, pcols, 1,num_met_levels , &
1723 0 : begchunk, endchunk, tmp_data, readvar, gridname='physgrid',timelevel=recnos(i))
1724 :
1725 0 : met_ti(i)%data(:,blev1:elev1,:) = tmp_data(:, blev2:elev2, :)
1726 :
1727 0 : met_qi(i)%data = 0._r8
1728 :
1729 : call infld(Qname, fids(i), dim1name, 'lev', dim2name, 1, pcols, 1,num_met_levels, &
1730 0 : begchunk, endchunk, tmp_data, readvar, gridname='physgrid',timelevel=recnos(i))
1731 :
1732 0 : met_qi(i)%data(:,blev1:elev1,:) = tmp_data(:, blev2:elev2, :)
1733 :
1734 0 : if (met_cell_wall_winds) then
1735 :
1736 0 : wrk3_xy = 0._r8
1737 0 : met_usi(i)%data(:,:,:) = 0._r8
1738 : call infld(Uname, fids(i), 'lon', 'slat', 'lev', ifirstxy, ilastxy, jfirstxy, jlastxy, &
1739 0 : 1,num_met_levels, wrk3_xy(:,:,blev4:elev4), readvar, gridname='fv_u_stagger',timelevel=recnos(i))
1740 :
1741 : ! transpose xy -> yz decomposition
1742 0 : call transpose_xy2yz_3d( wrk3_xy(:,:,blev3:elev3), met_usi(i)%data(:,:,:), grid )
1743 :
1744 0 : wrk3_xy = 0._r8
1745 0 : met_vsi(i)%data(:,:,:) = 0._r8
1746 : call infld(Vname, fids(i), 'slon', 'lat', 'lev', ifirstxy, ilastxy, jfirstxy, jlastxy, &
1747 0 : 1,num_met_levels, wrk3_xy(:,:,blev4:elev4), readvar, gridname='fv_v_stagger',timelevel=recnos(i))
1748 :
1749 : ! transpose xy -> yz decomposition
1750 0 : call transpose_xy2yz_3d( wrk3_xy(:,:,blev3:elev3), met_vsi(i)%data(:,:,:), grid )
1751 :
1752 : else
1753 :
1754 : ! read into lower portion of the array...
1755 :
1756 0 : wrk3_xy = 0._r8
1757 0 : met_ui(i)%data = 0._r8
1758 : call infld(Uname, fids(i), 'lon', 'lat', 'lev', ifirstxy, ilastxy, jfirstxy, jlastxy, &
1759 0 : 1,num_met_levels, wrk3_xy(:,:,blev4:elev4), readvar, gridname='fv_centers',timelevel=recnos(i))
1760 :
1761 : ! transpose xy -> yz decomposition
1762 0 : call transpose_xy2yz_3d( wrk3_xy(:,:,blev3:elev3), met_ui(i)%data(:,:,:), grid )
1763 :
1764 0 : wrk3_xy = 0._r8
1765 0 : met_vi(i)%data = 0._r8
1766 : call infld(Vname, fids(i), 'lon', 'lat', 'lev', ifirstxy, ilastxy, jfirstxy, jlastxy, &
1767 0 : 1,num_met_levels, wrk3_xy(:,:,blev4:elev4), readvar, gridname='fv_centers',timelevel=recnos(i))
1768 :
1769 : ! transpose xy -> yz decomposition
1770 0 : call transpose_xy2yz_3d( wrk3_xy(:,:,blev3:elev3), met_vi(i)%data(:,:,:), grid )
1771 :
1772 : endif ! met_cell_wall_winds
1773 :
1774 : call infld(PSname, fids(i), 'lon', 'lat', ifirstxy, ilastxy, jfirstxy, jlastxy, &
1775 0 : wrk2_xy, readvar, gridname='fv_centers', timelevel=recnos(i))
1776 :
1777 : ! transpose xy -> yz decomposition
1778 0 : call transpose_xy2yz_2d( wrk2_xy, met_psi_curr(i)%data, grid )
1779 :
1780 : enddo
1781 :
1782 0 : deallocate(tmp_data)
1783 0 : deallocate(wrk3_xy)
1784 0 : deallocate(wrk2_xy)
1785 :
1786 : ! 2-D feilds
1787 0 : call read_phys_srf_flds( )
1788 :
1789 0 : if(masterproc) write(iulog,*)'READ_NEXT_METDATA: Read meteorological data '
1790 :
1791 0 : call t_stopf('MET__read_next_metdata')
1792 :
1793 0 : end subroutine read_next_metdata
1794 :
1795 : !------------------------------------------------------------------------------
1796 : !------------------------------------------------------------------------------
1797 0 : subroutine read_phys_srf_flds( )
1798 0 : use ncdio_atm, only: infld
1799 :
1800 : integer :: i, recnos(2)
1801 0 : type(file_desc_t) :: fids(2)
1802 : logical :: readvar
1803 :
1804 0 : call find_times( recnos, fids, datatimem, datatimep, curr_mod_time )
1805 0 : do i=1,2
1806 :
1807 0 : call infld(met_shflx_name, fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1808 0 : met_shflxi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1809 :
1810 0 : if (has_lhflx) then
1811 : call infld('LHFLX', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1812 0 : met_lhflxi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1813 : end if
1814 :
1815 : call infld(met_qflx_name, fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1816 0 : met_qflxi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1817 : call infld('TAUX', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1818 0 : met_tauxi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1819 : call infld('TAUY', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1820 0 : met_tauyi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1821 :
1822 0 : if ( .not.met_srf_feedback ) then
1823 : call infld('SNOWH', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1824 0 : met_snowhi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1825 : endif
1826 :
1827 0 : if (has_ts) then
1828 : call infld('TS', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1829 0 : met_tsi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1830 : endif
1831 :
1832 0 : if (met_srf_rad) then
1833 : call infld('ASDIR', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1834 0 : met_asdiri(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1835 : call infld('ASDIF', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1836 0 : met_asdifi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1837 : call infld('ALDIR', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1838 0 : met_aldiri(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1839 : call infld('ALDIF', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1840 0 : met_aldifi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1841 : call infld('LWUP', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1842 0 : met_lwupi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1843 : endif
1844 :
1845 0 : if (met_srf_refs) then
1846 : call infld('QREF', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1847 0 : met_qrefi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1848 : call infld('TREF', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1849 0 : met_trefi(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1850 : call infld('U10', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1851 0 : met_u10i(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1852 : endif
1853 :
1854 0 : if (met_srf_sst) then
1855 : call infld('SST', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1856 0 : met_ssti(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1857 : call infld('ICEFRAC', fids(i), 'lon', 'lat', 1, pcols, begchunk, endchunk, &
1858 0 : met_icefraci(i)%data, readvar, gridname='physgrid',timelevel=recnos(i))
1859 : endif
1860 : enddo
1861 0 : end subroutine read_phys_srf_flds
1862 :
1863 : !------------------------------------------------------------------------------
1864 : !------------------------------------------------------------------------------
1865 0 : subroutine interp_phys_srf_flds( )
1866 0 : use phys_grid, only: get_ncols_p
1867 : real(r4) :: fact1, fact2
1868 : real(r8) :: deltat
1869 : integer :: i, c, ncol
1870 :
1871 0 : deltat = datatimep - datatimem
1872 0 : fact1 = (datatimep - curr_mod_time)/deltat
1873 0 : fact2 = D1_0-fact1
1874 :
1875 :
1876 0 : do c=begchunk,endchunk
1877 0 : ncol = get_ncols_p(c)
1878 0 : do i=1,ncol
1879 0 : if (has_lhflx) then
1880 0 : met_lhflx(i,c) = fact1*met_lhflxi(nm)%data(i,c) + fact2*met_lhflxi(np)%data(i,c)
1881 : end if
1882 0 : met_shflx(i,c) = fact1*met_shflxi(nm)%data(i,c) + fact2*met_shflxi(np)%data(i,c)
1883 0 : met_qflx(i,c) = fact1*met_qflxi(nm)%data(i,c) + fact2*met_qflxi(np)%data(i,c)
1884 0 : met_taux(i,c) = fact1*met_tauxi(nm)%data(i,c) + fact2*met_tauxi(np)%data(i,c)
1885 0 : met_tauy(i,c) = fact1*met_tauyi(nm)%data(i,c) + fact2*met_tauyi(np)%data(i,c)
1886 : enddo
1887 : enddo
1888 0 : if ( .not.met_srf_feedback ) then
1889 0 : do c=begchunk,endchunk
1890 0 : ncol = get_ncols_p(c)
1891 0 : do i=1,ncol
1892 0 : met_snowh(i,c) = fact1*met_snowhi(nm)%data(i,c) + fact2*met_snowhi(np)%data(i,c)
1893 : enddo
1894 : enddo
1895 : endif
1896 0 : if (has_ts) then
1897 0 : do c=begchunk,endchunk
1898 0 : ncol = get_ncols_p(c)
1899 0 : do i=1,ncol
1900 0 : met_ts(i,c) = fact1*met_tsi(nm)%data(i,c) + fact2*met_tsi(np)%data(i,c)
1901 : enddo
1902 : enddo
1903 : endif
1904 0 : if (met_srf_rad) then
1905 0 : do c=begchunk,endchunk
1906 0 : ncol = get_ncols_p(c)
1907 0 : do i=1,ncol
1908 : ! Albedo can be fillValue from the meteorology where the is no downwelling
1909 : ! solar. However, this changes slowly, so for interpolation use either end-point
1910 : ! if nothing is present. If there is no solar, then the albedo won't matter, so
1911 : ! should not cause problems.
1912 0 : if (met_asdiri(nm)%data(i,c) == srf_fill_value) then
1913 0 : met_asdir(i,c) = met_asdiri(np)%data(i,c)
1914 0 : else if (met_asdiri(np)%data(i,c) == srf_fill_value) then
1915 0 : met_asdir(i,c) = met_asdiri(nm)%data(i,c)
1916 : else
1917 0 : met_asdir(i,c) = fact1*met_asdiri(nm)%data(i,c) + fact2*met_asdiri(np)%data(i,c)
1918 : endif
1919 :
1920 0 : if (met_asdifi(nm)%data(i,c) == srf_fill_value) then
1921 0 : met_asdif(i,c) = met_asdifi(np)%data(i,c)
1922 0 : else if (met_asdifi(np)%data(i,c) == srf_fill_value) then
1923 0 : met_asdif(i,c) = met_asdifi(nm)%data(i,c)
1924 : else
1925 0 : met_asdif(i,c) = fact1*met_asdifi(nm)%data(i,c) + fact2*met_asdifi(np)%data(i,c)
1926 : endif
1927 :
1928 0 : if (met_aldiri(nm)%data(i,c) == srf_fill_value) then
1929 0 : met_aldir(i,c) = met_aldiri(np)%data(i,c)
1930 0 : else if (met_aldiri(np)%data(i,c) == srf_fill_value) then
1931 0 : met_aldir(i,c) = met_aldiri(nm)%data(i,c)
1932 : else
1933 0 : met_aldir(i,c) = fact1*met_aldiri(nm)%data(i,c) + fact2*met_aldiri(np)%data(i,c)
1934 : endif
1935 :
1936 0 : if (met_aldifi(nm)%data(i,c) == srf_fill_value) then
1937 0 : met_aldif(i,c) = met_aldifi(np)%data(i,c)
1938 0 : else if (met_aldifi(np)%data(i,c) == srf_fill_value) then
1939 0 : met_aldif(i,c) = met_aldifi(nm)%data(i,c)
1940 : else
1941 0 : met_aldif(i,c) = fact1*met_aldifi(nm)%data(i,c) + fact2*met_aldifi(np)%data(i,c)
1942 : endif
1943 :
1944 0 : met_lwup(i,c) = fact1*met_lwupi(nm)%data(i,c) + fact2*met_lwupi(np)%data(i,c)
1945 : enddo
1946 : enddo
1947 : endif
1948 0 : if (met_srf_refs) then
1949 0 : do c=begchunk,endchunk
1950 0 : ncol = get_ncols_p(c)
1951 0 : do i=1,ncol
1952 0 : met_qref(i,c) = fact1*met_qrefi(nm)%data(i,c) + fact2*met_qrefi(np)%data(i,c)
1953 0 : met_tref(i,c) = fact1*met_trefi(nm)%data(i,c) + fact2*met_trefi(np)%data(i,c)
1954 0 : met_u10(i,c) = fact1*met_u10i(nm)%data(i,c) + fact2*met_u10i(np)%data(i,c)
1955 : enddo
1956 : enddo
1957 : endif
1958 0 : if (met_srf_sst) then
1959 0 : do c=begchunk,endchunk
1960 0 : ncol = get_ncols_p(c)
1961 0 : do i=1,ncol
1962 : ! The sst is fill value over land, which should not change from timestep to
1963 : ! timestep, but just in case use the sst value if only one is present.
1964 0 : if (met_ssti(nm)%data(i,c) == srf_fill_value) then
1965 0 : met_sst(i,c) = met_ssti(np)%data(i,c)
1966 0 : else if (met_ssti(np)%data(i,c) == srf_fill_value) then
1967 0 : met_sst(i,c) = met_ssti(nm)%data(i,c)
1968 : else
1969 0 : met_sst(i,c) = fact1*met_ssti(nm)%data(i,c) + fact2*met_ssti(np)%data(i,c)
1970 : endif
1971 :
1972 0 : if (met_icefraci(nm)%data(i,c) == srf_fill_value) then
1973 0 : met_icefrac(i,c) = met_icefraci(np)%data(i,c)
1974 0 : else if (met_ssti(np)%data(i,c) == srf_fill_value) then
1975 0 : met_icefrac(i,c) = met_icefraci(nm)%data(i,c)
1976 : else
1977 0 : met_icefrac(i,c) = fact1*met_icefraci(nm)%data(i,c) + fact2*met_icefraci(np)%data(i,c)
1978 : endif
1979 : enddo
1980 : enddo
1981 : endif
1982 :
1983 0 : end subroutine interp_phys_srf_flds
1984 : !------------------------------------------------------------------------------
1985 : !------------------------------------------------------------------------------
1986 0 : subroutine interpolate_metdata(grid)
1987 0 : use phys_grid, only: get_ncols_p
1988 : #if defined( SPMD )
1989 : use mod_comm, only : mp_send4d_ns, mp_recv4d_ns
1990 : #endif
1991 :
1992 : implicit none
1993 : type (T_FVDYCORE_GRID), intent(in) :: grid
1994 :
1995 : real(r4) fact1, fact2
1996 : real(r4) nfact1, nfact2
1997 : real(r8) deltat,deltatn
1998 : integer :: i,c,k, ncol
1999 : integer :: im, jm, km, jfirst, jlast, kfirst, klast, ng_d, ng_s
2000 :
2001 0 : im = grid%im
2002 0 : jm = grid%jm
2003 0 : km = grid%km
2004 0 : jfirst = grid%jfirst
2005 0 : jlast = grid%jlast
2006 0 : kfirst = grid%kfirst
2007 0 : klast = grid%klast
2008 0 : ng_d = grid%ng_d
2009 0 : ng_s = grid%ng_s
2010 :
2011 0 : call t_startf('MET__interpolate_metdata')
2012 :
2013 0 : deltat = datatimep - datatimem
2014 0 : deltatn = datatimepn - datatimemn
2015 :
2016 0 : fact1 = (datatimep - curr_mod_time)/deltat
2017 : ! fact2 = (curr_mod_time - datatimem)/deltat
2018 0 : fact2 = D1_0-fact1
2019 :
2020 0 : nfact1 = (datatimepn - next_mod_time)/deltatn
2021 : ! nfact2 = (next_mod_time - datatimemn)/deltatn
2022 0 : nfact2 = D1_0-nfact1
2023 :
2024 0 : met_q = 0.0_r8
2025 0 : do c=begchunk,endchunk
2026 0 : ncol = get_ncols_p(c)
2027 0 : do k=1,pver
2028 0 : do i=1,ncol
2029 0 : met_t(i,k,c) = fact1*met_ti(nm)%data(i,k,c) + fact2*met_ti(np)%data(i,k,c)
2030 0 : met_q(i,k,c) = fact1*met_qi(nm)%data(i,k,c) + fact2*met_qi(np)%data(i,k,c)
2031 : enddo
2032 : enddo
2033 : enddo
2034 :
2035 0 : if (.not. online_test) where (met_q .lt. D0_0) met_q = D0_0
2036 :
2037 0 : met_ps_next(:,:) = nfact1*met_psi_next(nm)%data(:,:) + nfact2*met_psi_next(np)%data(:,:)
2038 0 : met_ps_curr(:,:) = fact1*met_psi_curr(nm)%data(:,:) + fact2*met_psi_curr(np)%data(:,:)
2039 :
2040 0 : call interp_phys_srf_flds()
2041 :
2042 0 : if (has_ts) then
2043 0 : do c=begchunk,endchunk
2044 0 : ncol = get_ncols_p(c)
2045 0 : do i=1,ncol
2046 0 : met_ts(i,c) = fact1*met_tsi(nm)%data(i,c) + fact2*met_tsi(np)%data(i,c)
2047 : enddo
2048 : enddo
2049 : endif
2050 :
2051 0 : if ( .not. met_cell_wall_winds ) then
2052 :
2053 0 : met_u(1:im,jfirst:jlast,kfirst:klast) = fact1*met_ui(nm)%data(1:im,jfirst:jlast,kfirst:klast) &
2054 0 : + fact2*met_ui(np)%data(1:im,jfirst:jlast,kfirst:klast)
2055 0 : met_v(1:im,jfirst:jlast,kfirst:klast) = fact1*met_vi(nm)%data(1:im,jfirst:jlast,kfirst:klast) &
2056 0 : + fact2*met_vi(np)%data(1:im,jfirst:jlast,kfirst:klast)
2057 :
2058 : ! ghost u,v
2059 : #if defined( SPMD )
2060 : call mp_send4d_ns( grid%commxy, im, jm, km, 1, jfirst, jlast, &
2061 0 : kfirst, klast, ng_d, ng_d, met_u )
2062 : call mp_send4d_ns( grid%commxy, im, jm, km, 1, jfirst, jlast, &
2063 0 : kfirst, klast, ng_d, ng_s, met_v )
2064 : call mp_recv4d_ns( grid%commxy, im, jm, km, 1, jfirst, jlast, &
2065 0 : kfirst, klast, ng_d, ng_d, met_u )
2066 : call mp_recv4d_ns( grid%commxy, im, jm, km, 1, jfirst, jlast, &
2067 0 : kfirst, klast, ng_d, ng_s, met_v )
2068 : #endif
2069 :
2070 : ! average to cell walls (vorticity winds)
2071 0 : call transfer_windsToWalls(grid)
2072 : else
2073 0 : met_us(:,jfirst:jlast,kfirst:klast) = fact1*met_usi(nm)%data(:,jfirst:jlast,kfirst:klast) + &
2074 0 : fact2*met_usi(np)%data(:,jfirst:jlast,kfirst:klast)
2075 0 : met_vs(:,jfirst:jlast,kfirst:klast) = fact1*met_vsi(nm)%data(:,jfirst:jlast,kfirst:klast) + &
2076 0 : fact2*met_vsi(np)%data(:,jfirst:jlast,kfirst:klast)
2077 :
2078 : endif
2079 :
2080 : ! ghost staggered u,v
2081 : ! WS 2006.04.11: not necessary here since it will be done in cd_core
2082 :
2083 : ! write(iulog,*)'INTERPOLATE_METDATA: complete.'
2084 :
2085 0 : call t_stopf('MET__interpolate_metdata')
2086 :
2087 0 : end subroutine interpolate_metdata
2088 :
2089 : !-----------------------------------------------------------------------
2090 : !-----------------------------------------------------------------------
2091 0 : subroutine get_dimension( fid, dname, dsize )
2092 : implicit none
2093 : type(file_desc_t), intent(in) :: fid
2094 : character(*), intent(in) :: dname
2095 : integer, intent(out) :: dsize
2096 :
2097 : integer :: dimid, ierr
2098 :
2099 0 : ierr = pio_inq_dimid( fid, dname, dimid )
2100 0 : ierr = pio_inq_dimlen( fid, dimid, dsize )
2101 :
2102 0 : end subroutine get_dimension
2103 :
2104 : !-----------------------------------------------------------------------
2105 : !-----------------------------------------------------------------------
2106 0 : subroutine open_met_datafile( fname, fileid, times, datapath, check_dims, grid )
2107 :
2108 : use ioFileMod, only: getfil
2109 : use pio, only: pio_seterrorhandling, PIO_INTERNAL_ERROR, PIO_BCAST_ERROR, PIO_NOERR
2110 :
2111 : implicit none
2112 :
2113 : character(*), intent(in) :: fname
2114 : type(file_desc_t), intent(inout) :: fileid
2115 : real(r8), pointer, intent(inout) :: times(:)
2116 : character(*), intent(in) :: datapath
2117 : logical, optional, intent(in) :: check_dims
2118 : type (T_FVDYCORE_GRID), optional, intent(in) :: grid
2119 :
2120 : character(len=256) :: filepath
2121 : character(len=256) :: filen
2122 : character(len=*), parameter :: subname = 'open_met_datafile'
2123 : integer :: year, month, day, dsize, i, timesize
2124 : integer :: dateid,secid
2125 0 : integer, allocatable , dimension(:) :: dates, datesecs
2126 : integer :: ierr
2127 : integer :: im, jm, km
2128 : integer :: varid
2129 :
2130 : !
2131 : ! open file and get fileid
2132 : !
2133 0 : if (len_trim( datapath )>0) then
2134 0 : filepath = trim(datapath)//'/'//trim(fname)
2135 : else
2136 0 : filepath = trim(fname)
2137 : endif
2138 0 : call getfil( filepath, filen, 0 )
2139 :
2140 0 : call cam_pio_openfile( fileid, filen, 0 )
2141 :
2142 0 : call pio_seterrorhandling(fileid, PIO_BCAST_ERROR)
2143 :
2144 0 : ierr = pio_inq_varid( fileid, 'TS', varid )
2145 0 : has_ts = ierr==PIO_NOERR
2146 :
2147 0 : ierr = pio_inq_varid( fileid, 'LHFLX', varid )
2148 0 : has_lhflx = ierr==PIO_NOERR
2149 :
2150 0 : call pio_seterrorhandling(fileid, PIO_INTERNAL_ERROR)
2151 :
2152 0 : if (masterproc) write(iulog,*) 'open_met_datafile: ',trim(filen)
2153 :
2154 0 : call get_dimension( fileid, 'time', timesize )
2155 :
2156 0 : if ( associated(times) ) deallocate(times)
2157 0 : allocate( times(timesize) )
2158 :
2159 0 : allocate( dates(timesize) )
2160 0 : allocate( datesecs(timesize) )
2161 :
2162 0 : ierr = pio_inq_varid( fileid, 'date', dateid )
2163 0 : ierr = pio_inq_varid( fileid, 'datesec', secid )
2164 :
2165 0 : ierr = pio_get_var( fileid, dateid, dates )
2166 0 : ierr = pio_get_var( fileid, secid, datesecs )
2167 :
2168 0 : do i=1,timesize
2169 0 : year = dates(i) / 10000
2170 0 : month = mod(dates(i),10000)/100
2171 0 : day = mod(dates(i),100)
2172 0 : times(i) = get_time_float( year, month, day, datesecs(i) )
2173 : enddo
2174 :
2175 0 : deallocate( dates )
2176 0 : deallocate( datesecs )
2177 :
2178 :
2179 : !
2180 : ! check that the data dim sizes match models dimensions
2181 : !
2182 0 : if (present(check_dims) .and. present(grid)) then
2183 0 : im = grid%im
2184 0 : jm = grid%jm
2185 0 : km = grid%km
2186 :
2187 0 : if (check_dims) then
2188 :
2189 0 : call get_dimension( fileid, 'lon', dsize )
2190 0 : if (dsize /= im) then
2191 0 : write(iulog,*)'open_met_datafile: lonsiz=',dsize,' must = ',im
2192 0 : call endrun (subname // ':: ERRROR - check atm.log for error message')
2193 : endif
2194 0 : call get_dimension( fileid, 'lat', dsize )
2195 0 : if (dsize /= jm) then
2196 0 : write(iulog,*)'open_met_datafile: latsiz=',dsize,' must = ',jm
2197 0 : call endrun (subname // ':: ERRROR - check atm.log for error message')
2198 : endif
2199 0 : call get_dimension( fileid, 'lev', dsize )
2200 0 : met_levels = min( dsize, km )
2201 0 : num_met_levels = dsize
2202 : endif
2203 : endif
2204 :
2205 0 : if (met_srf_rad) then
2206 0 : ierr = pio_inq_varid( fileid, 'ASDIR', varid )
2207 0 : ierr = pio_get_att( fileid, varid, '_FillValue', srf_fill_value)
2208 : endif
2209 :
2210 0 : end subroutine open_met_datafile
2211 :
2212 : !------------------------------------------------------------------------------
2213 : !------------------------------------------------------------------------------
2214 0 : function get_time_float( year, month, day, sec )
2215 :
2216 : ! returns float representation of time -- number of days
2217 : ! since 1 jan 0001 00:00:00.000
2218 :
2219 : implicit none
2220 :
2221 : integer, intent(in) :: year, month, day
2222 : integer, intent(in) :: sec
2223 : real(r8) :: get_time_float
2224 :
2225 : ! ref date is 1 jan 0001
2226 :
2227 : integer :: refyr, refmn, refdy
2228 : real(r8) :: refsc, fltdy
2229 : integer :: doy(12)
2230 :
2231 : ! jan feb mar apr may jun jul aug sep oct nov dec
2232 : ! 31 28 31 30 31 30 31 31 31 31 30 31
2233 : data doy / 1, 32, 60, 91,121,152,182,213,244,274,305,335 /
2234 :
2235 0 : refyr = 1
2236 0 : refmn = 1
2237 0 : refdy = 1
2238 0 : refsc = D0_0
2239 :
2240 0 : if ( timemgr_is_caltype(trim(shr_cal_gregorian))) then
2241 0 : fltdy = greg2jday(year, month, day) - greg2jday(refyr,refmn,refdy)
2242 : else ! assume no_leap (all years are 365 days)
2243 : fltdy = (year - refyr)*days_per_non_leapyear + &
2244 0 : (doy(month)-doy(refmn)) + &
2245 0 : (day-refdy)
2246 : endif
2247 :
2248 0 : get_time_float = fltdy + ((sec-refsc)/seconds_per_day)
2249 :
2250 0 : endfunction get_time_float
2251 :
2252 : !-----------------------------------------------------------------------
2253 : ! ... Return Julian day number given Gregorian date.
2254 : !
2255 : ! Algorithm from Hatcher,D.A., Simple Formulae for Julian Day Numbers
2256 : ! and Calendar Dates, Q.Jl.R.astr.Soc. (1984) v25, pp 53-55.
2257 : !-----------------------------------------------------------------------
2258 0 : function greg2jday( year, month, day )
2259 :
2260 : implicit none
2261 :
2262 : integer, intent(in) :: year, month, day
2263 : integer :: greg2jday
2264 :
2265 : !-----------------------------------------------------------------------
2266 : ! ... Local variables
2267 : !-----------------------------------------------------------------------
2268 : integer :: ap, mp
2269 : integer :: y, d, n, g
2270 :
2271 : !-----------------------------------------------------------------------
2272 : ! ... Modify year and month numbers
2273 : !-----------------------------------------------------------------------
2274 0 : ap = year - (12 - month)/10
2275 0 : mp = MOD( month-3,12 )
2276 0 : if( mp < 0 ) then
2277 0 : mp = mp + 12
2278 : end if
2279 :
2280 : !-----------------------------------------------------------------------
2281 : ! ... Julian day
2282 : !-----------------------------------------------------------------------
2283 0 : y = INT( days_per_year*( ap + 4712 ) )
2284 0 : d = INT( days_per_month*mp + D0_5 )
2285 0 : n = y + d + day + 59
2286 0 : g = INT( D0_75*INT( ap/100 + 49 ) ) - 38
2287 0 : greg2jday = n - g
2288 :
2289 0 : end function greg2jday
2290 :
2291 : !-----------------------------------------------------------------------
2292 : !-----------------------------------------------------------------------
2293 0 : subroutine set_met_rlx( )
2294 :
2295 : use pmgrid, only: plev
2296 : use hycoef, only: hypm, ps0
2297 :
2298 : integer :: k, k_cnt, k_top
2299 : real(r8), parameter :: h0 = 7._r8 ! scale height (km)
2300 : real(r8), parameter :: hsec = 3600._r8 ! seconds per hour
2301 : real(r8) :: p_top, p_bot
2302 : real(r8) :: p_bot_top, p_bot_bot
2303 : real(r8) :: met_max_rlxdt, dtime_hrs
2304 :
2305 : 996 format( 'set_met_rlx: ',a15, I10 )
2306 : 997 format( 'set_met_rlx: ',a15, E10.2 )
2307 : 998 format( 'set_met_rlx: ',a15, PLEV(E10.2))
2308 : 999 format( 'set_met_rlx: ',a15, PLEV(F10.5))
2309 : 993 format( 'set_met_rlx: ',a25, E10.2 )
2310 : 991 format( 'set_met_rlx: ',a25, I4, E10.2, E10.2 )
2311 :
2312 0 : met_rlx(:) = 999._r8
2313 :
2314 0 : dtime_hrs = get_step_size()/hsec ! hours
2315 :
2316 0 : if (met_rlx_time > dtime_hrs) then
2317 0 : met_max_rlxdt = dtime_hrs/met_rlx_time
2318 0 : elseif (met_rlx_time < 0._r8) then
2319 0 : met_max_rlxdt = 0._r8
2320 : else
2321 0 : met_max_rlxdt = 1._r8
2322 : endif
2323 :
2324 0 : if (masterproc) then
2325 0 : write(iulog,fmt=993) ' met_rlx_time in hrs= ', met_rlx_time
2326 0 : write(iulog,fmt=993) ' met_max_rlxdt in % = ', met_max_rlxdt*100._r8
2327 : endif
2328 :
2329 0 : p_top = ps0 * exp( - met_rlx_top/h0 )
2330 0 : p_bot = ps0 * exp( - met_rlx_bot/h0 )
2331 :
2332 0 : p_bot_top = ps0 * exp( - met_rlx_bot_top/h0 )
2333 0 : p_bot_bot = ps0 * exp( - met_rlx_bot_bot/h0 )
2334 :
2335 0 : if (masterproc) then
2336 0 : write(iulog,fmt=997) 'p_top = ',p_top
2337 0 : write(iulog,fmt=997) 'p_bot = ',p_bot
2338 0 : write(iulog,fmt=997) 'p_bot_top = ',p_bot_top
2339 0 : write(iulog,fmt=997) 'p_bot_bot = ',p_bot_bot
2340 : endif
2341 :
2342 0 : if ( p_bot < hypm( pver-met_levels+1 ) .and. ( met_levels < pver ) ) then
2343 0 : call endrun( 'set_met_rlx: met_rlx_bot is too high ' )
2344 : endif
2345 :
2346 0 : met_rlx = 0._r8
2347 :
2348 0 : where (p_top <= hypm .and. hypm < p_bot)
2349 : met_rlx = 999._r8
2350 : endwhere
2351 :
2352 0 : where( p_bot <= hypm .and. hypm <= p_bot_top )
2353 : met_rlx = met_max_rlxdt
2354 : endwhere
2355 :
2356 0 : where(p_bot_top < hypm .and. hypm <= p_bot_bot)
2357 : met_rlx = -999._r8
2358 : endwhere
2359 :
2360 0 : if ( any( met_rlx(:) /= met_max_rlxdt) ) then
2361 0 : k_top = max(plev - met_levels, 1)
2362 :
2363 : ! ramp region at model top
2364 0 : k_cnt = count(p_top < hypm .and. hypm < p_bot)
2365 0 : if (k_cnt > 0) then
2366 : k_top = max(plev - met_levels, 1)
2367 0 : do while ( met_rlx(k_top) /= 999._r8 )
2368 0 : k_top = k_top + 1
2369 0 : if ( k_top == pver ) then
2370 0 : call endrun ( 'set_met_rlx: cannot find ramped region ')
2371 : endif
2372 : enddo
2373 :
2374 0 : met_rlx(1:k_top) = 0._r8
2375 :
2376 0 : k_cnt = count( met_rlx == 999._r8 )
2377 :
2378 0 : if (masterproc) then
2379 0 : write(iulog,*) 'top of model ramped region:'
2380 0 : write(iulog,fmt=996) 'k_cnt = ',k_cnt
2381 0 : write(iulog,fmt=996) 'k_top = ',k_top
2382 : endif
2383 :
2384 0 : do k = k_top,k_top+k_cnt
2385 0 : met_rlx(k) = met_max_rlxdt*real( k - k_top ) / real(k_cnt)
2386 : enddo
2387 : endif
2388 :
2389 : ! ramp region at model bottom
2390 0 : k_cnt = count(p_bot_top < hypm .and. hypm < p_bot_bot)
2391 0 : if (k_cnt > 0) then
2392 0 : k_top = max(plev - met_levels, 1)
2393 0 : do while ( met_rlx(k_top) /= -999._r8 )
2394 0 : k_top = k_top + 1
2395 0 : if ( k_top == pver ) then
2396 0 : call endrun ( 'set_met_rlx: cannot find ramped region ')
2397 : endif
2398 : enddo
2399 :
2400 0 : if (masterproc) then
2401 0 : write(iulog,*) 'bottom of model ramped region:'
2402 0 : write(iulog,fmt=996) 'k_cnt = ',k_cnt
2403 0 : write(iulog,fmt=996) 'k_top = ',k_top
2404 : endif
2405 :
2406 0 : do k = k_top,k_top+k_cnt-1
2407 0 : met_rlx(k) = met_max_rlxdt*(1._r8 - real( k - k_top +1) / real(k_cnt))
2408 : enddo
2409 : endif
2410 :
2411 : endif
2412 :
2413 0 : if (masterproc) then
2414 0 : write(iulog,fmt=996) ' met_levels = ',met_levels
2415 0 : write(iulog,fmt=996) 'non-zero terms:',count( met_rlx /= 0._r8 )
2416 : endif
2417 :
2418 0 : if ( met_levels < count( met_rlx /= 0._r8 ) ) then
2419 0 : call endrun('set_met_rlx: met_rlx_top is too high for the meteorology data')
2420 : endif
2421 :
2422 0 : if (masterproc) then
2423 0 : do k = max(plev - met_levels, 1), pver
2424 0 : write(iulog,fmt=991) 'press levels, met_rlx = ', k, hypm(k), met_rlx(k)
2425 : end do
2426 : endif
2427 :
2428 0 : if ( any( (met_rlx > 1._r8) .or. (met_rlx < 0._r8) ) ) then
2429 0 : if (masterproc) then
2430 0 : write(iulog,fmt=993) 'set_met_rlx: dtime_hrs in hrs = ', dtime_hrs
2431 0 : write(iulog,fmt=993) 'set_met_rlx:met_rlx_time in hrs = ', met_rlx_time
2432 0 : write(iulog,fmt=993) 'set_met_rlx: met_max_rlxdt = ', met_max_rlxdt
2433 0 : write(iulog,fmt=999) 'set_met_rlx: met_rlx = ', met_rlx
2434 : endif
2435 0 : call endrun('Offline meteorology relaxation function not set correctly.')
2436 : endif
2437 :
2438 0 : end subroutine set_met_rlx
2439 :
2440 0 : end module metdata
|