Line data Source code
1 : module nudging
2 : !=====================================================================
3 : !
4 : ! Purpose: Implement Nudging of the model state of U,V,T,Q, and/or PS
5 : ! toward specified values from analyses.
6 : !
7 : ! Author: Patrick Callaghan
8 : !
9 : ! Description:
10 : !
11 : ! This module assumes that the user has {U,V,T,Q,PS} values from analyses
12 : ! which have been preprocessed onto the current model grid and adjusted
13 : ! for differences in topography. It is also assumed that these resulting
14 : ! values and are stored in individual files which are indexed with respect
15 : ! to year, month, day, and second of the day. When the model is inbetween
16 : ! the given begining and ending times, a relaxation forcing is added to
17 : ! nudge the model toward the analyses values determined from the forcing
18 : ! option specified. After the model passes the ending analyses time, the
19 : ! forcing discontinues.
20 : !
21 : ! Some analyses products can have gaps in the available data, where values
22 : ! are missing for some interval of time. When files are missing, the nudging
23 : ! force is switched off for that interval of time, so we effectively 'coast'
24 : ! thru the gap.
25 : !
26 : ! Currently, the nudging module is set up to accomodate nudging of PS
27 : ! values, however that functionality requires forcing that is applied in
28 : ! the selected dycore and is not yet implemented.
29 : !
30 : ! The nudging of the model toward the analyses data is controlled by
31 : ! the 'nudging_nl' namelist in 'user_nl_cam'; whose variables control the
32 : ! time interval over which nudging is applied, the strength of the nudging
33 : ! tendencies, and its spatial distribution.
34 : !
35 : ! FORCING:
36 : ! --------
37 : ! Nudging tendencies are applied as a relaxation force between the current
38 : ! model state values and target state values derived from the avalilable
39 : ! analyses. The form of the target values is selected by the 'Nudge_Force_Opt'
40 : ! option, the timescale of the forcing is determined from the given
41 : ! 'Nudge_TimeScale_Opt', and the nudging strength Alpha=[0.,1.] for each
42 : ! variable is specified by the 'Nudge_Xcoef' values. Where X={U,V,T,Q,PS}
43 : !
44 : ! F_nudge = Alpha*((Target-Model(t_curr))/TimeScale
45 : !
46 : !
47 : ! WINDOWING:
48 : ! ----------
49 : ! The region of applied nudging can be limited using Horizontal/Vertical
50 : ! window functions that are constructed using a parameterization of the
51 : ! Heaviside step function.
52 : !
53 : ! The Heaviside window function is the product of separate horizonal and vertical
54 : ! windows that are controled via 12 parameters:
55 : !
56 : ! Nudge_Hwin_lat0: Specify the horizontal center of the window in degrees.
57 : ! Nudge_Hwin_lon0: The longitude must be in the range [0,360] and the
58 : ! latitude should be [-90,+90].
59 : ! Nudge_Hwin_latWidth: Specify the lat and lon widths of the window as positive
60 : ! Nudge_Hwin_lonWidth: values in degrees.Setting a width to a large value (e.g. 999)
61 : ! renders the window a constant in that direction.
62 : ! Nudge_Hwin_latDelta: Controls the sharpness of the window transition with a
63 : ! Nudge_Hwin_lonDelta: length in degrees. Small non-zero values yeild a step
64 : ! function while a large value yeilds a smoother transition.
65 : ! Nudge_Hwin_Invert : A logical flag used to invert the horizontal window function
66 : ! to get its compliment.(e.g. to nudge outside a given window).
67 : !
68 : ! Nudge_Vwin_Lindex: In the vertical, the window is specified in terms of model
69 : ! Nudge_Vwin_Ldelta: level indcies. The High and Low transition levels should
70 : ! Nudge_Vwin_Hindex: range from [0,(NLEV+1)]. The transition lengths are also
71 : ! Nudge_Vwin_Hdelta: specified in terms of model indices. For a window function
72 : ! constant in the vertical, the Low index should be set to 0,
73 : ! the High index should be set to (NLEV+1), and the transition
74 : ! lengths should be set to 0.001
75 : ! Nudge_Vwin_Invert : A logical flag used to invert the vertical window function
76 : ! to get its compliment.
77 : !
78 : ! EXAMPLE: For a channel window function centered at the equator and independent
79 : ! of the vertical (30 levels):
80 : ! Nudge_Hwin_lat0 = 0. Nudge_Vwin_Lindex = 0.
81 : ! Nudge_Hwin_latWidth = 30. Nudge_Vwin_Ldelta = 0.001
82 : ! Nudge_Hwin_latDelta = 5.0 Nudge_Vwin_Hindex = 31.
83 : ! Nudge_Hwin_lon0 = 180. Nudge_Vwin_Hdelta = 0.001
84 : ! Nudge_Hwin_lonWidth = 999. Nudge_Vwin_Invert = .false.
85 : ! Nudge_Hwin_lonDelta = 1.0
86 : ! Nudge_Hwin_Invert = .false.
87 : !
88 : ! If on the other hand one wanted to apply nudging at the poles and
89 : ! not at the equator, the settings would be similar but with:
90 : ! Nudge_Hwin_Invert = .true.
91 : !
92 : ! A user can preview the window resulting from a given set of namelist values before
93 : ! running the model. Lookat_NudgeWindow.ncl is a script avalable in the tools directory
94 : ! which will read in the values for a given namelist and display the resulting window.
95 : !
96 : ! The module is currently configured for only 1 window function. It can readily be
97 : ! extended for multiple windows if the need arises.
98 : !
99 : !
100 : ! Input/Output Values:
101 : ! Forcing contributions are available for history file output by
102 : ! the names: {'Nudge_U','Nudge_V','Nudge_T',and 'Nudge_Q'}
103 : ! The target values that the model state is nudged toward are available for history
104 : ! file output via the variables: {'Target_U','Target_V','Target_T',and 'Target_Q'}
105 : !
106 : ! &nudging_nl
107 : ! Nudge_Model - LOGICAL toggle to activate nudging.
108 : ! TRUE -> Nudging is on.
109 : ! FALSE -> Nudging is off. [DEFAULT]
110 : !
111 : ! Nudge_Path - CHAR path to the analyses files.
112 : ! (e.g. '/glade/scratch/USER/inputdata/nudging/ERAI-Data/')
113 : !
114 : ! Nudge_File_Template - CHAR Analyses filename with year, month, day, and second
115 : ! values replaced by %y, %m, %d, and %s respectively.
116 : ! (e.g. '%y/ERAI_ne30np4_L30.cam2.i.%y-%m-%d-%s.nc')
117 : !
118 : ! Nudge_Times_Per_Day - INT Number of analyses files available per day.
119 : ! 1 --> daily analyses.
120 : ! 4 --> 6 hourly analyses.
121 : ! 8 --> 3 hourly.
122 : !
123 : ! Model_Times_Per_Day - INT Number of times to update the model state (used for nudging)
124 : ! each day. The value is restricted to be longer than the
125 : ! current model timestep and shorter than the analyses
126 : ! timestep. As this number is increased, the nudging
127 : ! force has the form of newtonian cooling.
128 : ! 48 --> 1800 Second timestep.
129 : ! 96 --> 900 Second timestep.
130 : !
131 : ! Nudge_Beg_Year - INT nudging begining year. [1979- ]
132 : ! Nudge_Beg_Month - INT nudging begining month. [1-12]
133 : ! Nudge_Beg_Day - INT nudging begining day. [1-31]
134 : ! Nudge_End_Year - INT nudging ending year. [1979-]
135 : ! Nudge_End_Month - INT nudging ending month. [1-12]
136 : ! Nudge_End_Day - INT nudging ending day. [1-31]
137 : !
138 : ! Nudge_Force_Opt - INT Index to select the nudging Target for a relaxation
139 : ! forcing of the form:
140 : ! where (t'==Analysis times ; t==Model Times)
141 : !
142 : ! 0 -> NEXT-OBS: Target=Anal(t'_next) [DEFAULT]
143 : ! 1 -> LINEAR: Target=(F*Anal(t'_curr) +(1-F)*Anal(t'_next))
144 : ! F =(t'_next - t_curr )/Tdlt_Anal
145 : !
146 : ! Nudge_TimeScale_Opt - INT Index to select the timescale for nudging.
147 : ! where (t'==Analysis times ; t==Model Times)
148 : !
149 : ! 0 --> TimeScale = 1/Tdlt_Anal [DEFAULT]
150 : ! 1 --> TimeScale = 1/(t'_next - t_curr )
151 : !
152 : ! Nudge_Uprof - INT index of profile structure to use for U. [0,1,2]
153 : ! Nudge_Vprof - INT index of profile structure to use for V. [0,1,2]
154 : ! Nudge_Tprof - INT index of profile structure to use for T. [0,1,2]
155 : ! Nudge_Qprof - INT index of profile structure to use for Q. [0,1,2]
156 : ! Nudge_PSprof - INT index of profile structure to use for PS. [0,N/A]
157 : !
158 : ! The spatial distribution is specified with a profile index.
159 : ! Where: 0 == OFF (No Nudging of this variable)
160 : ! 1 == CONSTANT (Spatially Uniform Nudging)
161 : ! 2 == HEAVISIDE WINDOW FUNCTION
162 : !
163 : ! Nudge_Ucoef - REAL fractional nudging coeffcient for U.
164 : ! Nudge_Vcoef - REAL fractional nudging coeffcient for V.
165 : ! Nudge_Tcoef - REAL fractional nudging coeffcient for T.
166 : ! Nudge_Qcoef - REAL fractional nudging coeffcient for Q.
167 : ! Nudge_PScoef - REAL fractional nudging coeffcient for PS.
168 : !
169 : ! The strength of the nudging is specified as a fractional
170 : ! coeffcient between [0,1].
171 : !
172 : ! Nudge_Hwin_lat0 - REAL latitudinal center of window in degrees.
173 : ! Nudge_Hwin_lon0 - REAL longitudinal center of window in degrees.
174 : ! Nudge_Hwin_latWidth - REAL latitudinal width of window in degrees.
175 : ! Nudge_Hwin_lonWidth - REAL longitudinal width of window in degrees.
176 : ! Nudge_Hwin_latDelta - REAL latitudinal transition length of window in degrees.
177 : ! Nudge_Hwin_lonDelta - REAL longitudinal transition length of window in degrees.
178 : ! Nudge_Hwin_Invert - LOGICAL FALSE= value=1 inside the specified window, 0 outside
179 : ! TRUE = value=0 inside the specified window, 1 outside
180 : ! Nudge_Vwin_Lindex - REAL LO model index of transition
181 : ! Nudge_Vwin_Hindex - REAL HI model index of transition
182 : ! Nudge_Vwin_Ldelta - REAL LO transition length
183 : ! Nudge_Vwin_Hdelta - REAL HI transition length
184 : ! Nudge_Vwin_Invert - LOGICAL FALSE= value=1 inside the specified window, 0 outside
185 : ! TRUE = value=0 inside the specified window, 1 outside
186 : ! /
187 : !
188 : !================
189 : !
190 : ! TO DO:
191 : ! -----------
192 : ! ** Implement Ps Nudging????
193 : !
194 : !=====================================================================
195 : ! Useful modules
196 : !------------------
197 : use shr_kind_mod, only: r8=>SHR_KIND_R8, cs=>SHR_KIND_CS, cl=>SHR_KIND_CL
198 : use time_manager, only: timemgr_time_ge, timemgr_time_inc, get_curr_date
199 : use time_manager, only: get_step_size
200 : use cam_abortutils, only: endrun
201 : use spmd_utils, only: masterproc, mstrid=>masterprocid, mpicom, mpi_success
202 : use spmd_utils, only: mpi_integer, mpi_real8, mpi_logical, mpi_character
203 : use cam_logfile, only: iulog
204 : use zonal_mean_mod, only: ZonalMean_t
205 :
206 : ! Set all Global values and routines to private by default
207 : ! and then explicitly set their exposure.
208 : !----------------------------------------------------------
209 : implicit none
210 : private
211 :
212 : public :: Nudge_Model,Nudge_ON
213 : public :: nudging_readnl
214 : public :: nudging_init
215 : public :: nudging_timestep_init
216 : public :: nudging_timestep_tend
217 : private :: nudging_update_analyses
218 : private :: nudging_set_PSprofile
219 : private :: nudging_set_profile
220 : private :: calc_DryStaticEnergy
221 : public :: nudging_final
222 :
223 : ! Nudging Parameters
224 : !--------------------
225 : logical :: Nudge_Model =.false.
226 : logical :: Nudge_ON =.false.
227 : logical :: Nudge_Initialized =.false.
228 : character(len=cl) :: Nudge_Path
229 : character(len=cs) :: Nudge_File,Nudge_File_Template
230 : integer :: Nudge_Force_Opt
231 : integer :: Nudge_TimeScale_Opt
232 : integer :: Nudge_TSmode
233 : integer :: Nudge_Times_Per_Day
234 : integer :: Model_Times_Per_Day
235 : real(r8) :: Nudge_Ucoef,Nudge_Vcoef
236 : integer :: Nudge_Uprof,Nudge_Vprof
237 : real(r8) :: Nudge_Qcoef,Nudge_Tcoef
238 : integer :: Nudge_Qprof,Nudge_Tprof
239 : real(r8) :: Nudge_PScoef
240 : integer :: Nudge_PSprof
241 : integer :: Nudge_Beg_Year ,Nudge_Beg_Month
242 : integer :: Nudge_Beg_Day ,Nudge_Beg_Sec
243 : integer :: Nudge_End_Year ,Nudge_End_Month
244 : integer :: Nudge_End_Day ,Nudge_End_Sec
245 : integer :: Nudge_Curr_Year,Nudge_Curr_Month
246 : integer :: Nudge_Curr_Day ,Nudge_Curr_Sec
247 : integer :: Nudge_Next_Year,Nudge_Next_Month
248 : integer :: Nudge_Next_Day ,Nudge_Next_Sec
249 : integer :: Nudge_Step
250 : integer :: Model_Curr_Year,Model_Curr_Month
251 : integer :: Model_Curr_Day ,Model_Curr_Sec
252 : integer :: Model_Next_Year,Model_Next_Month
253 : integer :: Model_Next_Day ,Model_Next_Sec
254 : integer :: Model_Step
255 : real(r8) :: Nudge_Hwin_lat0
256 : real(r8) :: Nudge_Hwin_latWidth
257 : real(r8) :: Nudge_Hwin_latDelta
258 : real(r8) :: Nudge_Hwin_lon0
259 : real(r8) :: Nudge_Hwin_lonWidth
260 : real(r8) :: Nudge_Hwin_lonDelta
261 : logical :: Nudge_Hwin_Invert = .false.
262 : real(r8) :: Nudge_Hwin_lo
263 : real(r8) :: Nudge_Hwin_hi
264 : real(r8) :: Nudge_Vwin_Hindex
265 : real(r8) :: Nudge_Vwin_Hdelta
266 : real(r8) :: Nudge_Vwin_Lindex
267 : real(r8) :: Nudge_Vwin_Ldelta
268 : logical :: Nudge_Vwin_Invert =.false.
269 : real(r8) :: Nudge_Vwin_lo
270 : real(r8) :: Nudge_Vwin_hi
271 : real(r8) :: Nudge_Hwin_latWidthH
272 : real(r8) :: Nudge_Hwin_lonWidthH
273 : real(r8) :: Nudge_Hwin_max
274 : real(r8) :: Nudge_Hwin_min
275 :
276 : ! Nudging Zonal Filter variables
277 : !---------------------------------
278 : logical :: Nudge_ZonalFilter =.false.
279 : integer :: Nudge_ZonalNbasis = -1
280 : type(ZonalMean_t) :: ZM
281 : real(r8),allocatable:: Zonal_Bamp2d(:)
282 : real(r8),allocatable:: Zonal_Bamp3d(:,:)
283 :
284 : ! Nudging State Arrays
285 : !-----------------------
286 : integer :: Nudge_nlon,Nudge_nlat,Nudge_ncol,Nudge_nlev
287 : real(r8),allocatable:: Target_U (:,:,:) !(pcols,pver,begchunk:endchunk)
288 : real(r8),allocatable:: Target_V (:,:,:) !(pcols,pver,begchunk:endchunk)
289 : real(r8),allocatable:: Target_T (:,:,:) !(pcols,pver,begchunk:endchunk)
290 : real(r8),allocatable:: Target_S (:,:,:) !(pcols,pver,begchunk:endchunk)
291 : real(r8),allocatable:: Target_Q (:,:,:) !(pcols,pver,begchunk:endchunk)
292 : real(r8),allocatable:: Target_PS (:,:) !(pcols,begchunk:endchunk)
293 : real(r8),allocatable:: Model_U (:,:,:) !(pcols,pver,begchunk:endchunk)
294 : real(r8),allocatable:: Model_V (:,:,:) !(pcols,pver,begchunk:endchunk)
295 : real(r8),allocatable:: Model_T (:,:,:) !(pcols,pver,begchunk:endchunk)
296 : real(r8),allocatable:: Model_S (:,:,:) !(pcols,pver,begchunk:endchunk)
297 : real(r8),allocatable:: Model_Q (:,:,:) !(pcols,pver,begchunk:endchunk)
298 : real(r8),allocatable:: Model_PS (:,:) !(pcols,begchunk:endchunk)
299 : real(r8),allocatable:: Nudge_Utau (:,:,:) !(pcols,pver,begchunk:endchunk)
300 : real(r8),allocatable:: Nudge_Vtau (:,:,:) !(pcols,pver,begchunk:endchunk)
301 : real(r8),allocatable:: Nudge_Stau (:,:,:) !(pcols,pver,begchunk:endchunk)
302 : real(r8),allocatable:: Nudge_Qtau (:,:,:) !(pcols,pver,begchunk:endchunk)
303 : real(r8),allocatable:: Nudge_PStau (:,:) !(pcols,begchunk:endchunk)
304 : real(r8),allocatable:: Nudge_Ustep (:,:,:) !(pcols,pver,begchunk:endchunk)
305 : real(r8),allocatable:: Nudge_Vstep (:,:,:) !(pcols,pver,begchunk:endchunk)
306 : real(r8),allocatable:: Nudge_Sstep (:,:,:) !(pcols,pver,begchunk:endchunk)
307 : real(r8),allocatable:: Nudge_Qstep (:,:,:) !(pcols,pver,begchunk:endchunk)
308 : real(r8),allocatable:: Nudge_PSstep(:,:) !(pcols,begchunk:endchunk)
309 :
310 : ! Nudging Observation Arrays
311 : !-----------------------------
312 : integer :: Nudge_NumObs
313 : integer,allocatable:: Nudge_ObsInd(:)
314 : logical ,allocatable::Nudge_File_Present(:)
315 : real(r8),allocatable::Nobs_U (:,:,:,:) !(pcols,pver,begchunk:endchunk,Nudge_NumObs)
316 : real(r8),allocatable::Nobs_V (:,:,:,:) !(pcols,pver,begchunk:endchunk,Nudge_NumObs)
317 : real(r8),allocatable::Nobs_T (:,:,:,:) !(pcols,pver,begchunk:endchunk,Nudge_NumObs)
318 : real(r8),allocatable::Nobs_Q (:,:,:,:) !(pcols,pver,begchunk:endchunk,Nudge_NumObs)
319 : real(r8),allocatable::Nobs_PS(:,:,:) !(pcols,begchunk:endchunk,Nudge_NumObs)
320 :
321 : contains
322 : !================================================================
323 1536 : subroutine nudging_readnl(nlfile)
324 : !
325 : ! NUDGING_READNL: Initialize default values controlling the Nudging
326 : ! process. Then read namelist values to override
327 : ! them.
328 : !===============================================================
329 : use ppgrid, only: pver
330 : use namelist_utils, only:find_group_name
331 : !
332 : ! Arguments
333 : !-------------
334 : character(len=*), intent(in) :: nlfile
335 : !
336 : ! Local Values
337 : !---------------
338 : integer :: ierr, unitn
339 :
340 : character(len=*), parameter :: prefix = 'nudging_readnl: '
341 :
342 : namelist /nudging_nl/ Nudge_Model, Nudge_Path, &
343 : Nudge_File_Template, Nudge_Force_Opt, &
344 : Nudge_TimeScale_Opt, &
345 : Nudge_Times_Per_Day, Model_Times_Per_Day, &
346 : Nudge_Ucoef , Nudge_Uprof, &
347 : Nudge_Vcoef , Nudge_Vprof, &
348 : Nudge_Qcoef , Nudge_Qprof, &
349 : Nudge_Tcoef , Nudge_Tprof, &
350 : Nudge_PScoef, Nudge_PSprof, &
351 : Nudge_Beg_Year, Nudge_Beg_Month, Nudge_Beg_Day, &
352 : Nudge_End_Year, Nudge_End_Month, Nudge_End_Day, &
353 : Nudge_Hwin_lat0, Nudge_Hwin_lon0, &
354 : Nudge_Hwin_latWidth, Nudge_Hwin_lonWidth, &
355 : Nudge_Hwin_latDelta, Nudge_Hwin_lonDelta, &
356 : Nudge_Hwin_Invert, &
357 : Nudge_Vwin_Lindex, Nudge_Vwin_Hindex, &
358 : Nudge_Vwin_Ldelta, Nudge_Vwin_Hdelta, &
359 : Nudge_Vwin_Invert
360 :
361 : ! For Zonal Mean Filtering
362 : namelist /nudging_nl/ Nudge_ZonalFilter, Nudge_ZonalNbasis
363 :
364 : ! Nudging is NOT initialized yet, For now
365 : ! Nudging will always begin/end at midnight.
366 : !--------------------------------------------
367 1536 : Nudge_Initialized =.false.
368 1536 : Nudge_ON =.false.
369 1536 : Nudge_Beg_Sec=0
370 1536 : Nudge_End_Sec=0
371 :
372 : ! Set Default Namelist values
373 : !-----------------------------
374 1536 : Nudge_Model = .false.
375 1536 : Nudge_Path = './Data/YOTC_ne30np4_001/'
376 1536 : Nudge_File_Template = 'YOTC_ne30np4_L30.cam2.i.%y-%m-%d-%s.nc'
377 1536 : Nudge_Force_Opt = 0
378 1536 : Nudge_TimeScale_Opt = 0
379 1536 : Nudge_TSmode = 0
380 1536 : Nudge_Times_Per_Day = 4
381 1536 : Model_Times_Per_Day = 4
382 1536 : Nudge_Ucoef = 0._r8
383 1536 : Nudge_Vcoef = 0._r8
384 1536 : Nudge_Qcoef = 0._r8
385 1536 : Nudge_Tcoef = 0._r8
386 1536 : Nudge_PScoef = 0._r8
387 1536 : Nudge_Uprof = 0
388 1536 : Nudge_Vprof = 0
389 1536 : Nudge_Qprof = 0
390 1536 : Nudge_Tprof = 0
391 1536 : Nudge_PSprof = 0
392 1536 : Nudge_Beg_Year = 2008
393 1536 : Nudge_Beg_Month = 5
394 1536 : Nudge_Beg_Day = 1
395 1536 : Nudge_End_Year = 2008
396 1536 : Nudge_End_Month = 9
397 1536 : Nudge_End_Day = 1
398 1536 : Nudge_Hwin_lat0 = 0._r8
399 1536 : Nudge_Hwin_latWidth = 9999._r8
400 1536 : Nudge_Hwin_latDelta = 1.0_r8
401 1536 : Nudge_Hwin_lon0 = 180._r8
402 1536 : Nudge_Hwin_lonWidth = 9999._r8
403 1536 : Nudge_Hwin_lonDelta = 1.0_r8
404 1536 : Nudge_Hwin_Invert = .false.
405 1536 : Nudge_Hwin_lo = 0.0_r8
406 1536 : Nudge_Hwin_hi = 1.0_r8
407 1536 : Nudge_Vwin_Hindex = float(pver+1)
408 1536 : Nudge_Vwin_Hdelta = 0.001_r8
409 1536 : Nudge_Vwin_Lindex = 0.0_r8
410 1536 : Nudge_Vwin_Ldelta = 0.001_r8
411 1536 : Nudge_Vwin_Invert = .false.
412 1536 : Nudge_Vwin_lo = 0.0_r8
413 1536 : Nudge_Vwin_hi = 1.0_r8
414 :
415 : ! Read in namelist values
416 : !------------------------
417 1536 : if(masterproc) then
418 2 : open(newunit=unitn, file=trim(nlfile), status='old')
419 2 : call find_group_name(unitn, 'nudging_nl', status=ierr)
420 2 : if(ierr == 0) then
421 0 : read(unitn,nudging_nl,iostat=ierr)
422 0 : if(ierr /= 0) then
423 0 : call endrun('nudging_readnl:: ERROR reading namelist')
424 : end if
425 : end if
426 2 : close(unitn)
427 : end if
428 :
429 : ! Set hi/lo values according to the given '_Invert' parameters
430 : !--------------------------------------------------------------
431 1536 : if(Nudge_Hwin_Invert) then
432 0 : Nudge_Hwin_lo = 1.0_r8
433 0 : Nudge_Hwin_hi = 0.0_r8
434 : else
435 1536 : Nudge_Hwin_lo = 0.0_r8
436 1536 : Nudge_Hwin_hi = 1.0_r8
437 : end if
438 :
439 1536 : if(Nudge_Vwin_Invert) then
440 0 : Nudge_Vwin_lo = 1.0_r8
441 0 : Nudge_Vwin_hi = 0.0_r8
442 : else
443 1536 : Nudge_Vwin_lo = 0.0_r8
444 1536 : Nudge_Vwin_hi = 1.0_r8
445 : end if
446 :
447 : ! Check for valid namelist values
448 : !----------------------------------
449 1536 : if((Nudge_Hwin_lat0 < -90._r8) .or. (Nudge_Hwin_lat0 > +90._r8)) then
450 0 : write(iulog,*) 'NUDGING: Window lat0 must be in [-90,+90]'
451 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_lat0=',Nudge_Hwin_lat0
452 0 : call endrun('nudging_readnl:: ERROR in namelist')
453 : endif
454 :
455 1536 : if((Nudge_Hwin_lon0 < 0._r8) .or. (Nudge_Hwin_lon0 >= 360._r8)) then
456 0 : write(iulog,*) 'NUDGING: Window lon0 must be in [0,+360)'
457 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_lon0=',Nudge_Hwin_lon0
458 0 : call endrun('nudging_readnl:: ERROR in namelist')
459 : endif
460 :
461 : if((Nudge_Vwin_Lindex > Nudge_Vwin_Hindex) .or. &
462 : (Nudge_Vwin_Hindex > float(pver+1)) .or. (Nudge_Vwin_Hindex < 0._r8) .or. &
463 1536 : (Nudge_Vwin_Lindex > float(pver+1)) .or. (Nudge_Vwin_Lindex < 0._r8) ) then
464 0 : write(iulog,*) 'NUDGING: Window Lindex must be in [0,pver+1]'
465 0 : write(iulog,*) 'NUDGING: Window Hindex must be in [0,pver+1]'
466 0 : write(iulog,*) 'NUDGING: Lindex must be LE than Hindex'
467 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_Lindex=',Nudge_Vwin_Lindex
468 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_Hindex=',Nudge_Vwin_Hindex
469 0 : call endrun('nudging_readnl:: ERROR in namelist')
470 : endif
471 :
472 : if((Nudge_Hwin_latDelta <= 0._r8) .or. (Nudge_Hwin_lonDelta <= 0._r8) .or. &
473 1536 : (Nudge_Vwin_Hdelta <= 0._r8) .or. (Nudge_Vwin_Ldelta <= 0._r8) ) then
474 0 : write(iulog,*) 'NUDGING: Window Deltas must be positive'
475 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_latDelta=',Nudge_Hwin_latDelta
476 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_lonDelta=',Nudge_Hwin_lonDelta
477 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_Hdelta=',Nudge_Vwin_Hdelta
478 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_Ldelta=',Nudge_Vwin_Ldelta
479 0 : call endrun('nudging_readnl:: ERROR in namelist')
480 :
481 : endif
482 :
483 1536 : if((Nudge_Hwin_latWidth <= 0._r8) .or. (Nudge_Hwin_lonWidth <= 0._r8)) then
484 0 : write(iulog,*) 'NUDGING: Window widths must be positive'
485 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_latWidth=',Nudge_Hwin_latWidth
486 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_lonWidth=',Nudge_Hwin_lonWidth
487 0 : call endrun('nudging_readnl:: ERROR in namelist')
488 : endif
489 :
490 : ! Broadcast namelist variables
491 : !------------------------------
492 : call MPI_bcast(Nudge_Path , len(Nudge_Path), &
493 1536 : mpi_character, mstrid, mpicom, ierr)
494 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Path ')
495 : call MPI_bcast(Nudge_File_Template,len(Nudge_File_Template), &
496 1536 : mpi_character, mstrid, mpicom, ierr)
497 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_File_Template')
498 1536 : call MPI_bcast(Nudge_Model , 1, mpi_logical, mstrid, mpicom, ierr)
499 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Model')
500 1536 : call MPI_bcast(Nudge_Initialized , 1, mpi_logical, mstrid, mpicom, ierr)
501 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Initialized')
502 1536 : call MPI_bcast(Nudge_ON , 1, mpi_logical, mstrid, mpicom, ierr)
503 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_ON')
504 1536 : call MPI_bcast(Nudge_Force_Opt , 1, mpi_integer, mstrid, mpicom, ierr)
505 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Force_Opt')
506 1536 : call MPI_bcast(Nudge_TimeScale_Opt, 1, mpi_integer, mstrid, mpicom, ierr)
507 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_TimeScale_Opt')
508 1536 : call MPI_bcast(Nudge_TSmode , 1, mpi_integer, mstrid, mpicom, ierr)
509 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_TSmode')
510 1536 : call MPI_bcast(Nudge_Times_Per_Day, 1, mpi_integer, mstrid, mpicom, ierr)
511 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Times_Per_Day')
512 1536 : call MPI_bcast(Model_Times_Per_Day, 1, mpi_integer, mstrid, mpicom, ierr)
513 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Model_Times_Per_Day')
514 1536 : call MPI_bcast(Nudge_Ucoef , 1, mpi_real8 , mstrid, mpicom, ierr)
515 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Ucoef')
516 1536 : call MPI_bcast(Nudge_Vcoef , 1, mpi_real8 , mstrid, mpicom, ierr)
517 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Vcoef')
518 1536 : call MPI_bcast(Nudge_Tcoef , 1, mpi_real8 , mstrid, mpicom, ierr)
519 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Tcoef')
520 1536 : call MPI_bcast(Nudge_Qcoef , 1, mpi_real8 , mstrid, mpicom, ierr)
521 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Qcoef')
522 1536 : call MPI_bcast(Nudge_PScoef , 1, mpi_real8 , mstrid, mpicom, ierr)
523 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_PScoef')
524 1536 : call MPI_bcast(Nudge_Uprof , 1, mpi_integer, mstrid, mpicom, ierr)
525 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Uprof')
526 1536 : call MPI_bcast(Nudge_Vprof , 1, mpi_integer, mstrid, mpicom, ierr)
527 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Vprof')
528 1536 : call MPI_bcast(Nudge_Tprof , 1, mpi_integer, mstrid, mpicom, ierr)
529 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Tprof')
530 1536 : call MPI_bcast(Nudge_Qprof , 1, mpi_integer, mstrid, mpicom, ierr)
531 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Qprof')
532 1536 : call MPI_bcast(Nudge_PSprof , 1, mpi_integer, mstrid, mpicom, ierr)
533 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_PSprof')
534 1536 : call MPI_bcast(Nudge_Beg_Year , 1, mpi_integer, mstrid, mpicom, ierr)
535 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Beg_Year')
536 1536 : call MPI_bcast(Nudge_Beg_Month , 1, mpi_integer, mstrid, mpicom, ierr)
537 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Beg_Month')
538 1536 : call MPI_bcast(Nudge_Beg_Day , 1, mpi_integer, mstrid, mpicom, ierr)
539 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Beg_Day')
540 1536 : call MPI_bcast(Nudge_Beg_Sec , 1, mpi_integer, mstrid, mpicom, ierr)
541 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Beg_Sec')
542 1536 : call MPI_bcast(Nudge_End_Year , 1, mpi_integer, mstrid, mpicom, ierr)
543 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_End_Year')
544 1536 : call MPI_bcast(Nudge_End_Month , 1, mpi_integer, mstrid, mpicom, ierr)
545 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_End_Month')
546 1536 : call MPI_bcast(Nudge_End_Day , 1, mpi_integer, mstrid, mpicom, ierr)
547 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_End_Day')
548 1536 : call MPI_bcast(Nudge_End_Sec , 1, mpi_integer, mstrid, mpicom, ierr)
549 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_End_Sec')
550 1536 : call MPI_bcast(Nudge_Hwin_lo , 1, mpi_real8 , mstrid, mpicom, ierr)
551 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_lo')
552 1536 : call MPI_bcast(Nudge_Hwin_hi , 1, mpi_real8 , mstrid, mpicom, ierr)
553 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_hi')
554 1536 : call MPI_bcast(Nudge_Hwin_lat0 , 1, mpi_real8 , mstrid, mpicom, ierr)
555 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_lat0')
556 1536 : call MPI_bcast(Nudge_Hwin_latWidth, 1, mpi_real8 , mstrid, mpicom, ierr)
557 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_latWidth')
558 1536 : call MPI_bcast(Nudge_Hwin_latDelta, 1, mpi_real8 , mstrid, mpicom, ierr)
559 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_latDelta')
560 1536 : call MPI_bcast(Nudge_Hwin_lon0 , 1, mpi_real8 , mstrid, mpicom, ierr)
561 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_lon0')
562 1536 : call MPI_bcast(Nudge_Hwin_lonWidth, 1, mpi_real8 , mstrid, mpicom, ierr)
563 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_lonWidth')
564 1536 : call MPI_bcast(Nudge_Hwin_lonDelta, 1, mpi_real8 , mstrid, mpicom, ierr)
565 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_lonDelta')
566 1536 : call MPI_bcast(Nudge_Hwin_Invert, 1, mpi_logical, mstrid, mpicom, ierr)
567 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_Invert')
568 1536 : call MPI_bcast(Nudge_Vwin_lo , 1, mpi_real8 , mstrid, mpicom, ierr)
569 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Vwin_lo')
570 1536 : call MPI_bcast(Nudge_Vwin_hi , 1, mpi_real8 , mstrid, mpicom, ierr)
571 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Vwin_hi')
572 1536 : call MPI_bcast(Nudge_Vwin_Hindex , 1, mpi_real8 , mstrid, mpicom, ierr)
573 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Vwin_Hindex')
574 1536 : call MPI_bcast(Nudge_Vwin_Hdelta , 1, mpi_real8 , mstrid, mpicom, ierr)
575 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Vwin_Hdelta')
576 1536 : call MPI_bcast(Nudge_Vwin_Lindex , 1, mpi_real8 , mstrid, mpicom, ierr)
577 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Vwin_Lindex')
578 1536 : call MPI_bcast(Nudge_Vwin_Ldelta , 1, mpi_real8 , mstrid, mpicom, ierr)
579 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Vwin_Ldelta')
580 1536 : call MPI_bcast(Nudge_Vwin_Invert, 1, mpi_logical, mstrid, mpicom, ierr)
581 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Vwin_Invert')
582 1536 : call MPI_bcast(Nudge_ZonalFilter, 1, mpi_logical, mstrid, mpicom, ierr)
583 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_ZonalFilter')
584 1536 : call MPI_bcast(Nudge_ZonalNbasis, 1, mpi_integer, mstrid, mpicom, ierr)
585 1536 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_ZonalNbasis')
586 :
587 : ! End Routine
588 : !------------
589 :
590 1536 : end subroutine nudging_readnl
591 : !================================================================
592 :
593 :
594 : !================================================================
595 0 : subroutine nudging_init
596 : !
597 : ! NUDGING_INIT: Allocate space and initialize Nudging values
598 : !===============================================================
599 : use ppgrid ,only: pver,pcols,begchunk,endchunk
600 : use error_messages,only: alloc_err
601 : use dycore ,only: dycore_is
602 : use dyn_grid ,only: get_horiz_grid_dim_d
603 : use phys_grid ,only: get_rlat_p,get_rlon_p,get_ncols_p
604 : use cam_history ,only: addfld
605 : use shr_const_mod ,only: SHR_CONST_PI
606 : use filenames ,only: interpret_filename_spec
607 :
608 : ! Local values
609 : !----------------
610 : integer :: Year,Month,Day,Sec
611 : integer :: YMD1,YMD
612 : logical :: After_Beg,Before_End
613 : integer :: istat,lchnk,ncol,icol,ilev
614 : integer :: hdim1_d,hdim2_d
615 : integer :: ierr
616 : integer :: dtime
617 : real(r8) :: rlat,rlon
618 : real(r8) :: Wprof(pver)
619 : real(r8) :: lonp,lon0,lonn,latp,lat0,latn
620 : real(r8) :: Val1_p,Val2_p,Val3_p,Val4_p
621 : real(r8) :: Val1_0,Val2_0,Val3_0,Val4_0
622 : real(r8) :: Val1_n,Val2_n,Val3_n,Val4_n
623 : integer :: nn
624 :
625 : character(len=*), parameter :: prefix = 'nudging_init: '
626 :
627 : ! Get the time step size
628 : !------------------------
629 0 : dtime = get_step_size()
630 :
631 : ! Allocate Space for Nudging data arrays
632 : !-----------------------------------------
633 0 : allocate(Target_U(pcols,pver,begchunk:endchunk),stat=istat)
634 0 : call alloc_err(istat,'nudging_init','Target_U',pcols*pver*((endchunk-begchunk)+1))
635 0 : allocate(Target_V(pcols,pver,begchunk:endchunk),stat=istat)
636 0 : call alloc_err(istat,'nudging_init','Target_V',pcols*pver*((endchunk-begchunk)+1))
637 0 : allocate(Target_T(pcols,pver,begchunk:endchunk),stat=istat)
638 0 : call alloc_err(istat,'nudging_init','Target_T',pcols*pver*((endchunk-begchunk)+1))
639 0 : allocate(Target_S(pcols,pver,begchunk:endchunk),stat=istat)
640 0 : call alloc_err(istat,'nudging_init','Target_S',pcols*pver*((endchunk-begchunk)+1))
641 0 : allocate(Target_Q(pcols,pver,begchunk:endchunk),stat=istat)
642 0 : call alloc_err(istat,'nudging_init','Target_Q',pcols*pver*((endchunk-begchunk)+1))
643 0 : allocate(Target_PS(pcols,begchunk:endchunk),stat=istat)
644 0 : call alloc_err(istat,'nudging_init','Target_PS',pcols*((endchunk-begchunk)+1))
645 :
646 0 : allocate(Model_U(pcols,pver,begchunk:endchunk),stat=istat)
647 0 : call alloc_err(istat,'nudging_init','Model_U',pcols*pver*((endchunk-begchunk)+1))
648 0 : allocate(Model_V(pcols,pver,begchunk:endchunk),stat=istat)
649 0 : call alloc_err(istat,'nudging_init','Model_V',pcols*pver*((endchunk-begchunk)+1))
650 0 : allocate(Model_T(pcols,pver,begchunk:endchunk),stat=istat)
651 0 : call alloc_err(istat,'nudging_init','Model_T',pcols*pver*((endchunk-begchunk)+1))
652 0 : allocate(Model_S(pcols,pver,begchunk:endchunk),stat=istat)
653 0 : call alloc_err(istat,'nudging_init','Model_S',pcols*pver*((endchunk-begchunk)+1))
654 0 : allocate(Model_Q(pcols,pver,begchunk:endchunk),stat=istat)
655 0 : call alloc_err(istat,'nudging_init','Model_Q',pcols*pver*((endchunk-begchunk)+1))
656 0 : allocate(Model_PS(pcols,begchunk:endchunk),stat=istat)
657 0 : call alloc_err(istat,'nudging_init','Model_PS',pcols*((endchunk-begchunk)+1))
658 :
659 : ! Allocate Space for spatial dependence of
660 : ! Nudging Coefs and Nudging Forcing.
661 : !-------------------------------------------
662 0 : allocate(Nudge_Utau(pcols,pver,begchunk:endchunk),stat=istat)
663 0 : call alloc_err(istat,'nudging_init','Nudge_Utau',pcols*pver*((endchunk-begchunk)+1))
664 0 : allocate(Nudge_Vtau(pcols,pver,begchunk:endchunk),stat=istat)
665 0 : call alloc_err(istat,'nudging_init','Nudge_Vtau',pcols*pver*((endchunk-begchunk)+1))
666 0 : allocate(Nudge_Stau(pcols,pver,begchunk:endchunk),stat=istat)
667 0 : call alloc_err(istat,'nudging_init','Nudge_Stau',pcols*pver*((endchunk-begchunk)+1))
668 0 : allocate(Nudge_Qtau(pcols,pver,begchunk:endchunk),stat=istat)
669 0 : call alloc_err(istat,'nudging_init','Nudge_Qtau',pcols*pver*((endchunk-begchunk)+1))
670 0 : allocate(Nudge_PStau(pcols,begchunk:endchunk),stat=istat)
671 0 : call alloc_err(istat,'nudging_init','Nudge_PStau',pcols*((endchunk-begchunk)+1))
672 :
673 0 : allocate(Nudge_Ustep(pcols,pver,begchunk:endchunk),stat=istat)
674 0 : call alloc_err(istat,'nudging_init','Nudge_Ustep',pcols*pver*((endchunk-begchunk)+1))
675 0 : allocate(Nudge_Vstep(pcols,pver,begchunk:endchunk),stat=istat)
676 0 : call alloc_err(istat,'nudging_init','Nudge_Vstep',pcols*pver*((endchunk-begchunk)+1))
677 0 : allocate(Nudge_Sstep(pcols,pver,begchunk:endchunk),stat=istat)
678 0 : call alloc_err(istat,'nudging_init','Nudge_Sstep',pcols*pver*((endchunk-begchunk)+1))
679 0 : allocate(Nudge_Qstep(pcols,pver,begchunk:endchunk),stat=istat)
680 0 : call alloc_err(istat,'nudging_init','Nudge_Qstep',pcols*pver*((endchunk-begchunk)+1))
681 0 : allocate(Nudge_PSstep(pcols,begchunk:endchunk),stat=istat)
682 0 : call alloc_err(istat,'nudging_init','Nudge_PSstep',pcols*((endchunk-begchunk)+1))
683 :
684 : ! Register output fields with the cam history module
685 : !-----------------------------------------------------
686 0 : call addfld( 'Nudge_U',(/ 'lev' /),'A','m/s/s' ,'U Nudging Tendency')
687 0 : call addfld( 'Nudge_V',(/ 'lev' /),'A','m/s/s' ,'V Nudging Tendency')
688 0 : call addfld( 'Nudge_T',(/ 'lev' /),'A','K/s' ,'T Nudging Tendency')
689 0 : call addfld( 'Nudge_Q',(/ 'lev' /),'A','kg/kg/s','Q Nudging Tendency')
690 0 : call addfld('Target_U',(/ 'lev' /),'A','m/s' ,'U Nudging Target' )
691 0 : call addfld('Target_V',(/ 'lev' /),'A','m/s' ,'V Nudging Target' )
692 0 : call addfld('Target_T',(/ 'lev' /),'A','K' ,'T Nudging Target' )
693 0 : call addfld('Target_Q',(/ 'lev' /),'A','kg/kg' ,'Q Nudging Target ')
694 :
695 : !-----------------------------------------
696 : ! Values initialized only by masterproc
697 : !-----------------------------------------
698 0 : if(masterproc) then
699 :
700 : ! Set the Stepping intervals for Model and Nudging values
701 : ! Ensure that the Model_Step is not smaller then one timestep
702 : ! and not larger then the Nudge_Step.
703 : !--------------------------------------------------------
704 0 : Model_Step=86400/Model_Times_Per_Day
705 0 : Nudge_Step=86400/Nudge_Times_Per_Day
706 0 : if(Model_Step < dtime) then
707 0 : write(iulog,*) ' '
708 0 : write(iulog,*) 'NUDGING: Model_Step cannot be less than a model timestep'
709 0 : write(iulog,*) 'NUDGING: Setting Model_Step=dtime , dtime=',dtime
710 0 : write(iulog,*) ' '
711 0 : Model_Step=dtime
712 : endif
713 0 : if(Model_Step > Nudge_Step) then
714 0 : write(iulog,*) ' '
715 0 : write(iulog,*) 'NUDGING: Model_Step cannot be more than Nudge_Step'
716 0 : write(iulog,*) 'NUDGING: Setting Model_Step=Nudge_Step, Nudge_Step=',Nudge_Step
717 0 : write(iulog,*) ' '
718 0 : Model_Step=Nudge_Step
719 : endif
720 :
721 : ! Initialize column and level dimensions
722 : !--------------------------------------------------------
723 0 : call get_horiz_grid_dim_d(hdim1_d,hdim2_d)
724 0 : Nudge_nlon=hdim1_d
725 0 : Nudge_nlat=hdim2_d
726 0 : Nudge_ncol=hdim1_d*hdim2_d
727 0 : Nudge_nlev=pver
728 :
729 : ! Check the time relative to the nudging window
730 : !------------------------------------------------
731 0 : call get_curr_date(Year,Month,Day,Sec)
732 0 : YMD=(Year*10000) + (Month*100) + Day
733 0 : YMD1=(Nudge_Beg_Year*10000) + (Nudge_Beg_Month*100) + Nudge_Beg_Day
734 : call timemgr_time_ge(YMD1,Nudge_Beg_Sec, &
735 0 : YMD ,Sec ,After_Beg)
736 0 : YMD1=(Nudge_End_Year*10000) + (Nudge_End_Month*100) + Nudge_End_Day
737 : call timemgr_time_ge(YMD ,Sec , &
738 0 : YMD1,Nudge_End_Sec,Before_End)
739 :
740 0 : if((After_Beg) .and. (Before_End)) then
741 : ! Set Time indicies so that the next call to
742 : ! timestep_init will initialize the data arrays.
743 : !--------------------------------------------
744 0 : Model_Next_Year =Year
745 0 : Model_Next_Month=Month
746 0 : Model_Next_Day =Day
747 0 : Model_Next_Sec =(Sec/Model_Step)*Model_Step
748 0 : Nudge_Next_Year =Year
749 0 : Nudge_Next_Month=Month
750 0 : Nudge_Next_Day =Day
751 0 : Nudge_Next_Sec =(Sec/Nudge_Step)*Nudge_Step
752 0 : elseif(.not.After_Beg) then
753 : ! Set Time indicies to Nudging start,
754 : ! timestep_init will initialize the data arrays.
755 : !--------------------------------------------
756 0 : Model_Next_Year =Nudge_Beg_Year
757 0 : Model_Next_Month=Nudge_Beg_Month
758 0 : Model_Next_Day =Nudge_Beg_Day
759 0 : Model_Next_Sec =Nudge_Beg_Sec
760 0 : Nudge_Next_Year =Nudge_Beg_Year
761 0 : Nudge_Next_Month=Nudge_Beg_Month
762 0 : Nudge_Next_Day =Nudge_Beg_Day
763 0 : Nudge_Next_Sec =Nudge_Beg_Sec
764 0 : elseif(.not.Before_End) then
765 : ! Nudging will never occur, so switch it off
766 : !--------------------------------------------
767 0 : Nudge_Model=.false.
768 0 : Nudge_ON =.false.
769 0 : write(iulog,*) ' '
770 0 : write(iulog,*) 'NUDGING: WARNING - Nudging has been requested by it will'
771 0 : write(iulog,*) 'NUDGING: never occur for the given time values'
772 0 : write(iulog,*) ' '
773 : endif
774 :
775 : ! Initialize values for window function
776 : !----------------------------------------
777 0 : lonp= 180._r8
778 0 : lon0= 0._r8
779 0 : lonn=-180._r8
780 0 : latp= 90._r8-Nudge_Hwin_lat0
781 0 : lat0= 0._r8
782 0 : latn= -90._r8-Nudge_Hwin_lat0
783 :
784 0 : Nudge_Hwin_lonWidthH=Nudge_Hwin_lonWidth/2._r8
785 0 : Nudge_Hwin_latWidthH=Nudge_Hwin_latWidth/2._r8
786 :
787 0 : Val1_p=(1._r8+tanh((Nudge_Hwin_lonWidthH+lonp)/Nudge_Hwin_lonDelta))/2._r8
788 0 : Val2_p=(1._r8+tanh((Nudge_Hwin_lonWidthH-lonp)/Nudge_Hwin_lonDelta))/2._r8
789 0 : Val3_p=(1._r8+tanh((Nudge_Hwin_latWidthH+latp)/Nudge_Hwin_latDelta))/2._r8
790 0 : Val4_p=(1._r8+tanh((Nudge_Hwin_latWidthH-latp)/Nudge_Hwin_latDelta))/2_r8
791 0 : Val1_0=(1._r8+tanh((Nudge_Hwin_lonWidthH+lon0)/Nudge_Hwin_lonDelta))/2._r8
792 0 : Val2_0=(1._r8+tanh((Nudge_Hwin_lonWidthH-lon0)/Nudge_Hwin_lonDelta))/2._r8
793 0 : Val3_0=(1._r8+tanh((Nudge_Hwin_latWidthH+lat0)/Nudge_Hwin_latDelta))/2._r8
794 0 : Val4_0=(1._r8+tanh((Nudge_Hwin_latWidthH-lat0)/Nudge_Hwin_latDelta))/2._r8
795 :
796 0 : Val1_n=(1._r8+tanh((Nudge_Hwin_lonWidthH+lonn)/Nudge_Hwin_lonDelta))/2._r8
797 0 : Val2_n=(1._r8+tanh((Nudge_Hwin_lonWidthH-lonn)/Nudge_Hwin_lonDelta))/2._r8
798 0 : Val3_n=(1._r8+tanh((Nudge_Hwin_latWidthH+latn)/Nudge_Hwin_latDelta))/2._r8
799 0 : Val4_n=(1._r8+tanh((Nudge_Hwin_latWidthH-latn)/Nudge_Hwin_latDelta))/2._r8
800 :
801 0 : Nudge_Hwin_max= Val1_0*Val2_0*Val3_0*Val4_0
802 : Nudge_Hwin_min=min((Val1_p*Val2_p*Val3_n*Val4_n), &
803 : (Val1_p*Val2_p*Val3_p*Val4_p), &
804 : (Val1_n*Val2_n*Val3_n*Val4_n), &
805 0 : (Val1_n*Val2_n*Val3_p*Val4_p))
806 :
807 : ! Initialize number of nudging observation values to keep track of.
808 : ! Allocate and initialize observation indices
809 : !-----------------------------------------------------------------
810 0 : if((Nudge_Force_Opt >= 0) .and. (Nudge_Force_Opt <= 1)) then
811 0 : Nudge_NumObs=2
812 : else
813 : ! Additional Options may need OBS values at more times.
814 : !------------------------------------------------------
815 0 : Nudge_NumObs=2
816 0 : write(iulog,*) 'NUDGING: Setting Nudge_NumObs=2'
817 0 : write(iulog,*) 'NUDGING: WARNING: Unknown Nudge_Force_Opt=',Nudge_Force_Opt
818 0 : call endrun('NUDGING: Unknown Forcing Option')
819 : endif
820 0 : allocate(Nudge_ObsInd(Nudge_NumObs),stat=istat)
821 0 : call alloc_err(istat,'nudging_init','Nudge_ObsInd',Nudge_NumObs)
822 0 : allocate(Nudge_File_Present(Nudge_NumObs),stat=istat)
823 0 : call alloc_err(istat,'nudging_init','Nudge_File_Present',Nudge_NumObs)
824 0 : do nn=1,Nudge_NumObs
825 0 : Nudge_ObsInd(nn) = Nudge_NumObs+1-nn
826 : end do
827 0 : Nudge_File_Present(:) = .false.
828 :
829 : ! Initialization is done,
830 : !--------------------------
831 0 : Nudge_Initialized = .true.
832 :
833 : ! Informational Output
834 : !---------------------------
835 0 : write(iulog,*) ' '
836 0 : write(iulog,*) '---------------------------------------------------------'
837 0 : write(iulog,*) ' MODEL NUDGING INITIALIZED WITH THE FOLLOWING SETTINGS: '
838 0 : write(iulog,*) '---------------------------------------------------------'
839 0 : write(iulog,*) 'NUDGING: Nudge_Model=',Nudge_Model
840 0 : write(iulog,*) 'NUDGING: Nudge_Path=',Nudge_Path
841 0 : write(iulog,*) 'NUDGING: Nudge_File_Template =',Nudge_File_Template
842 0 : write(iulog,*) 'NUDGING: Nudge_Force_Opt=',Nudge_Force_Opt
843 0 : write(iulog,*) 'NUDGING: Nudge_TimeScale_Opt=',Nudge_TimeScale_Opt
844 0 : write(iulog,*) 'NUDGING: Nudge_TSmode=',Nudge_TSmode
845 0 : write(iulog,*) 'NUDGING: Nudge_Times_Per_Day=',Nudge_Times_Per_Day
846 0 : write(iulog,*) 'NUDGING: Model_Times_Per_Day=',Model_Times_Per_Day
847 0 : write(iulog,*) 'NUDGING: Nudge_Step=',Nudge_Step
848 0 : write(iulog,*) 'NUDGING: Model_Step=',Model_Step
849 0 : write(iulog,*) 'NUDGING: Nudge_ZonalFilter=',Nudge_ZonalFilter
850 0 : write(iulog,*) 'NUDGING: Nudge_ZonalNbasis=',Nudge_ZonalNbasis
851 0 : write(iulog,*) 'NUDGING: Nudge_Ucoef =',Nudge_Ucoef
852 0 : write(iulog,*) 'NUDGING: Nudge_Vcoef =',Nudge_Vcoef
853 0 : write(iulog,*) 'NUDGING: Nudge_Qcoef =',Nudge_Qcoef
854 0 : write(iulog,*) 'NUDGING: Nudge_Tcoef =',Nudge_Tcoef
855 0 : write(iulog,*) 'NUDGING: Nudge_PScoef =',Nudge_PScoef
856 0 : write(iulog,*) 'NUDGING: Nudge_Uprof =',Nudge_Uprof
857 0 : write(iulog,*) 'NUDGING: Nudge_Vprof =',Nudge_Vprof
858 0 : write(iulog,*) 'NUDGING: Nudge_Qprof =',Nudge_Qprof
859 0 : write(iulog,*) 'NUDGING: Nudge_Tprof =',Nudge_Tprof
860 0 : write(iulog,*) 'NUDGING: Nudge_PSprof =',Nudge_PSprof
861 0 : write(iulog,*) 'NUDGING: Nudge_Beg_Year =',Nudge_Beg_Year
862 0 : write(iulog,*) 'NUDGING: Nudge_Beg_Month=',Nudge_Beg_Month
863 0 : write(iulog,*) 'NUDGING: Nudge_Beg_Day =',Nudge_Beg_Day
864 0 : write(iulog,*) 'NUDGING: Nudge_End_Year =',Nudge_End_Year
865 0 : write(iulog,*) 'NUDGING: Nudge_End_Month=',Nudge_End_Month
866 0 : write(iulog,*) 'NUDGING: Nudge_End_Day =',Nudge_End_Day
867 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_lat0 =',Nudge_Hwin_lat0
868 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_latWidth =',Nudge_Hwin_latWidth
869 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_latDelta =',Nudge_Hwin_latDelta
870 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_lon0 =',Nudge_Hwin_lon0
871 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_lonWidth =',Nudge_Hwin_lonWidth
872 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_lonDelta =',Nudge_Hwin_lonDelta
873 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_Invert =',Nudge_Hwin_Invert
874 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_lo =',Nudge_Hwin_lo
875 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_hi =',Nudge_Hwin_hi
876 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_Hindex =',Nudge_Vwin_Hindex
877 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_Hdelta =',Nudge_Vwin_Hdelta
878 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_Lindex =',Nudge_Vwin_Lindex
879 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_Ldelta =',Nudge_Vwin_Ldelta
880 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_Invert =',Nudge_Vwin_Invert
881 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_lo =',Nudge_Vwin_lo
882 0 : write(iulog,*) 'NUDGING: Nudge_Vwin_hi =',Nudge_Vwin_hi
883 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_latWidthH=',Nudge_Hwin_latWidthH
884 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_lonWidthH=',Nudge_Hwin_lonWidthH
885 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_max =',Nudge_Hwin_max
886 0 : write(iulog,*) 'NUDGING: Nudge_Hwin_min =',Nudge_Hwin_min
887 0 : write(iulog,*) 'NUDGING: Nudge_Initialized =',Nudge_Initialized
888 0 : write(iulog,*) ' '
889 0 : write(iulog,*) 'NUDGING: Nudge_NumObs=',Nudge_NumObs
890 0 : write(iulog,*) ' '
891 :
892 : endif ! (masterproc) then
893 :
894 : ! Broadcast other variables that have changed
895 : !---------------------------------------------
896 0 : call MPI_bcast(Model_Step , 1, mpi_real8 , mstrid, mpicom, ierr)
897 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Model_Step')
898 0 : call MPI_bcast(Nudge_Step , 1, mpi_real8 , mstrid, mpicom, ierr)
899 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Step')
900 0 : call MPI_bcast(Model_Next_Year , 1, mpi_integer, mstrid, mpicom, ierr)
901 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Model_Next_Year')
902 0 : call MPI_bcast(Model_Next_Month , 1, mpi_integer, mstrid, mpicom, ierr)
903 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Model_Next_Month')
904 0 : call MPI_bcast(Model_Next_Day , 1, mpi_integer, mstrid, mpicom, ierr)
905 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Model_Next_Day')
906 0 : call MPI_bcast(Model_Next_Sec , 1, mpi_integer, mstrid, mpicom, ierr)
907 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Model_Next_Sec')
908 0 : call MPI_bcast(Nudge_Next_Year , 1, mpi_integer, mstrid, mpicom, ierr)
909 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Next_Year')
910 0 : call MPI_bcast(Nudge_Next_Month , 1, mpi_integer, mstrid, mpicom, ierr)
911 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Next_Month')
912 0 : call MPI_bcast(Nudge_Next_Day , 1, mpi_integer, mstrid, mpicom, ierr)
913 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Next_Day')
914 0 : call MPI_bcast(Nudge_Next_Sec , 1, mpi_integer, mstrid, mpicom, ierr)
915 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Next_Sec')
916 0 : call MPI_bcast(Nudge_Model , 1, mpi_logical, mstrid, mpicom, ierr)
917 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Model')
918 0 : call MPI_bcast(Nudge_ON , 1, mpi_logical, mstrid, mpicom, ierr)
919 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_ON')
920 0 : call MPI_bcast(Nudge_Initialized , 1, mpi_logical, mstrid, mpicom, ierr)
921 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Initialized')
922 0 : call MPI_bcast(Nudge_ncol , 1, mpi_integer, mstrid, mpicom, ierr)
923 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_ncol')
924 0 : call MPI_bcast(Nudge_nlev , 1, mpi_integer, mstrid, mpicom, ierr)
925 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_nlev')
926 0 : call MPI_bcast(Nudge_nlon , 1, mpi_integer, mstrid, mpicom, ierr)
927 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_nlon')
928 0 : call MPI_bcast(Nudge_nlat , 1, mpi_integer, mstrid, mpicom, ierr)
929 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_nlat')
930 0 : call MPI_bcast(Nudge_Hwin_max , 1, mpi_real8 , mstrid, mpicom, ierr)
931 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_max')
932 0 : call MPI_bcast(Nudge_Hwin_min , 1, mpi_real8 , mstrid, mpicom, ierr)
933 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_min')
934 0 : call MPI_bcast(Nudge_Hwin_lonWidthH, 1, mpi_real8 , mstrid, mpicom, ierr)
935 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_lonWidthH')
936 0 : call MPI_bcast(Nudge_Hwin_latWidthH, 1, mpi_real8 , mstrid, mpicom, ierr)
937 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_Hwin_latWidthH')
938 0 : call MPI_bcast(Nudge_NumObs , 1, mpi_integer, mstrid, mpicom, ierr)
939 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_NumObs')
940 :
941 : ! All non-masterproc processes also need to allocate space
942 : ! before the broadcast of Nudge_NumObs dependent data.
943 : !------------------------------------------------------------
944 0 : if(.not.masterproc) then
945 0 : allocate(Nudge_ObsInd(Nudge_NumObs),stat=istat)
946 0 : call alloc_err(istat,'nudging_init','Nudge_ObsInd',Nudge_NumObs)
947 0 : allocate(Nudge_File_Present(Nudge_NumObs),stat=istat)
948 0 : call alloc_err(istat,'nudging_init','Nudge_File_Present',Nudge_NumObs)
949 : endif
950 :
951 0 : call MPI_bcast(Nudge_ObsInd , Nudge_NumObs, mpi_integer, mstrid, mpicom, ierr)
952 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: ')
953 0 : call MPI_bcast(Nudge_File_Present, Nudge_NumObs, mpi_logical, mstrid, mpicom, ierr)
954 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: ')
955 :
956 : ! Allocate Space for Nudging observation arrays, initialize with 0's
957 : !---------------------------------------------------------------------
958 0 : allocate(Nobs_U(pcols,pver,begchunk:endchunk,Nudge_NumObs),stat=istat)
959 0 : call alloc_err(istat,'nudging_init','Nobs_U',pcols*pver*((endchunk-begchunk)+1)*Nudge_NumObs)
960 0 : allocate(Nobs_V(pcols,pver,begchunk:endchunk,Nudge_NumObs),stat=istat)
961 0 : call alloc_err(istat,'nudging_init','Nobs_V',pcols*pver*((endchunk-begchunk)+1)*Nudge_NumObs)
962 0 : allocate(Nobs_T(pcols,pver,begchunk:endchunk,Nudge_NumObs),stat=istat)
963 0 : call alloc_err(istat,'nudging_init','Nobs_T',pcols*pver*((endchunk-begchunk)+1)*Nudge_NumObs)
964 0 : allocate(Nobs_Q(pcols,pver,begchunk:endchunk,Nudge_NumObs),stat=istat)
965 0 : call alloc_err(istat,'nudging_init','Nobs_Q',pcols*pver*((endchunk-begchunk)+1)*Nudge_NumObs)
966 0 : allocate(Nobs_PS(pcols,begchunk:endchunk,Nudge_NumObs),stat=istat)
967 0 : call alloc_err(istat,'nudging_init','Nobs_PS',pcols*((endchunk-begchunk)+1)*Nudge_NumObs)
968 :
969 0 : Nobs_U(:pcols,:pver,begchunk:endchunk,:Nudge_NumObs)=0._r8
970 0 : Nobs_V(:pcols,:pver,begchunk:endchunk,:Nudge_NumObs)=0._r8
971 0 : Nobs_T(:pcols,:pver,begchunk:endchunk,:Nudge_NumObs)=0._r8
972 0 : Nobs_Q(:pcols,:pver,begchunk:endchunk,:Nudge_NumObs)=0._r8
973 0 : Nobs_PS(:pcols ,begchunk:endchunk,:Nudge_NumObs)=0._r8
974 :
975 : !!DIAG
976 0 : if(masterproc) then
977 0 : write(iulog,*) 'NUDGING: nudging_init() OBS arrays allocated and initialized'
978 0 : write(iulog,*) 'NUDGING: nudging_init() SIZE#',(9*pcols*pver*((endchunk-begchunk)+1)*Nudge_NumObs)
979 0 : write(iulog,*) 'NUDGING: nudging_init() MB:',float(8*9*pcols*pver*((endchunk-begchunk)+1)*Nudge_NumObs)/(1024._r8*1024._r8)
980 0 : write(iulog,*) 'NUDGING: nudging_init() pcols=',pcols,' pver=',pver
981 0 : write(iulog,*) 'NUDGING: nudging_init() begchunk:',begchunk,' endchunk=',endchunk
982 0 : write(iulog,*) 'NUDGING: nudging_init() chunk:',(endchunk-begchunk+1),' Nudge_NumObs=',Nudge_NumObs
983 0 : write(iulog,*) 'NUDGING: nudging_init() Nudge_ObsInd=',Nudge_ObsInd
984 0 : write(iulog,*) 'NUDGING: nudging_init() Nudge_File_Present=',Nudge_File_Present
985 : endif
986 : !!DIAG
987 :
988 : ! Initialize the Zonal Mean type if needed
989 : !------------------------------------------
990 0 : if(Nudge_ZonalFilter) then
991 0 : call ZM%init(Nudge_ZonalNbasis)
992 0 : allocate(Zonal_Bamp2d(Nudge_ZonalNbasis),stat=istat)
993 0 : call alloc_err(istat,'nudging_init','Zonal_Bamp2d',Nudge_ZonalNbasis)
994 0 : allocate(Zonal_Bamp3d(Nudge_ZonalNbasis,pver),stat=istat)
995 0 : call alloc_err(istat,'nudging_init','Zonal_Bamp3d',Nudge_ZonalNbasis*pver)
996 : endif
997 :
998 : ! Initialize the analysis filename at the NEXT time for startup.
999 : !---------------------------------------------------------------
1000 : Nudge_File=interpret_filename_spec(Nudge_File_Template , &
1001 : yr_spec=Nudge_Next_Year , &
1002 : mon_spec=Nudge_Next_Month, &
1003 : day_spec=Nudge_Next_Day , &
1004 0 : sec_spec=Nudge_Next_Sec )
1005 0 : if(masterproc) then
1006 0 : write(iulog,*) 'NUDGING: Reading analyses:',trim(Nudge_Path)//trim(Nudge_File)
1007 : endif
1008 :
1009 : ! Rotate Nudge_ObsInd() indices for new data, then update
1010 : ! the Nudge observation arrays with analysis data at the
1011 : ! NEXT==Nudge_ObsInd(1) time.
1012 : !----------------------------------------------------------
1013 0 : call nudging_update_analyses (trim(Nudge_Path)//trim(Nudge_File))
1014 :
1015 : ! Initialize Nudging Coeffcient profiles in local arrays
1016 : ! Load zeros into nudging arrays
1017 : !------------------------------------------------------
1018 0 : do lchnk=begchunk,endchunk
1019 0 : ncol=get_ncols_p(lchnk)
1020 0 : do icol=1,ncol
1021 0 : rlat=get_rlat_p(lchnk,icol)*180._r8/SHR_CONST_PI
1022 0 : rlon=get_rlon_p(lchnk,icol)*180._r8/SHR_CONST_PI
1023 :
1024 0 : call nudging_set_profile(rlat,rlon,Nudge_Uprof,Wprof,pver)
1025 0 : Nudge_Utau(icol,:,lchnk)=Wprof(:)
1026 0 : call nudging_set_profile(rlat,rlon,Nudge_Vprof,Wprof,pver)
1027 0 : Nudge_Vtau(icol,:,lchnk)=Wprof(:)
1028 0 : call nudging_set_profile(rlat,rlon,Nudge_Tprof,Wprof,pver)
1029 0 : Nudge_Stau(icol,:,lchnk)=Wprof(:)
1030 0 : call nudging_set_profile(rlat,rlon,Nudge_Qprof,Wprof,pver)
1031 0 : Nudge_Qtau(icol,:,lchnk)=Wprof(:)
1032 :
1033 0 : Nudge_PStau(icol,lchnk)=nudging_set_PSprofile(rlat,rlon,Nudge_PSprof)
1034 : end do
1035 : Nudge_Utau(:ncol,:pver,lchnk) = &
1036 0 : Nudge_Utau(:ncol,:pver,lchnk) * Nudge_Ucoef/float(Nudge_Step)
1037 : Nudge_Vtau(:ncol,:pver,lchnk) = &
1038 0 : Nudge_Vtau(:ncol,:pver,lchnk) * Nudge_Vcoef/float(Nudge_Step)
1039 : Nudge_Stau(:ncol,:pver,lchnk) = &
1040 0 : Nudge_Stau(:ncol,:pver,lchnk) * Nudge_Tcoef/float(Nudge_Step)
1041 : Nudge_Qtau(:ncol,:pver,lchnk) = &
1042 0 : Nudge_Qtau(:ncol,:pver,lchnk) * Nudge_Qcoef/float(Nudge_Step)
1043 : Nudge_PStau(:ncol,lchnk)= &
1044 0 : Nudge_PStau(:ncol,lchnk)* Nudge_PScoef/float(Nudge_Step)
1045 :
1046 0 : Nudge_Ustep(:pcols,:pver,lchnk)=0._r8
1047 0 : Nudge_Vstep(:pcols,:pver,lchnk)=0._r8
1048 0 : Nudge_Sstep(:pcols,:pver,lchnk)=0._r8
1049 0 : Nudge_Qstep(:pcols,:pver,lchnk)=0._r8
1050 0 : Nudge_PSstep(:pcols,lchnk)=0._r8
1051 0 : Target_U(:pcols,:pver,lchnk)=0._r8
1052 0 : Target_V(:pcols,:pver,lchnk)=0._r8
1053 0 : Target_T(:pcols,:pver,lchnk)=0._r8
1054 0 : Target_S(:pcols,:pver,lchnk)=0._r8
1055 0 : Target_Q(:pcols,:pver,lchnk)=0._r8
1056 0 : Target_PS(:pcols,lchnk)=0._r8
1057 : end do
1058 :
1059 : ! End Routine
1060 : !------------
1061 :
1062 0 : end subroutine nudging_init
1063 : !================================================================
1064 :
1065 :
1066 : !================================================================
1067 0 : subroutine nudging_timestep_init(phys_state)
1068 : !
1069 : ! NUDGING_TIMESTEP_INIT:
1070 : ! Check the current time and update Model/Nudging
1071 : ! arrays when necessary. Toggle the Nudging flag
1072 : ! when the time is withing the nudging window.
1073 : !===============================================================
1074 0 : use physconst ,only: cpair
1075 : use physics_types,only: physics_state
1076 : use constituents ,only: cnst_get_ind
1077 : use dycore ,only: dycore_is
1078 : use ppgrid ,only: pver,pcols,begchunk,endchunk
1079 : use filenames ,only: interpret_filename_spec
1080 : use ESMF
1081 :
1082 : ! Arguments
1083 : !-----------
1084 : type(physics_state),intent(in):: phys_state(begchunk:endchunk)
1085 :
1086 : ! Local values
1087 : !----------------
1088 : integer :: Year,Month,Day,Sec
1089 : integer :: YMD1,YMD2,YMD
1090 : logical :: Update_Model,Update_Nudge,Sync_Error
1091 : logical :: After_Beg ,Before_End
1092 : integer :: lchnk,ncol,indw
1093 :
1094 : type(ESMF_Time) :: Date1,Date2
1095 : type(ESMF_TimeInterval) :: DateDiff
1096 : integer :: DeltaT
1097 : real(r8) :: Tscale
1098 : real(r8) :: Tfrac
1099 : integer :: rc
1100 : integer :: nn
1101 : integer :: kk
1102 : real(r8) :: Sbar,Qbar,Wsum
1103 : integer :: dtime
1104 :
1105 : ! Check if Nudging is initialized
1106 : !---------------------------------
1107 0 : if(.not.Nudge_Initialized) then
1108 0 : call endrun('nudging_timestep_init:: Nudging NOT Initialized')
1109 : endif
1110 :
1111 : ! Get time step size
1112 : !--------------------
1113 0 : dtime = get_step_size()
1114 :
1115 : ! Get Current time
1116 : !--------------------
1117 0 : call get_curr_date(Year,Month,Day,Sec)
1118 0 : YMD=(Year*10000) + (Month*100) + Day
1119 :
1120 : !-------------------------------------------------------
1121 : ! Determine if the current time is AFTER the begining time
1122 : ! and if it is BEFORE the ending time.
1123 : !-------------------------------------------------------
1124 0 : YMD1=(Nudge_Beg_Year*10000) + (Nudge_Beg_Month*100) + Nudge_Beg_Day
1125 : call timemgr_time_ge(YMD1,Nudge_Beg_Sec, &
1126 0 : YMD ,Sec ,After_Beg)
1127 :
1128 0 : YMD1=(Nudge_End_Year*10000) + (Nudge_End_Month*100) + Nudge_End_Day
1129 : call timemgr_time_ge(YMD ,Sec, &
1130 0 : YMD1,Nudge_End_Sec,Before_End)
1131 :
1132 : !--------------------------------------------------------------
1133 : ! When past the NEXT time, Update Model Arrays and time indices
1134 : !--------------------------------------------------------------
1135 0 : YMD1=(Model_Next_Year*10000) + (Model_Next_Month*100) + Model_Next_Day
1136 : call timemgr_time_ge(YMD1,Model_Next_Sec, &
1137 0 : YMD ,Sec ,Update_Model)
1138 :
1139 0 : if((Before_End) .and. (Update_Model)) then
1140 : ! Increment the Model times by the current interval
1141 : !---------------------------------------------------
1142 0 : Model_Curr_Year =Model_Next_Year
1143 0 : Model_Curr_Month=Model_Next_Month
1144 0 : Model_Curr_Day =Model_Next_Day
1145 0 : Model_Curr_Sec =Model_Next_Sec
1146 0 : YMD1=(Model_Curr_Year*10000) + (Model_Curr_Month*100) + Model_Curr_Day
1147 : call timemgr_time_inc(YMD1,Model_Curr_Sec, &
1148 0 : YMD2,Model_Next_Sec,Model_Step,0,0)
1149 :
1150 : ! Check for Sync Error where NEXT model time after the update
1151 : ! is before the current time. If so, reset the next model
1152 : ! time to a Model_Step after the current time.
1153 : !--------------------------------------------------------------
1154 : call timemgr_time_ge(YMD2,Model_Next_Sec, &
1155 0 : YMD ,Sec ,Sync_Error)
1156 0 : if(Sync_Error) then
1157 0 : Model_Curr_Year =Year
1158 0 : Model_Curr_Month=Month
1159 0 : Model_Curr_Day =Day
1160 0 : Model_Curr_Sec =Sec
1161 : call timemgr_time_inc(YMD ,Model_Curr_Sec, &
1162 0 : YMD2,Model_Next_Sec,Model_Step,0,0)
1163 0 : write(iulog,*) 'NUDGING: WARNING - Model_Time Sync ERROR... CORRECTED'
1164 : endif
1165 0 : Model_Next_Year =(YMD2/10000)
1166 0 : YMD2 = YMD2-(Model_Next_Year*10000)
1167 0 : Model_Next_Month=(YMD2/100)
1168 0 : Model_Next_Day = YMD2-(Model_Next_Month*100)
1169 :
1170 : ! Load values at Current into the Model arrays
1171 : !-----------------------------------------------
1172 0 : call cnst_get_ind('Q',indw)
1173 0 : do lchnk=begchunk,endchunk
1174 0 : ncol=phys_state(lchnk)%ncol
1175 0 : Model_U(:ncol,:pver,lchnk)=phys_state(lchnk)%u(:ncol,:pver)
1176 0 : Model_V(:ncol,:pver,lchnk)=phys_state(lchnk)%v(:ncol,:pver)
1177 0 : Model_T(:ncol,:pver,lchnk)=phys_state(lchnk)%t(:ncol,:pver)
1178 0 : Model_Q(:ncol,:pver,lchnk)=phys_state(lchnk)%q(:ncol,:pver,indw)
1179 0 : Model_PS(:ncol,lchnk)=phys_state(lchnk)%ps(:ncol)
1180 : end do
1181 :
1182 : ! Load Dry Static Energy values for Model
1183 : !-----------------------------------------
1184 0 : if(Nudge_TSmode == 0) then
1185 : ! DSE tendencies from Temperature only
1186 : !---------------------------------------
1187 0 : do lchnk=begchunk,endchunk
1188 0 : ncol=phys_state(lchnk)%ncol
1189 0 : Model_S(:ncol,:pver,lchnk)=cpair*Model_T(:ncol,:pver,lchnk)
1190 : end do
1191 0 : elseif(Nudge_TSmode == 1) then
1192 : ! Caluculate DSE tendencies from Temperature, Water Vapor, and Surface Pressure
1193 : !------------------------------------------------------------------------------
1194 0 : do lchnk=begchunk,endchunk
1195 0 : ncol=phys_state(lchnk)%ncol
1196 : call calc_DryStaticEnergy(Model_T(:,:,lchnk) , Model_Q(:,:,lchnk), &
1197 : phys_state(lchnk)%phis, Model_PS(:,lchnk), &
1198 0 : Model_S(:,:,lchnk), ncol)
1199 : end do
1200 : endif
1201 :
1202 : ! Optionally: Apply Zonal Filtering to Model state data
1203 : !-------------------------------------------------------
1204 0 : if(Nudge_ZonalFilter) then
1205 0 : call ZM%calc_amps(Model_U,Zonal_Bamp3d)
1206 0 : call ZM%eval_grid(Zonal_Bamp3d,Model_U)
1207 :
1208 0 : call ZM%calc_amps(Model_V,Zonal_Bamp3d)
1209 0 : call ZM%eval_grid(Zonal_Bamp3d,Model_V)
1210 :
1211 0 : call ZM%calc_amps(Model_T,Zonal_Bamp3d)
1212 0 : call ZM%eval_grid(Zonal_Bamp3d,Model_T)
1213 :
1214 0 : call ZM%calc_amps(Model_S,Zonal_Bamp3d)
1215 0 : call ZM%eval_grid(Zonal_Bamp3d,Model_S)
1216 :
1217 0 : call ZM%calc_amps(Model_Q,Zonal_Bamp3d)
1218 0 : call ZM%eval_grid(Zonal_Bamp3d,Model_Q)
1219 :
1220 0 : call ZM%calc_amps(Model_PS,Zonal_Bamp2d)
1221 0 : call ZM%eval_grid(Zonal_Bamp2d,Model_PS)
1222 : endif
1223 : endif ! ((Before_End) .and. (Update_Model)) then
1224 :
1225 : !----------------------------------------------------------------
1226 : ! When past the NEXT time, Update Nudging Arrays and time indices
1227 : !----------------------------------------------------------------
1228 0 : YMD1=(Nudge_Next_Year*10000) + (Nudge_Next_Month*100) + Nudge_Next_Day
1229 : call timemgr_time_ge(YMD1,Nudge_Next_Sec, &
1230 0 : YMD ,Sec ,Update_Nudge)
1231 :
1232 0 : if((Before_End) .and. (Update_Nudge)) then
1233 : ! Increment the Nudge times by the current interval
1234 : !---------------------------------------------------
1235 0 : Nudge_Curr_Year =Nudge_Next_Year
1236 0 : Nudge_Curr_Month=Nudge_Next_Month
1237 0 : Nudge_Curr_Day =Nudge_Next_Day
1238 0 : Nudge_Curr_Sec =Nudge_Next_Sec
1239 0 : YMD1=(Nudge_Curr_Year*10000) + (Nudge_Curr_Month*100) + Nudge_Curr_Day
1240 : call timemgr_time_inc(YMD1,Nudge_Curr_Sec, &
1241 0 : YMD2,Nudge_Next_Sec,Nudge_Step,0,0)
1242 0 : Nudge_Next_Year =(YMD2/10000)
1243 0 : YMD2 = YMD2-(Nudge_Next_Year*10000)
1244 0 : Nudge_Next_Month=(YMD2/100)
1245 0 : Nudge_Next_Day = YMD2-(Nudge_Next_Month*100)
1246 :
1247 : ! Set the analysis filename at the NEXT time.
1248 : !---------------------------------------------------------------
1249 : Nudge_File=interpret_filename_spec(Nudge_File_Template , &
1250 : yr_spec=Nudge_Next_Year , &
1251 : mon_spec=Nudge_Next_Month, &
1252 : day_spec=Nudge_Next_Day , &
1253 0 : sec_spec=Nudge_Next_Sec )
1254 0 : if(masterproc) then
1255 0 : write(iulog,*) 'NUDGING: Reading analyses:',trim(Nudge_Path)//trim(Nudge_File)
1256 : endif
1257 :
1258 : ! Rotate Nudge_ObsInd() indices for new data, then update
1259 : ! the Nudge observation arrays with analysis data at the
1260 : ! NEXT==Nudge_ObsInd(1) time.
1261 : !----------------------------------------------------------
1262 0 : call nudging_update_analyses (trim(Nudge_Path)//trim(Nudge_File))
1263 : endif ! ((Before_End) .and. (Update_Nudge)) then
1264 :
1265 : !----------------------------------------------------------------
1266 : ! Toggle Nudging flag when the time interval is between
1267 : ! beginning and ending times, and all of the analyses files exist.
1268 : !----------------------------------------------------------------
1269 0 : if((After_Beg) .and. (Before_End)) then
1270 0 : if(Nudge_Force_Opt == 0) then
1271 : ! Verify that the NEXT analyses are available
1272 : !---------------------------------------------
1273 0 : Nudge_ON=Nudge_File_Present(Nudge_ObsInd(1))
1274 0 : elseif(Nudge_Force_Opt == 1) then
1275 : ! Verify that the CURR and NEXT analyses are available
1276 : !-----------------------------------------------------
1277 0 : Nudge_ON=(Nudge_File_Present(Nudge_ObsInd(1)) .and. &
1278 0 : Nudge_File_Present(Nudge_ObsInd(2)) )
1279 : else
1280 : ! Verify that the ALL analyses are available
1281 : !---------------------------------------------
1282 0 : Nudge_ON=.true.
1283 0 : do nn=1,Nudge_NumObs
1284 0 : if(.not.Nudge_File_Present(nn)) Nudge_ON=.false.
1285 : end do
1286 : endif
1287 0 : if(.not.Nudge_ON) then
1288 0 : if(masterproc) then
1289 0 : write(iulog,*) 'NUDGING: WARNING - analyses file NOT FOUND. Switching '
1290 0 : write(iulog,*) 'NUDGING: nudging OFF to coast thru the gap. '
1291 : endif
1292 : endif
1293 : else
1294 0 : Nudge_ON=.false.
1295 : endif
1296 :
1297 : !-------------------------------------------------------
1298 : ! HERE Implement time dependence of Nudging Coefs HERE
1299 : !-------------------------------------------------------
1300 :
1301 :
1302 : !---------------------------------------------------
1303 : ! If Data arrays have changed update stepping arrays
1304 : !---------------------------------------------------
1305 0 : if((Before_End) .and. ((Update_Nudge) .or. (Update_Model))) then
1306 :
1307 : ! Now Load the Target values for nudging tendencies
1308 : !---------------------------------------------------
1309 0 : if(Nudge_Force_Opt == 0) then
1310 : ! Target is OBS data at NEXT time
1311 : !----------------------------------
1312 0 : do lchnk=begchunk,endchunk
1313 0 : ncol=phys_state(lchnk)%ncol
1314 0 : Target_U(:ncol,:pver,lchnk)=Nobs_U(:ncol,:pver,lchnk,Nudge_ObsInd(1))
1315 0 : Target_V(:ncol,:pver,lchnk)=Nobs_V(:ncol,:pver,lchnk,Nudge_ObsInd(1))
1316 0 : Target_T(:ncol,:pver,lchnk)=Nobs_T(:ncol,:pver,lchnk,Nudge_ObsInd(1))
1317 0 : Target_Q(:ncol,:pver,lchnk)=Nobs_Q(:ncol,:pver,lchnk,Nudge_ObsInd(1))
1318 0 : Target_PS(:ncol ,lchnk)=Nobs_PS(:ncol ,lchnk,Nudge_ObsInd(1))
1319 : end do
1320 0 : elseif(Nudge_Force_Opt == 1) then
1321 : ! Target is linear interpolation of OBS data CURR<-->NEXT time
1322 : !---------------------------------------------------------------
1323 0 : call ESMF_TimeSet(Date1,YY=Year,MM=Month,DD=Day,S=Sec)
1324 : call ESMF_TimeSet(Date2,YY=Nudge_Next_Year,MM=Nudge_Next_Month, &
1325 0 : DD=Nudge_Next_Day , S=Nudge_Next_Sec )
1326 0 : DateDiff =Date2-Date1
1327 0 : call ESMF_TimeIntervalGet(DateDiff,S=DeltaT,rc=rc)
1328 0 : Tfrac= float(DeltaT)/float(Nudge_Step)
1329 0 : do lchnk=begchunk,endchunk
1330 0 : ncol=phys_state(lchnk)%ncol
1331 0 : Target_U(:ncol,:pver,lchnk)=(1._r8-Tfrac)*Nobs_U(:ncol,:pver,lchnk,Nudge_ObsInd(1)) &
1332 0 : +Tfrac *Nobs_U(:ncol,:pver,lchnk,Nudge_ObsInd(2))
1333 0 : Target_V(:ncol,:pver,lchnk)=(1._r8-Tfrac)*Nobs_V(:ncol,:pver,lchnk,Nudge_ObsInd(1)) &
1334 0 : +Tfrac *Nobs_V(:ncol,:pver,lchnk,Nudge_ObsInd(2))
1335 0 : Target_T(:ncol,:pver,lchnk)=(1._r8-Tfrac)*Nobs_T(:ncol,:pver,lchnk,Nudge_ObsInd(1)) &
1336 0 : +Tfrac *Nobs_T(:ncol,:pver,lchnk,Nudge_ObsInd(2))
1337 0 : Target_Q(:ncol,:pver,lchnk)=(1._r8-Tfrac)*Nobs_Q(:ncol,:pver,lchnk,Nudge_ObsInd(1)) &
1338 0 : +Tfrac *Nobs_Q(:ncol,:pver,lchnk,Nudge_ObsInd(2))
1339 0 : Target_PS(:ncol ,lchnk)=(1._r8-Tfrac)*Nobs_PS(:ncol ,lchnk,Nudge_ObsInd(1)) &
1340 0 : +Tfrac *Nobs_PS(:ncol ,lchnk,Nudge_ObsInd(2))
1341 : end do
1342 : else
1343 0 : write(iulog,*) 'NUDGING: Unknown Nudge_Force_Opt=',Nudge_Force_Opt
1344 0 : call endrun('nudging_timestep_init:: ERROR unknown Nudging_Force_Opt')
1345 : endif
1346 :
1347 : ! Now load Dry Static Energy values for Target
1348 : !---------------------------------------------
1349 0 : if(Nudge_TSmode == 0) then
1350 : ! DSE tendencies from Temperature only
1351 : !---------------------------------------
1352 0 : do lchnk=begchunk,endchunk
1353 0 : ncol=phys_state(lchnk)%ncol
1354 0 : Target_S(:ncol,:pver,lchnk)=cpair*Target_T(:ncol,:pver,lchnk)
1355 : end do
1356 0 : elseif(Nudge_TSmode == 1) then
1357 : ! Caluculate DSE tendencies from Temperature, Water Vapor, and Surface Pressure
1358 : !------------------------------------------------------------------------------
1359 0 : do lchnk=begchunk,endchunk
1360 0 : ncol=phys_state(lchnk)%ncol
1361 : call calc_DryStaticEnergy(Target_T(:,:,lchnk), Target_Q(:,:,lchnk), &
1362 : phys_state(lchnk)%phis, Target_PS(:,lchnk), &
1363 0 : Target_S(:,:,lchnk), ncol)
1364 : end do
1365 : endif
1366 :
1367 : ! Set Tscale for the specified Forcing Option
1368 : !-----------------------------------------------
1369 0 : if(Nudge_TimeScale_Opt == 0) then
1370 : Tscale=1._r8
1371 0 : elseif(Nudge_TimeScale_Opt == 1) then
1372 0 : call ESMF_TimeSet(Date1,YY=Year,MM=Month,DD=Day,S=Sec)
1373 : call ESMF_TimeSet(Date2,YY=Nudge_Next_Year,MM=Nudge_Next_Month, &
1374 0 : DD=Nudge_Next_Day , S=Nudge_Next_Sec )
1375 0 : DateDiff =Date2-Date1
1376 0 : call ESMF_TimeIntervalGet(DateDiff,S=DeltaT,rc=rc)
1377 0 : Tscale=float(Nudge_Step)/float(DeltaT)
1378 : else
1379 0 : write(iulog,*) 'NUDGING: Unknown Nudge_TimeScale_Opt=',Nudge_TimeScale_Opt
1380 0 : call endrun('nudging_timestep_init:: ERROR unknown Nudging_TimeScale_Opt')
1381 : endif
1382 :
1383 : ! Update the nudging tendencies
1384 : !--------------------------------
1385 0 : do lchnk=begchunk,endchunk
1386 0 : ncol=phys_state(lchnk)%ncol
1387 0 : Nudge_Ustep(:ncol,:pver,lchnk)=( Target_U(:ncol,:pver,lchnk) &
1388 0 : -Model_U(:ncol,:pver,lchnk)) &
1389 0 : *Tscale*Nudge_Utau(:ncol,:pver,lchnk)
1390 0 : Nudge_Vstep(:ncol,:pver,lchnk)=( Target_V(:ncol,:pver,lchnk) &
1391 0 : -Model_V(:ncol,:pver,lchnk)) &
1392 0 : *Tscale*Nudge_Vtau(:ncol,:pver,lchnk)
1393 0 : Nudge_Sstep(:ncol,:pver,lchnk)=( Target_S(:ncol,:pver,lchnk) &
1394 0 : -Model_S(:ncol,:pver,lchnk)) &
1395 0 : *Tscale*Nudge_Stau(:ncol,:pver,lchnk)
1396 0 : Nudge_Qstep(:ncol,:pver,lchnk)=( Target_Q(:ncol,:pver,lchnk) &
1397 0 : -Model_Q(:ncol,:pver,lchnk)) &
1398 0 : *Tscale*Nudge_Qtau(:ncol,:pver,lchnk)
1399 0 : Nudge_PSstep(:ncol, lchnk)=( Target_PS(:ncol,lchnk) &
1400 0 : -Model_PS(:ncol,lchnk)) &
1401 0 : *Tscale*Nudge_PStau(:ncol,lchnk)
1402 : end do
1403 :
1404 : !******************
1405 : ! DIAG
1406 : !******************
1407 : ! if(masterproc) then
1408 : ! write(iulog,*) 'PFC: Target_T(1,:pver,begchunk)=',Target_T(1,:pver,begchunk)
1409 : ! write(iulog,*) 'PFC: Model_T(1,:pver,begchunk)=',Model_T(1,:pver,begchunk)
1410 : ! write(iulog,*) 'PFC: Target_S(1,:pver,begchunk)=',Target_S(1,:pver,begchunk)
1411 : ! write(iulog,*) 'PFC: Model_S(1,:pver,begchunk)=',Model_S(1,:pver,begchunk)
1412 : ! write(iulog,*) 'PFC: Target_PS(1,begchunk)=',Target_PS(1,begchunk)
1413 : ! write(iulog,*) 'PFC: Model_PS(1,begchunk)=',Model_PS(1,begchunk)
1414 : ! write(iulog,*) 'PFC: Nudge_Sstep(1,:pver,begchunk)=',Nudge_Sstep(1,:pver,begchunk)
1415 : ! write(iulog,*) 'PFC: Nudge_Xstep arrays updated:'
1416 : ! endif
1417 : endif ! ((Before_End) .and. ((Update_Nudge) .or. (Update_Model))) then
1418 :
1419 : ! End Routine
1420 : !------------
1421 :
1422 0 : end subroutine nudging_timestep_init
1423 : !================================================================
1424 :
1425 :
1426 : !================================================================
1427 0 : subroutine nudging_timestep_tend(phys_state,phys_tend)
1428 : !
1429 : ! NUDGING_TIMESTEP_TEND:
1430 : ! If Nudging is ON, return the Nudging contributions
1431 : ! to forcing using the current contents of the Nudge
1432 : ! arrays. Send output to the cam history module as well.
1433 : !===============================================================
1434 0 : use physconst ,only: cpair
1435 : use physics_types,only: physics_state,physics_ptend,physics_ptend_init
1436 : use constituents ,only: cnst_get_ind,pcnst
1437 : use ppgrid ,only: pver,pcols,begchunk,endchunk
1438 : use cam_history ,only: outfld
1439 :
1440 : ! Arguments
1441 : !-------------
1442 : type(physics_state), intent(in) :: phys_state
1443 : type(physics_ptend), intent(out):: phys_tend
1444 :
1445 : ! Local values
1446 : !--------------------
1447 : integer :: indw,ncol,lchnk
1448 : logical :: lq(pcnst)
1449 :
1450 0 : call cnst_get_ind('Q',indw)
1451 0 : lq(:) =.false.
1452 0 : lq(indw)=.true.
1453 0 : call physics_ptend_init(phys_tend,phys_state%psetcols,'nudging',lu=.true.,lv=.true.,ls=.true.,lq=lq)
1454 :
1455 0 : if(Nudge_ON) then
1456 0 : lchnk=phys_state%lchnk
1457 0 : ncol =phys_state%ncol
1458 0 : phys_tend%u(:ncol,:pver) =Nudge_Ustep(:ncol,:pver,lchnk)
1459 0 : phys_tend%v(:ncol,:pver) =Nudge_Vstep(:ncol,:pver,lchnk)
1460 0 : phys_tend%s(:ncol,:pver) =Nudge_Sstep(:ncol,:pver,lchnk)
1461 0 : phys_tend%q(:ncol,:pver,indw)=Nudge_Qstep(:ncol,:pver,lchnk)
1462 :
1463 0 : call outfld( 'Nudge_U',phys_tend%u ,pcols,lchnk)
1464 0 : call outfld( 'Nudge_V',phys_tend%v ,pcols,lchnk)
1465 0 : call outfld( 'Nudge_T',phys_tend%s/cpair ,pcols,lchnk)
1466 0 : call outfld( 'Nudge_Q',phys_tend%q(1,1,indw) ,pcols,lchnk)
1467 0 : call outfld('Target_U',Target_U(:,:,lchnk),pcols,lchnk)
1468 0 : call outfld('Target_V',Target_V(:,:,lchnk),pcols,lchnk)
1469 0 : call outfld('Target_T',Target_T(:,:,lchnk),pcols,lchnk)
1470 0 : call outfld('Target_Q',Target_Q(:,:,lchnk),pcols,lchnk)
1471 : endif
1472 :
1473 : ! End Routine
1474 : !------------
1475 :
1476 0 : end subroutine nudging_timestep_tend
1477 : !================================================================
1478 :
1479 : !================================================================
1480 0 : subroutine nudging_update_analyses(anal_file)
1481 : !
1482 : ! NUDGING_UPDATE_ANALYSES:
1483 : ! Open the given analyses data file, read in
1484 : ! U,V,T,Q, and PS values and then distribute
1485 : ! the values to all of the chunks.
1486 : !===============================================================
1487 0 : use ppgrid ,only: pcols,pver,begchunk,endchunk
1488 : use cam_pio_utils ,only: cam_pio_openfile
1489 : use pio ,only: PIO_BCAST_ERROR,PIO_INTERNAL_ERROR
1490 : use pio ,only: pio_closefile,pio_seterrorhandling,file_desc_t
1491 : use ncdio_atm ,only: infld
1492 : use cam_grid_support,only: cam_grid_id,cam_grid_get_dim_names,DLEN=>max_hcoordname_len
1493 :
1494 : ! Arguments
1495 : !-------------
1496 : character(len=*),intent(in):: anal_file
1497 :
1498 : ! Local values
1499 : !-------------
1500 : type(file_desc_t) :: fileID
1501 : integer :: nn,Nindex
1502 : logical :: VARflag
1503 : integer :: grid_id
1504 : integer :: ierr
1505 : character(len=DLEN):: dim1name,dim2name
1506 : integer :: err_handling
1507 :
1508 0 : real(r8),allocatable:: Tmp3D(:,:,:)
1509 0 : real(r8),allocatable:: Tmp2D(:,:)
1510 :
1511 : character(len=*), parameter :: prefix = 'nudging_update_analyses: '
1512 :
1513 : ! Rotate Nudge_ObsInd() indices, then check the existence of the analyses
1514 : ! file; broadcast the updated indices and file status to all the other MPI nodes.
1515 : ! If the file is not there, then just return.
1516 : !------------------------------------------------------------------------
1517 0 : if(masterproc) then
1518 0 : Nindex=Nudge_ObsInd(Nudge_NumObs)
1519 0 : do nn=Nudge_NumObs,2,-1
1520 0 : Nudge_ObsInd(nn)=Nudge_ObsInd(nn-1)
1521 : end do
1522 0 : Nudge_ObsInd(1)=Nindex
1523 0 : inquire(FILE=trim(anal_file),EXIST=Nudge_File_Present(Nudge_ObsInd(1)))
1524 0 : write(iulog,*)'NUDGING: Nudge_ObsInd=',Nudge_ObsInd
1525 0 : write(iulog,*)'NUDGING: Nudge_File_Present=',Nudge_File_Present
1526 : endif
1527 :
1528 0 : call MPI_bcast(Nudge_File_Present, Nudge_NumObs, mpi_logical, mstrid, mpicom, ierr)
1529 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_File_Present')
1530 0 : call MPI_bcast(Nudge_ObsInd , Nudge_NumObs, mpi_integer, mstrid, mpicom, ierr)
1531 0 : if (ierr /= mpi_success) call endrun(prefix//'FATAL: mpi_bcast: Nudge_ObsInd')
1532 :
1533 0 : if(.not. Nudge_File_Present(Nudge_ObsInd(1))) then
1534 0 : return
1535 : end if
1536 :
1537 : ! Open the file and get the fileID.
1538 : !-------------------------------------
1539 0 : call cam_pio_openfile(fileID,trim(anal_file),0)
1540 0 : call pio_seterrorhandling(fileID,PIO_BCAST_ERROR,oldmethod=err_handling)
1541 0 : if(masterproc) write(iulog,*)'PIO_OPEN: file=',trim(anal_file)
1542 :
1543 0 : grid_id = cam_grid_id('physgrid')
1544 0 : call cam_grid_get_dim_names(grid_id,dim1name,dim2name)
1545 :
1546 0 : allocate(Tmp3D(pcols,pver,begchunk:endchunk))
1547 0 : allocate(Tmp2D(pcols,begchunk:endchunk))
1548 :
1549 : ! Read in, U,V,T,Q, and PS
1550 : !----------------------------------
1551 : call infld('U',fileID,dim1name,'lev',dim2name, &
1552 : 1,pcols,1,pver,begchunk,endchunk,Tmp3D, &
1553 0 : VARflag,gridname='physgrid',timelevel=1 )
1554 0 : if(VARflag) then
1555 0 : if(Nudge_ZonalFilter) then
1556 0 : call ZM%calc_amps(Tmp3D,Zonal_Bamp3d)
1557 0 : call ZM%eval_grid(Zonal_Bamp3d,Tmp3D)
1558 : endif
1559 0 : Nobs_U(:,:,begchunk:endchunk,Nudge_ObsInd(1)) = Tmp3D(:,:,begchunk:endchunk)
1560 : else
1561 0 : call endrun('Variable "U" is missing in '//trim(anal_file))
1562 : endif
1563 :
1564 : call infld('V',fileID,dim1name,'lev',dim2name, &
1565 : 1,pcols,1,pver,begchunk,endchunk,Tmp3D, &
1566 0 : VARflag,gridname='physgrid',timelevel=1 )
1567 0 : if(VARflag) then
1568 0 : if(Nudge_ZonalFilter) then
1569 0 : call ZM%calc_amps(Tmp3D,Zonal_Bamp3d)
1570 0 : call ZM%eval_grid(Zonal_Bamp3d,Tmp3D)
1571 : endif
1572 0 : Nobs_V(:,:,begchunk:endchunk,Nudge_ObsInd(1)) = Tmp3D(:,:,begchunk:endchunk)
1573 : else
1574 0 : call endrun('Variable "V" is missing in '//trim(anal_file))
1575 : endif
1576 :
1577 : call infld('T',fileID,dim1name,'lev',dim2name, &
1578 : 1,pcols,1,pver,begchunk,endchunk,Tmp3D, &
1579 0 : VARflag,gridname='physgrid',timelevel=1 )
1580 0 : if(VARflag) then
1581 0 : if(Nudge_ZonalFilter) then
1582 0 : call ZM%calc_amps(Tmp3D,Zonal_Bamp3d)
1583 0 : call ZM%eval_grid(Zonal_Bamp3d,Tmp3D)
1584 : endif
1585 0 : Nobs_T(:,:,begchunk:endchunk,Nudge_ObsInd(1)) = Tmp3D(:,:,begchunk:endchunk)
1586 : else
1587 0 : call endrun('Variable "T" is missing in '//trim(anal_file))
1588 : endif
1589 :
1590 : call infld('Q',fileID,dim1name,'lev',dim2name, &
1591 : 1,pcols,1,pver,begchunk,endchunk,Tmp3D, &
1592 0 : VARflag,gridname='physgrid',timelevel=1 )
1593 0 : if(VARflag) then
1594 0 : if(Nudge_ZonalFilter) then
1595 0 : call ZM%calc_amps(Tmp3D,Zonal_Bamp3d)
1596 0 : call ZM%eval_grid(Zonal_Bamp3d,Tmp3D)
1597 : endif
1598 0 : Nobs_Q(:,:,begchunk:endchunk,Nudge_ObsInd(1)) = Tmp3D(:,:,begchunk:endchunk)
1599 : else
1600 0 : call endrun('Variable "Q" is missing in '//trim(anal_file))
1601 : endif
1602 :
1603 : call infld('PS',fileID,dim1name,dim2name, &
1604 : 1,pcols,begchunk,endchunk,Tmp2D, &
1605 0 : VARflag,gridname='physgrid',timelevel=1 )
1606 0 : if(VARflag) then
1607 0 : if(Nudge_ZonalFilter) then
1608 0 : call ZM%calc_amps(Tmp2D,Zonal_Bamp2d)
1609 0 : call ZM%eval_grid(Zonal_Bamp2d,Tmp2D)
1610 : endif
1611 0 : Nobs_PS(:,begchunk:endchunk,Nudge_ObsInd(1)) = Tmp2D(:,begchunk:endchunk)
1612 : else
1613 0 : call endrun('Variable "PS" is missing in '//trim(anal_file))
1614 : endif
1615 :
1616 : ! Restore old error handling
1617 : !----------------------------
1618 0 : call pio_seterrorhandling(fileID,err_handling)
1619 :
1620 : ! Close the analyses file
1621 : !-----------------------
1622 0 : deallocate(Tmp3D)
1623 0 : deallocate(Tmp2D)
1624 0 : call pio_closefile(fileID)
1625 :
1626 : ! End Routine
1627 : !------------
1628 :
1629 0 : end subroutine nudging_update_analyses
1630 : !================================================================
1631 :
1632 :
1633 : !================================================================
1634 0 : subroutine nudging_set_profile(rlat,rlon,Nudge_prof,Wprof,nlev)
1635 : !
1636 : ! NUDGING_SET_PROFILE: for the given lat,lon, and Nudging_prof, set
1637 : ! the verical profile of window coeffcients.
1638 : ! Values range from 0. to 1. to affect spatial
1639 : ! variations on nudging strength.
1640 : !===============================================================
1641 :
1642 : ! Arguments
1643 : !--------------
1644 : integer :: nlev,Nudge_prof
1645 : real(r8) :: rlat,rlon
1646 : real(r8) :: Wprof(nlev)
1647 :
1648 : ! Local values
1649 : !----------------
1650 : integer :: ilev
1651 : real(r8) :: Hcoef,latx,lonx,Vmax,Vmin
1652 : real(r8) :: lon_lo,lon_hi,lat_lo,lat_hi,lev_lo,lev_hi
1653 :
1654 : !---------------
1655 : ! set coeffcient
1656 : !---------------
1657 0 : if(Nudge_prof == 0) then
1658 : ! No Nudging
1659 : !-------------
1660 0 : Wprof(:)=0.0_r8
1661 0 : elseif(Nudge_prof == 1) then
1662 : ! Uniform Nudging
1663 : !-----------------
1664 0 : Wprof(:)=1.0_r8
1665 0 : elseif(Nudge_prof == 2) then
1666 : ! Localized Nudging with specified Heaviside window function
1667 : !------------------------------------------------------------
1668 0 : if(Nudge_Hwin_max <= Nudge_Hwin_min) then
1669 : ! For a constant Horizontal window function,
1670 : ! just set Hcoef to the maximum of Hlo/Hhi.
1671 : !--------------------------------------------
1672 0 : Hcoef=max(Nudge_Hwin_lo,Nudge_Hwin_hi)
1673 : else
1674 : ! get lat/lon relative to window center
1675 : !------------------------------------------
1676 0 : latx=rlat-Nudge_Hwin_lat0
1677 0 : lonx=rlon-Nudge_Hwin_lon0
1678 0 : if(lonx > 180._r8) lonx=lonx-360._r8
1679 0 : if(lonx <= -180._r8) lonx=lonx+360._r8
1680 :
1681 : ! Calcualte RAW window value
1682 : !-------------------------------
1683 0 : lon_lo=(Nudge_Hwin_lonWidthH+lonx)/Nudge_Hwin_lonDelta
1684 0 : lon_hi=(Nudge_Hwin_lonWidthH-lonx)/Nudge_Hwin_lonDelta
1685 0 : lat_lo=(Nudge_Hwin_latWidthH+latx)/Nudge_Hwin_latDelta
1686 0 : lat_hi=(Nudge_Hwin_latWidthH-latx)/Nudge_Hwin_latDelta
1687 : Hcoef=((1._r8+tanh(lon_lo))/2._r8)*((1._r8+tanh(lon_hi))/2._r8) &
1688 0 : *((1._r8+tanh(lat_lo))/2._r8)*((1._r8+tanh(lat_hi))/2._r8)
1689 :
1690 : ! Scale the horizontal window coef for specfied range of values.
1691 : !--------------------------------------------------------
1692 0 : Hcoef=(Hcoef-Nudge_Hwin_min)/(Nudge_Hwin_max-Nudge_Hwin_min)
1693 0 : Hcoef=(1._r8-Hcoef)*Nudge_Hwin_lo + Hcoef*Nudge_Hwin_hi
1694 : endif
1695 :
1696 : ! Load the RAW vertical window
1697 : !------------------------------
1698 0 : do ilev=1,nlev
1699 0 : lev_lo=(float(ilev)-Nudge_Vwin_Lindex)/Nudge_Vwin_Ldelta
1700 0 : lev_hi=(Nudge_Vwin_Hindex-float(ilev))/Nudge_Vwin_Hdelta
1701 0 : Wprof(ilev)=((1._r8+tanh(lev_lo))/2._r8)*((1._r8+tanh(lev_hi))/2._r8)
1702 : end do
1703 :
1704 : ! Scale the Window function to span the values between Vlo and Vhi:
1705 : !-----------------------------------------------------------------
1706 0 : Vmax=maxval(Wprof)
1707 0 : Vmin=minval(Wprof)
1708 0 : if((Vmax <= Vmin) .or. ((Nudge_Vwin_Hindex >= (nlev+1)) .and. &
1709 : (Nudge_Vwin_Lindex <= 0 ) )) then
1710 : ! For a constant Vertical window function,
1711 : ! load maximum of Vlo/Vhi into Wprof()
1712 : !--------------------------------------------
1713 0 : Vmax=max(Nudge_Vwin_lo,Nudge_Vwin_hi)
1714 0 : Wprof(:)=Vmax
1715 : else
1716 : ! Scale the RAW vertical window for specfied range of values.
1717 : !--------------------------------------------------------
1718 0 : Wprof(:)=(Wprof(:)-Vmin)/(Vmax-Vmin)
1719 0 : Wprof(:)=Nudge_Vwin_lo + Wprof(:)*(Nudge_Vwin_hi-Nudge_Vwin_lo)
1720 : endif
1721 :
1722 : ! The desired result is the product of the vertical profile
1723 : ! and the horizontal window coeffcient.
1724 : !----------------------------------------------------
1725 0 : Wprof(:)=Hcoef*Wprof(:)
1726 : else
1727 0 : call endrun('nudging_set_profile:: Unknown Nudge_prof value')
1728 : endif
1729 :
1730 : ! End Routine
1731 : !------------
1732 :
1733 0 : end subroutine nudging_set_profile
1734 : !================================================================
1735 :
1736 : !================================================================
1737 0 : subroutine nudging_final
1738 :
1739 0 : if (allocated(Target_U)) deallocate(Target_U)
1740 0 : if (allocated(Target_V)) deallocate(Target_V)
1741 0 : if (allocated(Target_T)) deallocate(Target_T)
1742 0 : if (allocated(Target_S)) deallocate(Target_S)
1743 0 : if (allocated(Target_Q)) deallocate(Target_Q)
1744 0 : if (allocated(Target_PS)) deallocate(Target_PS)
1745 0 : if (allocated(Model_U)) deallocate(Model_U)
1746 0 : if (allocated(Model_V)) deallocate(Model_V)
1747 0 : if (allocated(Model_T)) deallocate(Model_T)
1748 0 : if (allocated(Model_S)) deallocate(Model_S)
1749 0 : if (allocated(Model_Q)) deallocate(Model_Q)
1750 0 : if (allocated(Model_PS)) deallocate(Model_PS)
1751 0 : if (allocated(Nudge_Utau)) deallocate(Nudge_Utau)
1752 0 : if (allocated(Nudge_Vtau)) deallocate(Nudge_Vtau)
1753 0 : if (allocated(Nudge_Stau)) deallocate(Nudge_Stau)
1754 0 : if (allocated(Nudge_Qtau)) deallocate(Nudge_Qtau)
1755 0 : if (allocated(Nudge_PStau)) deallocate(Nudge_PStau)
1756 0 : if (allocated(Nudge_Ustep)) deallocate(Nudge_Ustep)
1757 0 : if (allocated(Nudge_Vstep)) deallocate(Nudge_Vstep)
1758 0 : if (allocated(Nudge_Sstep)) deallocate(Nudge_Sstep)
1759 0 : if (allocated(Nudge_Qstep)) deallocate(Nudge_Qstep)
1760 0 : if (allocated(Nudge_PSstep)) deallocate(Nudge_PSstep)
1761 :
1762 0 : if (allocated(Nudge_ObsInd)) deallocate(Nudge_ObsInd)
1763 0 : if (allocated(Nudge_File_Present)) deallocate(Nudge_File_Present)
1764 0 : if (allocated(Nobs_U)) deallocate(Nobs_U)
1765 0 : if (allocated(Nobs_V)) deallocate(Nobs_V)
1766 0 : if (allocated(Nobs_T)) deallocate(Nobs_T)
1767 0 : if (allocated(Nobs_Q)) deallocate(Nobs_Q)
1768 0 : if (allocated(Nobs_PS)) deallocate(Nobs_PS)
1769 0 : if (allocated(Zonal_Bamp2d)) deallocate(Zonal_Bamp2d)
1770 0 : if (allocated(Zonal_Bamp3d)) deallocate(Zonal_Bamp3d)
1771 :
1772 0 : call ZM%final()
1773 :
1774 0 : end subroutine nudging_final
1775 : !================================================================
1776 :
1777 : !================================================================
1778 0 : real(r8) function nudging_set_PSprofile(rlat,rlon,Nudge_PSprof)
1779 : !
1780 : ! NUDGING_SET_PSPROFILE: for the given lat and lon set the surface
1781 : ! pressure profile value for the specified index.
1782 : ! Values range from 0. to 1. to affect spatial
1783 : ! variations on nudging strength.
1784 : !===============================================================
1785 :
1786 : ! Arguments
1787 : !--------------
1788 : real(r8) :: rlat,rlon
1789 : integer :: Nudge_PSprof
1790 :
1791 : ! Local values
1792 : !----------------
1793 :
1794 : !---------------
1795 : ! set coeffcient
1796 : !---------------
1797 0 : if(Nudge_PSprof == 0) then
1798 : ! No Nudging
1799 : !-------------
1800 : nudging_set_PSprofile=0.0_r8
1801 0 : elseif(Nudge_PSprof == 1) then
1802 : ! Uniform Nudging
1803 : !-----------------
1804 : nudging_set_PSprofile=1.0_r8
1805 : else
1806 0 : call endrun('nudging_set_PSprofile:: Unknown Nudge_prof value')
1807 : endif
1808 :
1809 : ! End Routine
1810 : !------------
1811 :
1812 0 : end function nudging_set_PSprofile
1813 : !================================================================
1814 :
1815 :
1816 : !================================================================
1817 0 : subroutine calc_DryStaticEnergy(t, q, phis, ps, dse, ncol)
1818 : !
1819 : ! calc_DryStaticEnergy: Given the temperature, specific humidity, surface pressure,
1820 : ! and surface geopotential for a chunk containing 'ncol' columns,
1821 : ! calculate and return the corresponding dry static energy values.
1822 : !--------------------------------------------------------------------------------------
1823 : use shr_kind_mod, only: r8 => shr_kind_r8
1824 : use ppgrid, only: pver, pverp
1825 : use dycore, only: dycore_is
1826 : use hycoef, only: hyai, hybi, ps0, hyam, hybm
1827 : use physconst, only: zvir, gravit, cpair, rair
1828 : !
1829 : ! Input/Output arguments
1830 : !-----------------------
1831 : integer , intent(in) :: ncol ! Number of columns in chunk
1832 : real(r8), intent(in) :: t(:,:) ! (pcols,pver) - temperature
1833 : real(r8), intent(in) :: q(:,:) ! (pcols,pver) - specific humidity
1834 : real(r8), intent(in) :: ps(:) ! (pcols) - surface pressure
1835 : real(r8), intent(in) :: phis(:) ! (pcols) - surface geopotential
1836 : real(r8), intent(out):: dse(:,:) ! (pcols,pver) - dry static energy
1837 : !
1838 : ! Local variables
1839 : !------------------
1840 : logical :: fvdyn ! finite volume dynamics
1841 : integer :: ii,kk ! Lon, level, level indices
1842 : real(r8) :: tvfac ! Virtual temperature factor
1843 0 : real(r8) :: hkk(ncol) ! diagonal element of hydrostatic matrix
1844 0 : real(r8) :: hkl(ncol) ! off-diagonal element
1845 0 : real(r8) :: pint(ncol,pverp) ! Interface pressures
1846 0 : real(r8) :: pmid(ncol,pver ) ! Midpoint pressures
1847 0 : real(r8) :: zi(ncol,pverp) ! Height above surface at interfaces
1848 0 : real(r8) :: zm(ncol,pver ) ! Geopotential height at mid level
1849 :
1850 : ! Set dynamics flag
1851 : !-------------------
1852 0 : fvdyn = dycore_is ('LR')
1853 :
1854 : ! Load Pressure values and midpoint pressures
1855 : !----------------------------------------------
1856 0 : do kk=1,pverp
1857 0 : do ii=1,ncol
1858 0 : pint(ii,kk)=(hyai(kk)*ps0)+(hybi(kk)*ps(ii))
1859 : end do
1860 : end do
1861 0 : do kk=1,pver
1862 0 : do ii=1,ncol
1863 0 : pmid(ii,kk)=(hyam(kk)*ps0)+(hybm(kk)*ps(ii))
1864 : end do
1865 : end do
1866 :
1867 : ! The surface height is zero by definition.
1868 : !-------------------------------------------
1869 0 : do ii = 1,ncol
1870 0 : zi(ii,pverp) = 0.0_r8
1871 : end do
1872 :
1873 : ! Compute the dry static energy, zi, zm from bottom up
1874 : ! Note, zi(i,k) is the interface above zm(i,k)
1875 : !---------------------------------------------------------
1876 0 : do kk=pver,1,-1
1877 :
1878 : ! First set hydrostatic elements consistent with dynamics
1879 : !--------------------------------------------------------
1880 0 : if(fvdyn) then
1881 0 : do ii=1,ncol
1882 0 : hkl(ii)=log(pint(ii,kk+1))-log(pint(ii,kk))
1883 0 : hkk(ii)=1._r8-(hkl(ii)*pint(ii,kk)/(pint(ii,kk+1)-pint(ii,kk)))
1884 : end do
1885 : else
1886 0 : do ii=1,ncol
1887 0 : hkl(ii)=(pint(ii,kk+1)-pint(ii,kk))/pmid(ii,kk)
1888 0 : hkk(ii)=0.5_r8*hkl(ii)
1889 : end do
1890 : endif
1891 :
1892 : ! Now compute zm, zi, and dse (WACCM-X vars rairv/zairv/cpairv not used!)
1893 : !------------------------------------------------------------------------
1894 0 : do ii=1,ncol
1895 0 : tvfac=t(ii,kk)*rair*(1._r8+(zvir*q(ii,kk)))/gravit
1896 0 : zm (ii,kk)=zi(ii,kk+1) + (tvfac*hkk(ii))
1897 0 : zi (ii,kk)=zi(ii,kk+1) + (tvfac*hkl(ii))
1898 0 : dse(ii,kk)=(t(ii,kk)*cpair) + phis(ii) + (gravit*zm(ii,kk))
1899 : end do
1900 :
1901 : end do ! kk=pver,1,-1
1902 :
1903 : ! End Routine
1904 : !-----------
1905 :
1906 0 : end subroutine calc_DryStaticEnergy
1907 : !================================================================
1908 :
1909 : end module nudging
|
