ROMS
Loading...
Searching...
No Matches
state_initialize.F
Go to the documentation of this file.
1#include "cppdefs.h"
3!
4!git $Id$
5!================================================== Hernan G. Arango ===
6! Copyright (c) 2002-2025 The ROMS Group !
7! Licensed under a MIT/X style license !
8! See License_ROMS.md !
9!=======================================================================
10! !
11! This routine initialize the state variable as follows: !
12! !
13! s_var(...,Lout) = fac !
14! !
15! where fac is scalar usually zero. !
16! !
17!=======================================================================
18!
19 implicit none
20!
21 PUBLIC :: state_initialize
22!
23 CONTAINS
24!
25!***********************************************************************
26 SUBROUTINE state_initialize (ng, tile, &
27 & LBi, UBi, LBj, UBj, LBij, UBij, &
28 & Lout, fac, &
29#ifdef MASKING
30 & rmask, umask, vmask, &
31#endif
32#ifdef ADJUST_BOUNDARY
33# ifdef SOLVE3D
34 & s_t_obc, &
35 & s_u_obc, s_v_obc, &
36# endif
37 & s_ubar_obc, s_vbar_obc, &
38 & s_zeta_obc, &
39#endif
40#ifdef ADJUST_WSTRESS
41 & s_sustr, s_svstr, &
42#endif
43#ifdef SOLVE3D
44# ifdef ADJUST_STFLUX
45 & s_tflux, &
46# endif
47 & s_t, s_u, s_v, &
48#else
49 & s_ubar, s_vbar, &
50#endif
51 & s_zeta)
52!***********************************************************************
53!
54 USE mod_param
55#if defined ADJUST_BOUNDARY || defined ADJUST_STFLUX || \
56 defined adjust_wstress
57 USE mod_ncparam
58 USE mod_scalars
59#endif
60!
61! Imported variable declarations.
62!
63 integer, intent(in) :: ng, tile
64 integer, intent(in) :: lbi, ubi, lbj, ubj, lbij, ubij
65 integer, intent(in) :: lout
66!
67 real(r8), intent(in) :: fac
68!
69#ifdef ASSUMED_SHAPE
70# ifdef MASKING
71 real(r8), intent(in) :: rmask(lbi:,lbj:)
72 real(r8), intent(in) :: umask(lbi:,lbj:)
73 real(r8), intent(in) :: vmask(lbi:,lbj:)
74# endif
75# ifdef ADJUST_BOUNDARY
76# ifdef SOLVE3D
77 real(r8), intent(inout) :: s_t_obc(lbij:,:,:,:,:,:)
78 real(r8), intent(inout) :: s_u_obc(lbij:,:,:,:,:)
79 real(r8), intent(inout) :: s_v_obc(lbij:,:,:,:,:)
80# endif
81 real(r8), intent(inout) :: s_ubar_obc(lbij:,:,:,:)
82 real(r8), intent(inout) :: s_vbar_obc(lbij:,:,:,:)
83 real(r8), intent(inout) :: s_zeta_obc(lbij:,:,:,:)
84# endif
85# ifdef ADJUST_WSTRESS
86 real(r8), intent(inout) :: s_sustr(lbi:,lbj:,:,:)
87 real(r8), intent(inout) :: s_svstr(lbi:,lbj:,:,:)
88# endif
89# ifdef SOLVE3D
90# ifdef ADJUST_STFLUX
91 real(r8), intent(inout) :: s_tflux(lbi:,lbj:,:,:,:)
92# endif
93 real(r8), intent(inout) :: s_t(lbi:,lbj:,:,:,:)
94 real(r8), intent(inout) :: s_u(lbi:,lbj:,:,:)
95 real(r8), intent(inout) :: s_v(lbi:,lbj:,:,:)
96# else
97 real(r8), intent(inout) :: s_ubar(lbi:,lbj:,:)
98 real(r8), intent(inout) :: s_vbar(lbi:,lbj:,:)
99# endif
100 real(r8), intent(inout) :: s_zeta(lbi:,lbj:,:)
101
102#else
103
104# ifdef MASKING
105 real(r8), intent(in) :: rmask(lbi:ubi,lbj:ubj)
106 real(r8), intent(in) :: umask(lbi:ubi,lbj:ubj)
107 real(r8), intent(in) :: vmask(lbi:ubi,lbj:ubj)
108# endif
109# ifdef ADJUST_BOUNDARY
110# ifdef SOLVE3D
111 real(r8), intent(inout) :: s_t_obc(lbij:ubij,n(ng),4, &
112 & Nbrec(ng),2,NT(ng))
113 real(r8), intent(inout) :: s_u_obc(lbij:ubij,n(ng),4,nbrec(ng),2)
114 real(r8), intent(inout) :: s_v_obc(lbij:ubij,n(ng),4,nbrec(ng),2)
115# endif
116 real(r8), intent(inout) :: s_ubar_obc(lbij:ubij,4,nbrec(ng),2)
117 real(r8), intent(inout) :: s_vbar_obc(lbij:ubij,4,nbrec(ng),2)
118 real(r8), intent(inout) :: s_zeta_obc(lbij:ubij,4,nbrec(ng),2)
119# endif
120# ifdef ADJUST_WSTRESS
121 real(r8), intent(inout) :: s_sustr(lbi:ubi,lbj:ubj,nfrec(ng),2)
122 real(r8), intent(inout) :: s_svstr(lbi:ubi,lbj:ubj,nfrec(ng),2)
123# endif
124# ifdef SOLVE3D
125# ifdef ADJUST_STFLUX
126 real(r8), intent(inout) :: s_tflux(lbi:ubi,lbj:ubj, &
127 & Nfrec(ng),2,NT(ng))
128# endif
129 real(r8), intent(inout) :: s_t(lbi:ubi,lbj:ubj,n(ng),3,nt(ng))
130 real(r8), intent(inout) :: s_u(lbi:ubi,lbj:ubj,n(ng),2)
131 real(r8), intent(inout) :: s_v(lbi:ubi,lbj:ubj,n(ng),2)
132# else
133 real(r8), intent(inout) :: s_ubar(lbi:ubi,lbj:ubj,:)
134 real(r8), intent(inout) :: s_vbar(lbi:ubi,lbj:ubj,:)
135# endif
136 real(r8), intent(inout) :: s_zeta(lbi:ubi,lbj:ubj,:)
137#endif
138!
139! Local variable declarations.
140!
141 integer :: i, j, k
142 integer :: ib, ir, it
143
144#include "set_bounds.h"
145!
146!-----------------------------------------------------------------------
147! Scale model state variable by a constant.
148!-----------------------------------------------------------------------
149!
150! Free-surface.
151!
152 DO j=jstrt,jendt
153 DO i=istrt,iendt
154 s_zeta(i,j,lout)=fac
155 END DO
156 END DO
157
158#ifdef ADJUST_BOUNDARY
159!
160! Free-surface open boundaries.
161!
162 IF (any(lobc(:,isfsur,ng))) THEN
163 DO ir=1,nbrec(ng)
164 IF ((lobc(iwest,isfsur,ng)).and. &
165 & domain(ng)%Western_Edge(tile)) THEN
166 DO j=jstr,jend
167 s_zeta_obc(j,iwest,ir,lout)=fac
168 END DO
169 END IF
170 IF ((lobc(ieast,isfsur,ng)).and. &
171 & domain(ng)%Eastern_Edge(tile)) THEN
172 DO j=jstr,jend
173 s_zeta_obc(j,ieast,ir,lout)=fac
174 END DO
175 END IF
176 IF ((lobc(isouth,isfsur,ng)).and. &
177 & domain(ng)%Southern_Edge(tile)) THEN
178 ib=isouth
179 DO i=istr,iend
180 s_zeta_obc(i,isouth,ir,lout)=fac
181 END DO
182 END IF
183 IF ((lobc(inorth,isfsur,ng)).and. &
184 & domain(ng)%Northern_Edge(tile)) THEN
185 DO i=istr,iend
186 s_zeta_obc(i,inorth,ir,lout)=fac
187 END DO
188 END IF
189 END DO
190 END IF
191#endif
192
193#ifndef SOLVE3D
194!
195! 2D U-momentum component.
196!
197 DO j=jstrt,jendt
198 DO i=istrp,iendt
199 s_ubar(i,j,lout)=fac
200 END DO
201 END DO
202#endif
203
204#ifdef ADJUST_BOUNDARY
205!
206! 2D U-momentum open boundaries.
207!
208 IF (any(lobc(:,isubar,ng))) THEN
209 DO ir=1,nbrec(ng)
210 IF ((lobc(iwest,isubar,ng)).and. &
211 & domain(ng)%Western_Edge(tile)) THEN
212 DO j=jstr,jend
213 s_ubar_obc(j,iwest,ir,lout)=fac
214 END DO
215 END IF
216 IF ((lobc(ieast,isubar,ng)).and. &
217 & domain(ng)%Eastern_Edge(tile)) THEN
218 DO j=jstr,jend
219 s_ubar_obc(j,ieast,ir,lout)=fac
220 END DO
221 END IF
222 IF ((lobc(isouth,isubar,ng)).and. &
223 & domain(ng)%Southern_Edge(tile)) THEN
224 DO i=istru,iend
225 s_ubar_obc(i,isouth,ir,lout)=fac
226 END DO
227 END IF
228 IF ((lobc(inorth,isubar,ng)).and. &
229 & domain(ng)%Northern_Edge(tile)) THEN
230 DO i=istru,iend
231 s_ubar_obc(i,inorth,ir,lout)=fac
232 END DO
233 END IF
234 END DO
235 END IF
236#endif
237
238#ifndef SOLVE3D
239!
240! 2D V-momentum component.
241!
242 DO j=jstrp,jendt
243 DO i=istrt,iendt
244 s_vbar(i,j,lout)=fac
245 END DO
246 END DO
247#endif
248
249#ifdef ADJUST_BOUNDARY
250!
251! 2D V-momentum open boundaries.
252!
253 IF (any(lobc(:,isvbar,ng))) THEN
254 DO ir=1,nbrec(ng)
255 IF ((lobc(iwest,isvbar,ng)).and. &
256 & domain(ng)%Western_Edge(tile)) THEN
257 DO j=jstrv,jend
258 s_vbar_obc(j,iwest,ir,lout)=fac
259 END DO
260 END IF
261 IF ((lobc(ieast,isvbar,ng)).and. &
262 & domain(ng)%Eastern_Edge(tile)) THEN
263 DO j=jstrv,jend
264 s_vbar_obc(j,ieast,ir,lout)=fac
265 END DO
266 END IF
267 IF ((lobc(isouth,isvbar,ng)).and. &
268 & domain(ng)%Southern_Edge(tile)) THEN
269 DO i=istr,iend
270 s_vbar_obc(i,isouth,ir,lout)=fac
271 END DO
272 END IF
273 IF ((lobc(inorth,isvbar,ng)).and. &
274 & domain(ng)%Northern_Edge(tile)) THEN
275 DO i=istr,iend
276 s_vbar_obc(i,inorth,ir,lout)=fac
277 END DO
278 END IF
279 END DO
280 END IF
281#endif
282
283#ifdef ADJUST_WSTRESS
284!
285! Surface momentum stress.
286!
287 DO ir=1,nfrec(ng)
288 DO j=jstrt,jendt
289 DO i=istrp,iendt
290 s_sustr(i,j,ir,lout)=fac
291 END DO
292 END DO
293 DO j=jstrp,jendt
294 DO i=istrt,iendt
295 s_svstr(i,j,ir,lout)=fac
296 END DO
297 END DO
298 END DO
299#endif
300
301#ifdef SOLVE3D
302!
303! 3D U-momentum component.
304!
305 DO k=1,n(ng)
306 DO j=jstrt,jendt
307 DO i=istrp,iendt
308 s_u(i,j,k,lout)=fac
309 END DO
310 END DO
311 END DO
312
313# ifdef ADJUST_BOUNDARY
314!
315! 3D U-momentum open boundaries.
316!
317 IF (any(lobc(:,isuvel,ng))) THEN
318 DO ir=1,nbrec(ng)
319 IF ((lobc(iwest,isuvel,ng)).and. &
320 & domain(ng)%Western_Edge(tile)) THEN
321 DO k=1,n(ng)
322 DO j=jstr,jend
323 s_u_obc(j,k,iwest,ir,lout)=fac
324 END DO
325 END DO
326 END IF
327 IF ((lobc(ieast,isuvel,ng)).and. &
328 & domain(ng)%Eastern_Edge(tile)) THEN
329 DO k=1,n(ng)
330 DO j=jstr,jend
331 s_u_obc(j,k,ieast,ir,lout)=fac
332 END DO
333 END DO
334 END IF
335 IF ((lobc(isouth,isuvel,ng)).and. &
336 & domain(ng)%Southern_Edge(tile)) THEN
337 DO k=1,n(ng)
338 DO i=istru,iend
339 s_u_obc(i,k,isouth,ir,lout)=fac
340 END DO
341 END DO
342 END IF
343 IF ((lobc(inorth,isuvel,ng)).and. &
344 & domain(ng)%Northern_Edge(tile)) THEN
345 DO k=1,n(ng)
346 DO i=istru,iend
347 s_u_obc(i,k,inorth,ir,lout)=fac
348 END DO
349 END DO
350 END IF
351 END DO
352 END IF
353# endif
354!
355! 3D V-momentum component.
356!
357 DO k=1,n(ng)
358 DO j=jstrp,jendt
359 DO i=istrt,iendt
360 s_v(i,j,k,lout)=fac
361 END DO
362 END DO
363 END DO
364
365# ifdef ADJUST_BOUNDARY
366!
367! 3D V-momentum open boundaries.
368!
369 IF (any(lobc(:,isvvel,ng))) THEN
370 DO ir=1,nbrec(ng)
371 IF ((lobc(iwest,isvvel,ng)).and. &
372 & domain(ng)%Western_Edge(tile)) THEN
373 DO k=1,n(ng)
374 DO j=jstrv,jend
375 s_v_obc(j,k,iwest,ir,lout)=fac
376 END DO
377 END DO
378 END IF
379 IF ((lobc(ieast,isvvel,ng)).and. &
380 & domain(ng)%Eastern_Edge(tile)) THEN
381 DO k=1,n(ng)
382 DO j=jstrv,jend
383 s_v_obc(j,k,ieast,ir,lout)=fac
384 END DO
385 END DO
386 END IF
387 IF ((lobc(isouth,isvvel,ng)).and. &
388 & domain(ng)%Southern_Edge(tile)) THEN
389 DO k=1,n(ng)
390 DO i=istr,iend
391 s_v_obc(i,k,isouth,ir,lout)=fac
392 END DO
393 END DO
394 END IF
395 IF ((lobc(inorth,isvvel,ng)).and. &
396 & domain(ng)%Northern_Edge(tile)) THEN
397 DO k=1,n(ng)
398 DO i=istr,iend
399 s_v_obc(i,k,inorth,ir,lout)=fac
400 END DO
401 END DO
402 END IF
403 END DO
404 END IF
405# endif
406!
407! Tracers.
408!
409 DO it=1,nt(ng)
410 DO k=1,n(ng)
411 DO j=jstrt,jendt
412 DO i=istrt,iendt
413 s_t(i,j,k,lout,it)=fac
414 END DO
415 END DO
416 END DO
417 END DO
418
419# ifdef ADJUST_BOUNDARY
420!
421! Tracers open boundaries.
422!
423 DO it=1,nt(ng)
424 IF (any(lobc(:,istvar(it),ng))) THEN
425 DO ir=1,nbrec(ng)
426 IF ((lobc(iwest,istvar(it),ng)).and. &
427 & domain(ng)%Western_Edge(tile)) THEN
428 DO k=1,n(ng)
429 DO j=jstr,jend
430 s_t_obc(j,k,iwest,ir,lout,it)=fac
431 END DO
432 END DO
433 END IF
434 IF ((lobc(ieast,istvar(it),ng)).and. &
435 & domain(ng)%Eastern_Edge(tile)) THEN
436 DO k=1,n(ng)
437 DO j=jstr,jend
438 s_t_obc(j,k,ieast,ir,lout,it)=fac
439 END DO
440 END DO
441 END IF
442 IF ((lobc(isouth,istvar(it),ng)).and. &
443 & domain(ng)%Southern_Edge(tile)) THEN
444 DO k=1,n(ng)
445 DO i=istr,iend
446 s_t_obc(i,k,isouth,ir,lout,it)=fac
447 END DO
448 END DO
449 END IF
450 IF ((lobc(inorth,istvar(it),ng)).and. &
451 & domain(ng)%Northern_Edge(tile)) THEN
452 DO k=1,n(ng)
453 DO i=istr,iend
454 s_t_obc(i,k,inorth,ir,lout,it)=fac
455 END DO
456 END DO
457 END IF
458 END DO
459 END IF
460 END DO
461# endif
462
463# ifdef ADJUST_STFLUX
464!
465! Surface tracers flux.
466!
467 DO it=1,nt(ng)
468 IF (lstflux(it,ng)) THEN
469 DO ir=1,nfrec(ng)
470 DO j=jstrt,jendt
471 DO i=istrt,iendt
472 s_tflux(i,j,ir,lout,it)=fac
473 END DO
474 END DO
475 END DO
476 END IF
477 END DO
478# endif
479
480#endif
481!
482 RETURN
483 END SUBROUTINE state_initialize
484
485 END MODULE state_initialize_mod
integer isvvel
integer isvbar
integer, dimension(:), allocatable istvar
integer isuvel
integer isfsur
integer isubar
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable nt
Definition mod_param.F:489
logical, dimension(:,:,:), allocatable lobc
integer, parameter iwest
logical, dimension(:,:), allocatable lstflux
integer, dimension(:), allocatable nfrec
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
integer, dimension(:), allocatable nbrec
subroutine, public state_initialize(ng, tile, lbi, ubi, lbj, ubj, lbij, ubij, lout, fac, rmask, umask, vmask, s_t_obc, s_u_obc, s_v_obc, s_ubar_obc, s_vbar_obc, s_zeta_obc, s_sustr, s_svstr, s_tflux, s_t, s_u, s_v, s_zeta)