ROMS
Loading...
Searching...
No Matches
ana_scope.h
Go to the documentation of this file.
1!!
2 SUBROUTINE ana_scope (ng, tile, model)
3!
4!! git $Id$
5!!======================================================================
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 subroutine sets analytical adjoint sensitivity spatial scope !
12! masking arrays. !
13! !
14!=======================================================================
15!
16 USE mod_param
17 USE mod_grid
18 USE mod_ncparam
19!
20! Imported variable declarations.
21!
22 integer, intent(in) :: ng, tile, model
23!
24! Local variable declarations.
25!
26 character (len=*), parameter :: MyFile = &
27 & __FILE__
28!
29#include "tile.h"
30!
31 CALL ana_scope_tile (ng, tile, model, &
32 & lbi, ubi, lbj, ubj, &
33 & imins, imaxs, jmins, jmaxs, &
34#ifdef MASKING
35 & grid(ng) % rmask, &
36 & grid(ng) % umask, &
37 & grid(ng) % vmask, &
38#endif
39 & grid(ng) % Rscope, &
40 & grid(ng) % Uscope, &
41 & grid(ng) % Vscope)
42!
43! Set analytical header file name used.
44!
45#ifdef DISTRIBUTE
46 IF (lanafile) THEN
47#else
48 IF (lanafile.and.(tile.eq.0)) THEN
49#endif
50 ananame(22)=myfile
51 END IF
52!
53 RETURN
54 END SUBROUTINE ana_scope
55!
56!***********************************************************************
57 SUBROUTINE ana_scope_tile (ng, tile, model, &
58 & LBi, UBi, LBj, UBj, &
59 & IminS, ImaxS, JminS, JmaxS, &
60#ifdef MASKING
61 & rmask, umask, vmask, &
62#endif
63 & Rscope, Uscope, Vscope)
64!***********************************************************************
65!
66 USE mod_param
67 USE mod_scalars
68!
70#ifdef DISTRIBUTE
72#endif
73!
74! Imported variable declarations.
75!
76 integer, intent(in) :: ng, tile, model
77 integer, intent(in) :: LBi, UBi, LBj, UBj
78 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
79!
80#ifdef ASSUMED_SHAPE
81# ifdef MASKING
82 real(r8), intent(in) :: rmask(LBi:,LBj:)
83 real(r8), intent(in) :: umask(LBi:,LBj:)
84 real(r8), intent(in) :: vmask(LBi:,LBj:)
85# endif
86 real(r8), intent(out) :: Rscope(LBi:,LBj:)
87 real(r8), intent(out) :: Uscope(LBi:,LBj:)
88 real(r8), intent(out) :: Vscope(LBi:,LBj:)
89#else
90# ifdef MASKING
91 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
92 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
93 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
94# endif
95 real(r8), intent(out) :: Rscope(LBi:UBi,LBj:UBj)
96 real(r8), intent(out) :: Uscope(LBi:UBi,LBj:UBj)
97 real(r8), intent(out) :: Vscope(LBi:UBi,LBj:UBj)
98#endif
99!
100! Local variable declarations.
101!
102 integer :: Imin, Imax, Jmin, Jmax, i, j
103!
104 real(r8) :: scope(IminS:ImaxS,JminS:JmaxS)
105
106#include "set_bounds.h"
107!
108!-----------------------------------------------------------------------
109! Set Land/Sea mask of RHO-points: Land=0, Sea=1.
110!-----------------------------------------------------------------------
111!
112! Notice that private scratch array "mask" is used to allow
113! computation within a parallel loop.
114!
115#ifdef DOUBLE_GYRE
116 imin=-5+(lm(ng)+1)/2
117 imax=imin+10
118 jmin=-5+(mm(ng)+1)/2
119 jmax=jmin+10
120
121 DO j=jstrm2,jendp2
122 DO i=istrm2,iendp2
123 scope(i,j)=0.0_r8
124 IF (((imin.le.i).and.(i.le.imax)).and. &
125 & ((jmin.le.j).and.(j.le.jmax))) THEN
126 scope(i,j)=1.0_r8
127 END IF
128 END DO
129 END DO
130#else
131 ana_scope.h: no values provided for spatial scope masking.
132#endif
133!
134 DO j=jstrt,jendt
135 DO i=istrt,iendt
136 rscope(i,j)=scope(i,j)
137#ifdef MASKING
138 rscope(i,j)=rscope(i,j)*rmask(i,j)
139#endif
140 END DO
141 END DO
142!
143!-----------------------------------------------------------------------
144! Compute Land/Sea mask of U- and V-points.
145!-----------------------------------------------------------------------
146!
147 DO j=jstrt,jendt
148 DO i=istrp,iendt
149 uscope(i,j)=scope(i-1,j)*scope(i,j)
150#ifdef MASKING
151 uscope(i,j)=uscope(i,j)*umask(i,j)
152#endif
153 END DO
154 END DO
155 DO j=jstrp,jendt
156 DO i=istrt,iendt
157 vscope(i,j)=scope(i,j-1)*scope(i,j)
158#ifdef MASKING
159 vscope(i,j)=vscope(i,j)*vmask(i,j)
160#endif
161 END DO
162 END DO
163!
164!-----------------------------------------------------------------------
165! Exchange boundary edges.
166!-----------------------------------------------------------------------
167!
168 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
169 CALL exchange_r2d_tile (ng, tile, &
170 & lbi, ubi, lbj, ubj, &
171 & rscope)
172 CALL exchange_u2d_tile (ng, tile, &
173 & lbi, ubi, lbj, ubj, &
174 & uscope)
175 CALL exchange_v2d_tile (ng, tile, &
176 & lbi, ubi, lbj, ubj, &
177 & vscope)
178 END IF
179
180#ifdef DISTRIBUTE
181 CALL mp_exchange2d (ng, tile, model, 3, &
182 & lbi, ubi, lbj, ubj, &
183 & nghostpoints, &
184 & ewperiodic(ng), nsperiodic(ng), &
185 & rscope, uscope, vscope)
186#endif
187!
188 RETURN
189 END SUBROUTINE ana_scope_tile
subroutine ana_scope_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, rmask, umask, vmask, rscope, uscope, vscope)
Definition ana_scope.h:64
subroutine ana_scope(ng, tile, model)
Definition ana_scope.h:3
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
logical lanafile
character(len=256), dimension(39) ananame
integer nghostpoints
Definition mod_param.F:710
integer, dimension(:), allocatable lm
Definition mod_param.F:455
integer, dimension(:), allocatable mm
Definition mod_param.F:456
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)