ROMS
Loading...
Searching...
No Matches
ice_set_avg.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined ICE_MODEL && defined SOLVE3D
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This module accumulates and computes output time-averaged ice !
13! model fields. Due to synchronization, the time-averaged fields !
14! are computed in delayed mode. All averages are accumulated at !
15! the beggining of the next time-step. !
16! !
17!=======================================================================
18!
19 USE mod_param
20 USE mod_grid
21 USE mod_ice
22 USE mod_ncparam
23 USE mod_scalars
24 USE mod_stepping
25!
27# ifdef DISTRIBUTE
29# endif
30!
31 implicit none
32!
33 PRIVATE
34 PUBLIC :: ice_set_avg
35 PUBLIC :: ice_set_avg_tile
36!
37 CONTAINS
38!
39!***********************************************************************
40 SUBROUTINE ice_set_avg (ng, tile, model)
41!***********************************************************************
42!
43! Imported variable declarations.
44!
45 integer, intent(in) :: ng, tile, model
46!
47! Local variable declarations.
48!
49 character (len=*), parameter :: MyFile = &
50 & __FILE__
51!
52# include "tile.h"
53!
54# ifdef PROFILE
55 CALL wclock_on (ng, model, 5, __line__, myfile)
56# endif
57 CALL ice_set_avg_tile (ng, tile, model, &
58 & lbi, ubi, lbj, ubj, &
59 & imins, imaxs, jmins, jmaxs, &
60 & iout)
61
62# ifdef PROFILE
63 CALL wclock_off (ng, model, 5, __line__, myfile)
64# endif
65!
66 RETURN
67 END SUBROUTINE ice_set_avg
68!
69!***********************************************************************
70 SUBROUTINE ice_set_avg_tile (ng, tile, model, &
71 & LBi, UBi, LBj, UBj, &
72 & IminS, ImaxS, JminS, JmaxS, &
73 & Iout)
74!***********************************************************************
75!
76!
77! Imported variable declarations.
78!
79 integer, intent(in) :: ng, tile, model
80 integer, intent(in) :: LBi, UBi, LBj, UBj
81 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
82 integer, intent(in) :: Iout
83!
84! Local variable declarations.
85!
86 integer :: i, ifield, j, nf, ns
87!
88 real(r8) :: fac
89!
90 real(r8) :: rfac(IminS:ImaxS,JminS:JmaxS)
91 real(r8) :: ufac(IminS:ImaxS,JminS:JmaxS)
92 real(r8) :: vfac(IminS:ImaxS,JminS:JmaxS)
93
94# include "set_bounds.h"
95!
96!-----------------------------------------------------------------------
97! Return if time-averaging window is zero.
98!-----------------------------------------------------------------------
99!
100 IF (navg(ng).eq.0) RETURN
101!
102!-----------------------------------------------------------------------
103! Initialize time-averaged arrays when appropriate. Notice that
104! fields are initilized twice during re-start. However, the time-
105! averaged fields are computed correctly.
106!-----------------------------------------------------------------------
107!
108 IF (((iic(ng).gt.ntsavg(ng)).and. &
109 & (mod(iic(ng)-1,navg(ng)).eq.1)).or. &
110 & ((iic(ng).ge.ntsavg(ng)).and.(navg(ng).eq.1)).or. &
111 & ((nrrec(ng).gt.0).and.(iic(ng).eq.ntstart(ng)))) THEN
112!
113! Initialize ice model state variables.
114!
115 DO ns=1,nices
116 IF (isice(ns).gt.0) THEN
117 ifield=isice(ns)
118 IF (aout(ifield,ng)) THEN
119 DO j=jstrr,jendr
120 DO i=istrr,iendr
121 ice_savg(ns,ng)%var(i,j)=ice(ng)%Si(i,j,iout,ns)
122# ifdef MASKING
123 ice_savg(ns,ng)%var(i,j)=ice_savg(ns,ng)%var(i,j)* &
124 & grid(ng)%rmask_full(i,j)
125# endif
126 END DO
127 END DO
128 END IF
129 END IF
130 END DO
131!
132! Initialize ice model internal variables.
133!
134 DO nf=1,nicef
135 IF (ifice(nf).gt.0) THEN
136 ifield=ifice(nf)
137 IF (aout(ifield,ng)) THEN
138 DO j=jstrr,jendr
139 DO i=istrr,iendr
140 ice_favg(nf,ng)%var(i,j)=ice(ng)%Fi(i,j,nf)
141# ifdef MASKING
142 ice_favg(nf,ng)%var(i,j)=ice_favg(nf,ng)%var(i,j)* &
143 & grid(ng)%rmask_full(i,j)
144# endif
145 END DO
146 END DO
147 END IF
148 END IF
149 END DO
150!
151!-----------------------------------------------------------------------
152! Accumulate time-averaged fields.
153!-----------------------------------------------------------------------
154!
155 ELSE IF (iic(ng).gt.ntsavg(ng)) THEN
156!
157! Accumulate ice model state variables.
158!
159 DO ns=1,nices
160 IF (isice(ns).gt.0) THEN
161 ifield=isice(ns)
162 IF (aout(ifield,ng)) THEN
163 DO j=jstrr,jendr
164 DO i=istrr,iendr
165 ice_savg(ns,ng)%var(i,j)=ice_savg(ns,ng)%var(i,j)+ &
166# ifdef MASKING
167 & grid(ng)%rmask_full(i,j)* &
168# endif
169 & ice(ng)%Si(i,j,iout,ns)
170 END DO
171 END DO
172 END IF
173 END IF
174 END DO
175!
176! Accumulate ice model internal variables.
177!
178 DO nf=1,nicef
179 IF (ifice(nf).gt.0) THEN
180 ifield=ifice(nf)
181 IF (aout(ifield,ng)) THEN
182 DO j=jstrr,jendr
183 DO i=istrr,iendr
184 ice_favg(nf,ng)%var(i,j)=ice_favg(nf,ng)%var(i,j)+ &
185# ifdef MASKING
186 & grid(ng)%rmask_full(i,j)* &
187# endif
188 & ice(ng)%Fi(i,j,nf)
189 END DO
190 END DO
191 END IF
192 END IF
193 END DO
194 END IF
195!
196!-----------------------------------------------------------------------
197! Convert accumulated sums into time-averages, if appropriate.
198! Notice that we need to apply periodic conditions, if any, since
199! the full I- and J-ranges are different.
200!-----------------------------------------------------------------------
201!
202 IF (((iic(ng).gt.ntsavg(ng)).and. &
203 & (mod(iic(ng)-1,navg(ng)).eq.0).and. &
204 & ((iic(ng).ne.ntstart(ng)).or.(nrrec(ng).eq.0))).or. &
205 & ((iic(ng).ge.ntsavg(ng)).and.(navg(ng).eq.1))) THEN
206!
207! Set time-averaged factors for each C-grid variable type. Notice that
208! the I- and J-ranges are all grid types are the same for convinience.
209# ifdef WET_DRY
210! In wetting and drying, the sums are devided by the number of times
211! that each qrid point is wet.
212# endif
213!
214# ifdef WET_DRY
215 DO j=jstrr,jendr
216 DO i=istrr,iendr
217 rfac(i,j)=1.0_r8/max(1.0_r8, grid(ng)%rmask_avg(i,j))
218 ufac(i,j)=1.0_r8/max(1.0_r8, grid(ng)%umask_avg(i,j))
219 vfac(i,j)=1.0_r8/max(1.0_r8, grid(ng)%vmask_avg(i,j))
220 END DO
221 END DO
222# else
223 fac=1.0_r8/real(navg(ng),r8)
224 DO j=jstrr,jendr
225 DO i=istrr,iendr
226 rfac(i,j)=fac
227 ufac(i,j)=fac
228 vfac(i,j)=fac
229 END DO
230 END DO
231# endif
232!
233! Process ice model state variables.
234!
235 DO ns=1,nices
236 IF (isice(ns).gt.0) THEN
237 ifield=isice(ns)
238 IF (aout(ifield,ng)) THEN
239 DO j=jstrr,jendr
240 DO i=istrr,iendr
241 ice_savg(ns,ng)%var(i,j)=rfac(i,j)* &
242 & ice_savg(ns,ng)%var(i,j)
243 END DO
244 END DO
245 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
246 CALL exchange_r2d_tile (ng, tile, &
247 & lbi, ubi, lbj, ubj, &
248 & ice_savg(ns,ng)%var)
249
250# ifdef DISTRIBUTE
251 CALL mp_exchange2d (ng, tile, model, 1, &
252 & lbi, ubi, lbj, ubj, &
253 & nghostpoints, &
254 & ewperiodic(ng), nsperiodic(ng), &
255 & ice_savg(ns,ng)%var)
256# endif
257 END IF
258 END IF
259 END IF
260 END DO
261!
262! Process ice model internal variables.
263!
264 DO nf=1,nicef
265 IF (ifice(nf).gt.0) THEN
266 ifield=ifice(nf)
267 IF (aout(ifield,ng)) THEN
268 DO j=jstrr,jendr
269 DO i=istrr,iendr
270 ice_favg(nf,ng)%var(i,j)=rfac(i,j)* &
271 & ice_favg(nf,ng)%var(i,j)
272 END DO
273 END DO
274 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
275 CALL exchange_r2d_tile (ng, tile, &
276 & lbi, ubi, lbj, ubj, &
277 & ice_favg(ns,ng)%var)
278
279# ifdef DISTRIBUTE
280 CALL mp_exchange2d (ng, tile, model, 1, &
281 & lbi, ubi, lbj, ubj, &
282 & nghostpoints, &
283 & ewperiodic(ng), nsperiodic(ng), &
284 & ice_favg(ns,ng)%var)
285# endif
286 END IF
287 END IF
288 END IF
289 END DO
290 END IF
291!
292 RETURN
293 END SUBROUTINE ice_set_avg_tile
294#endif
295 END MODULE ice_set_avg_mod
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer, dimension(nicef) ifice
Definition ice_mod.h:169
type(t_ice), dimension(:), allocatable ice
Definition ice_mod.h:283
integer, parameter nicef
Definition ice_mod.h:167
integer, dimension(nices) isice
Definition ice_mod.h:135
type(t_ice_avg), dimension(:,:), allocatable ice_savg
Definition ice_mod.h:320
type(t_ice_avg), dimension(:,:), allocatable ice_favg
Definition ice_mod.h:319
integer, parameter nices
Definition ice_mod.h:130
logical, dimension(:,:), allocatable aout
integer nghostpoints
Definition mod_param.F:710
integer, dimension(:), allocatable nrrec
integer, dimension(:), allocatable iic
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
integer, dimension(:), allocatable navg
integer, dimension(:), allocatable ntstart
integer, dimension(:), allocatable ntsavg
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3