New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
obc_oce.F90 in trunk/NEMO/OPA_SRC/OBC – NEMO

source: trunk/NEMO/OPA_SRC/OBC/obc_oce.F90 @ 699

Last change on this file since 699 was 699, checked in by smasson, 17 years ago

insert revision Id

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 18.6 KB
Line 
1MODULE obc_oce
2   !!==============================================================================
3   !!                       ***  MODULE obc_oce   ***
4   !! Open Boundary Cond. :   define related variables
5   !!==============================================================================
6   !!  OPA 9.0 , LOCEAN-IPSL (2005)
7   !! $Id$
8   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
9   !!----------------------------------------------------------------------
10#if defined key_obc
11   !!----------------------------------------------------------------------
12   !!   'key_obc' :                                 Open Boundary Condition
13   !!----------------------------------------------------------------------
14   !! history :
15   !!  8.0   01/91   (CLIPPER)  Original code
16   !!  8.5   06/02   (C. Talandier)  modules
17   !!        06/04   (F. Durand) ORCA2_ZIND config
18   !!        06/04   (F. Durand) Dimensions of arrays vsdta, tsdta, ssdta,
19   !!                vndta, tndta, sndta, uwdta, twdta, swdta, uedta, tedta, sedta
20   !!                are defined to the actual dimensions of the OBs i.e.
21   !!         (jpisd:jpisf,jpk,jptobc) for the South OB
22   !!            (jpind:jpinf,jpk,jptobc) for the North OB
23   !!         (jpjwd:jpjwf,jpk,jptobc) for the West OB
24   !!            (jpjed:jpjef,jpk,jptobc) for the East OB
25   !!                 
26   !!----------------------------------------------------------------------
27   !! * Modules used
28   USE par_oce         ! ocean parameters
29   USE obc_par         ! open boundary condition parameters
30
31   IMPLICIT NONE
32   PUBLIC
33
34   !!----------------------------------------------------------------------
35   !! open boundary variables
36   !!----------------------------------------------------------------------
37   !!
38   !!General variables for open boundaries:
39   !!--------------------------------------
40   INTEGER ::              & !: * namelist ??? *
41      nbobc    = 2  ,      & !: number of open boundaries ( 1=< nbobc =< 4 )
42      nobc_dta = 0           !:  = 0 use the initial state as obc data
43      !                      !   = 1 read obc data in obcxxx.dta files
44
45   LOGICAL ::  ln_obc_clim = .true.  !:  obc data files are climatological
46   LOGICAL ::  ln_obc_fla  = .false. !:  Flather open boundary condition not used
47   LOGICAL ::  ln_vol_cst  = .true.  !:  Conservation of the whole volume
48
49   REAL(wp) ::             & !!: open boundary namelist (namobc)
50      rdpein =  1.  ,      &  !: damping time scale for inflow at East open boundary
51      rdpwin =  1.  ,      &  !:    "                      "   at West open boundary
52      rdpsin =  1.  ,      &  !:    "                      "   at South open boundary
53      rdpnin =  1.  ,      &  !:    "                      "   at North open boundary
54      rdpeob = 15.  ,      &  !: damping time scale for the climatology at East open boundary
55      rdpwob = 15.  ,      &  !:    "                           "       at West open boundary
56      rdpsob = 15.  ,      &  !:    "                           "       at South open boundary
57      rdpnob = 15.  ,      &  !:    "                           "       at North open boundary
58      volemp =  1.            !: = 0 the total volume will have the variability of the
59                              !      surface Flux E-P else (volemp = 1) the volume will be constant
60                              !  = 1 the volume will be constant during all the integration.
61
62   LOGICAL ::              &  !:
63      lfbceast, lfbcwest,  &  !: logical flag for a fixed East and West open boundaries       
64      lfbcnorth, lfbcsouth    !: logical flag for a fixed North and South open boundaries       
65      !                       !  These logical flags are set to 'true' if damping time
66      !                       !  scale are set to 0 in the namelist, for both inflow and outflow).
67
68   REAL(wp), PUBLIC ::    &  !:
69      obcsurftot       !: Total lateral surface of open boundaries
70   
71   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: &  !:
72      obctmsk,            &  !: mask array identical to tmask, execpt along OBC where it is set to 0
73      !                      !  it used to calculate the cumulate flux E-P in the obcvol.F90 routine
74      obcumask, obcvmask     !: u-, v- Force filtering mask for the open
75      !                      !  boundary condition on grad D
76
77   !!----------------
78   !! Rigid lid case:
79   !!----------------
80   INTEGER ::   nbic !: number of isolated coastlines ( 0 <= nbic <= 3 )
81         
82   INTEGER, DIMENSION(jpnic,0:4,3) ::   &  !:
83      miic, mjic     !: position of isolated coastlines points
84
85   INTEGER, DIMENSION(0:4,3) ::   &  !:
86      mnic           !: number of points on isolated coastlines
87
88   REAL(wp), DIMENSION(jpi,jpj) ::   &  !:
89      gcbob          !: right hand side of the barotropic elliptic equation associated
90      !              !  with the OBC
91                                             
92   REAL(wp), DIMENSION(jpi,jpj,3) ::   &  !:
93      gcfobc         !: coef. associated with the contribution of isolated coastlines
94      !              !  to the right hand side of the barotropic elliptic equation
95
96   REAL(wp), DIMENSION(3) ::   &  !:
97      gcbic          !: time variation of the barotropic stream function along the
98      !              !  isolated coastlines
99
100   REAL(wp), DIMENSION(1) ::   &  !:
101      bsfic0         !: barotropic stream function on isolated coastline
102         
103   REAL(wp), DIMENSION(3) ::   &  !:
104      bsfic          !: barotropic stream function on isolated coastline
105         
106   !!--------------------
107   !! East open boundary:
108   !!--------------------
109   INTEGER ::   nie0  , nie1      !: do loop index in mpp case for jpieob
110   INTEGER ::   nie0p1, nie1p1    !: do loop index in mpp case for jpieob+1
111   INTEGER ::   nie0m1, nie1m1    !: do loop index in mpp case for jpieob-1
112   INTEGER ::   nje0  , nje1      !: do loop index in mpp case for jpjed, jpjef
113   INTEGER ::   nje0p1, nje1m1    !: do loop index in mpp case for jpjedp1,jpjefm1
114   INTEGER ::   nje1m2, nje0m1    !: do loop index in mpp case for jpjefm1-1,jpjed
115
116   REAL(wp), DIMENSION(jpj) ::    &  !:
117      bsfeob              !: now barotropic stream fuction computed at the OBC. The corres-
118      !                   ! ponding bsfn will be computed by the forward time step in dynspg.
119
120   REAL(wp), DIMENSION(jpj,3,3) ::   &  !:
121      bebnd               !: east boundary barotropic streamfunction over 3 rows
122      !                   ! and 3 time step (now, before, and before before)
123
124   REAL(wp), DIMENSION(jpjed:jpjef) ::   &  !:
125      bfoe,             & !: now climatology of the east boundary barotropic stream function
126      sshfoe,           & !: now climatology of the east boundary sea surface height
127      ubtfoe,vbtfoe       !: now climatology of the east boundary barotropic transport
128     
129   REAL(wp), DIMENSION(jpj,jpk) ::   &  !:
130      ufoe, vfoe,       & !: now climatology of the east boundary velocities
131      tfoe, sfoe,       & !: now climatology of the east boundary temperature and salinity
132      uclie               !: baroclinic componant of the zonal velocity after radiation
133      !                   ! in the obcdyn.F90 routine
134
135   REAL(wp), DIMENSION(jpjed:jpjef,jpj) ::   &  !:
136      sshfoe_b            !: east boundary ssh correction averaged over the barotropic loop
137                          !: (if Flather's algoritm applied at open boundary)
138
139   REAL(wp), DIMENSION(jpjed:jpjef,0:jptobc+1) ::   &  !:
140      sshedta, ubtedta    !: array used for interpolating monthly data on the east boundary
141
142   REAL(wp), DIMENSION(jpjed:jpjef,jpk,jptobc) ::   &  !:
143      uedta, tedta, sedta !: array used for interpolating monthly data on the east boundary
144
145   !!-------------------------------
146   !! Arrays for radiative East OBC:
147   !!-------------------------------
148   REAL(wp), DIMENSION(jpj,jpk,3,3) ::   &  !:
149      uebnd, vebnd                  !: baroclinic u & v component of the velocity over 3 rows
150                                    ! and 3 time step (now, before, and before before)
151
152   REAL(wp), DIMENSION(jpj,jpk,2,2) ::   &  !:
153      tebnd, sebnd                  !: East boundary temperature and salinity over 2 rows
154                                    ! and 2 time step (now and before)
155
156   REAL(wp), DIMENSION(jpj,jpk) ::   &  !:
157      u_cxebnd, v_cxebnd            !: Zonal component of the phase speed ratio computed with
158                                    ! radiation of u and v velocity (respectively) at the
159                                    ! east open boundary (u_cxebnd = cx rdt )
160
161   REAL(wp), DIMENSION(jpj,jpk) ::   &  !:
162      uemsk, vemsk, temsk           !: 2D mask for the East OB
163
164   ! Note that those arrays are optimized for mpp case
165   ! (hence the dimension jpj is the size of one processor subdomain)
166
167   !!--------------------
168   !! West open boundary
169   !!--------------------
170   INTEGER ::   niw0  , niw1       !: do loop index in mpp case for jpiwob
171   INTEGER ::   niw0p1, niw1p1     !: do loop index in mpp case for jpiwob+1
172   INTEGER ::   njw0  , njw1       !: do loop index in mpp case for jpjwd, jpjwf
173   INTEGER ::   njw0p1, njw1m1     !: do loop index in mpp case for jpjwdp1,jpjwfm1
174   INTEGER ::   njw1m2, njw0m1     !: do loop index in mpp case for jpjwfm2,jpjwd
175
176   REAL(wp), DIMENSION(jpj) ::   &  !:
177      bsfwob              !: now barotropic stream fuction computed at the OBC. The corres-
178      !                   !  ponding bsfn will be computed by the forward time step in dynspg.
179
180   REAL(wp), DIMENSION(jpj,3,3) ::   &  !:
181      bwbnd               !: West boundary barotropic streamfunction over
182      !                   !  3 rows and 3 time step (now, before, and before before)
183
184   REAL(wp), DIMENSION(jpjwd:jpjwf) ::   &  !:
185      bfow,             & !: now climatology of the west boundary barotropic stream function
186      sshfow,           & !: now climatology of the west boundary sea surface height
187      ubtfow,vbtfow       !: now climatology of the west boundary barotropic transport
188
189   REAL(wp), DIMENSION(jpj,jpk) ::   &  !:
190      ufow, vfow,       & !: now climatology of the west velocities
191      tfow, sfow,       & !: now climatology of the west temperature and salinity
192      ucliw               !: baroclinic componant of the zonal velocity after the radiation
193      !                   !  in the obcdyn.F90 routine
194
195   REAL(wp), DIMENSION(jpjwd:jpjwf,jpj) ::   &  !:
196      sshfow_b            !: west boundary ssh correction averaged over the barotropic loop
197                          !: (if Flather's algoritm applied at open boundary)
198
199   REAL(wp), DIMENSION(jpjwd:jpjwf,0:jptobc+1) ::   &  !:
200      sshwdta, ubtwdta    !: array used for interpolating monthly data on the west boundary
201
202   REAL(wp), DIMENSION(jpjwd:jpjwf,jpk,jptobc) ::   &  !:
203      uwdta, twdta, swdta !: array used for interpolating monthly data on the west boundary
204
205   !!-------------------------------
206   !! Arrays for radiative West OBC
207   !!-------------------------------
208   REAL(wp), DIMENSION(jpj,jpk,3,3) ::   &  !:
209      uwbnd, vwbnd                  !: baroclinic u & v components of the velocity over 3 rows
210      !                             !  and 3 time step (now, before, and before before)
211
212   REAL(wp), DIMENSION(jpj,jpk,2,2) ::   &  !:
213      twbnd, swbnd                  !: west boundary temperature and salinity over 2 rows and
214      !                             !  2 time step (now and before)
215
216   REAL(wp), DIMENSION(jpj,jpk) ::    &  !:
217      u_cxwbnd, v_cxwbnd            !: Zonal component of the phase speed ratio computed with
218      !                             !  radiation of zonal and meridional velocity (respectively)
219      !                             !  at the west open boundary (u_cxwbnd = cx rdt )
220
221   REAL(wp), DIMENSION(jpj,jpk) ::    &  !:
222      uwmsk, vwmsk, twmsk           !: 2D mask for the West OB
223
224   ! Note that those arrays are optimized for mpp case
225   ! (hence the dimension jpj is the size of one processor subdomain)
226
227   !!---------------------
228   !! North open boundary
229   !!---------------------
230   INTEGER ::   nin0  , nin1       !: do loop index in mpp case for jpind, jpinf
231   INTEGER ::   nin0p1, nin1m1     !: do loop index in mpp case for jpindp1, jpinfm1
232   INTEGER ::   nin1m2, nin0m1     !: do loop index in mpp case for jpinfm1-1,jpind
233   INTEGER ::   njn0  , njn1       !: do loop index in mpp case for jpnob
234   INTEGER ::   njn0p1, njn1p1     !: do loop index in mpp case for jpnob+1
235   INTEGER ::   njn0m1, njn1m1     !: do loop index in mpp case for jpnob-1
236
237   REAL(wp), DIMENSION(jpi) ::   &  !:
238      bsfnob              !: now barotropic stream fuction computed at the OBC. The corres-
239      !                   !  ponding bsfn will be computed by the forward time step in dynspg.
240
241   REAL(wp), DIMENSION(jpi,3,3) ::   &  !:
242      bnbnd               !: north boundary barotropic streamfunction over
243      !                   !  3 rows and 3 time step (now, before, and before before)
244
245   REAL(wp), DIMENSION(jpind:jpinf) ::   &  !:
246      bfon,             & !: now climatology of the north boundary barotropic stream function
247      sshfon,           & !: now climatology of the north boundary sea surface height
248      ubtfon,vbtfon       !: now climatology of the north boundary barotropic transport
249
250   REAL(wp), DIMENSION(jpi,jpk) ::   &    !:
251      ufon, vfon,       & !: now climatology of the north boundary velocities
252      tfon, sfon,       & !: now climatology of the north boundary temperature and salinity
253      vclin               !: baroclinic componant of the meridian velocity after the radiation
254      !                   !  in yhe obcdyn.F90 routine
255
256   REAL(wp), DIMENSION(jpind:jpinf,jpj) ::   &  !:
257      sshfon_b            !: north boundary ssh correction averaged over the barotropic loop
258                          !: (if Flather's algoritm applied at open boundary)
259
260   REAL(wp), DIMENSION(jpind:jpinf,0:jptobc+1) ::   &  !:
261      sshndta, vbtndta   !: array used for interpolating monthly data on the north boundary
262
263   REAL(wp), DIMENSION(jpind:jpinf,jpk,jptobc) ::   &  !:
264      vndta, tndta, sndta !: array used for interpolating monthly data on the north boundary
265
266   !!--------------------------------
267   !! Arrays for radiative North OBC
268   !!--------------------------------
269   !!   
270   REAL(wp), DIMENSION(jpi,jpk,3,3) ::   &   !:
271      unbnd, vnbnd                  !: baroclinic u & v components of the velocity over 3
272      !                             !  rows and 3 time step (now, before, and before before)
273
274   REAL(wp), DIMENSION(jpi,jpk,2,2) ::   &   !:
275      tnbnd, snbnd                  !: north boundary temperature and salinity over
276      !                             !  2 rows and 2 time step (now and before)
277
278   REAL(wp), DIMENSION(jpi,jpk) ::   &     !:
279      u_cynbnd, v_cynbnd            !: Meridional component of the phase speed ratio compu-
280      !                             !  ted with radiation of zonal and meridional velocity
281      !                             !  (respectively) at the north OB (u_cynbnd = cx rdt )
282
283   REAL(wp), DIMENSION(jpi,jpk) ::   &  !:
284      unmsk, vnmsk, tnmsk           !: 2D mask for the North OB
285
286   ! Note that those arrays are optimized for mpp case
287   ! (hence the dimension jpj is the size of one processor subdomain)
288   
289   !!---------------------
290   !! South open boundary
291   !!---------------------
292   INTEGER ::   nis0  , nis1       !: do loop index in mpp case for jpisd, jpisf
293   INTEGER ::   nis0p1, nis1m1     !: do loop index in mpp case for jpisdp1, jpisfm1
294   INTEGER ::   nis1m2, nis0m1     !: do loop index in mpp case for jpisfm1-1,jpisd
295   INTEGER ::   njs0  , njs1       !: do loop index in mpp case for jpsob
296   INTEGER ::   njs0p1, njs1p1     !: do loop index in mpp case for jpsob+1
297
298   REAL(wp), DIMENSION(jpi) ::    &   !:
299      bsfsob              !: now barotropic stream fuction computed at the OBC.The corres-
300      !                   !  ponding bsfn will be computed by the forward time step in dynspg.
301   REAL(wp), DIMENSION(jpi,3,3) ::   &   !:
302      bsbnd               !: south boundary barotropic stream function over
303      !                   !  3 rows and 3 time step (now, before, and before before)
304
305   REAL(wp), DIMENSION(jpisd:jpisf) ::    &   !:
306      bfos,             & !: now climatology of the south boundary barotropic stream function
307      sshfos,           & !: now climatology of the south boundary sea surface height
308      ubtfos,vbtfos       !: now climatology of the south boundary barotropic transport
309
310   REAL(wp), DIMENSION(jpi,jpk) ::    &   !:
311      ufos, vfos,       & !: now climatology of the south boundary velocities
312      tfos, sfos,       & !: now climatology of the south boundary temperature and salinity
313      vclis               !: baroclinic componant of the meridian velocity after the radiation
314      !                   !  in the obcdyn.F90 routine
315
316   REAL(wp), DIMENSION(jpisd:jpisf,jpj) ::   &  !:
317      sshfos_b            !: south boundary ssh correction averaged over the barotropic loop
318                          !: (if Flather's algoritm applied at open boundary)
319
320   REAL(wp), DIMENSION(jpisd:jpisf,0:jptobc+1) ::    &    !:
321      sshsdta, vbtsdta    !: array used for interpolating monthly data on the south boundary
322
323   REAL(wp), DIMENSION(jpisd:jpisf,jpk,jptobc) ::    &    !:
324      vsdta, tsdta, ssdta   !: array used for interpolating monthly data on the south boundary
325
326   !!--------------------------------
327   !! Arrays for radiative South OBC
328   !!--------------------------------
329   !!                        computed by the forward time step in dynspg.
330   REAL(wp), DIMENSION(jpi,jpk,3,3) ::   &   !:
331      usbnd, vsbnd                  !: baroclinic u & v components of the velocity over 3
332      !                             !  rows and 3 time step (now, before, and before before)
333
334   REAL(wp), DIMENSION(jpi,jpk,2,2) ::   &  !:
335      tsbnd, ssbnd                  !: south boundary temperature and salinity over
336      !                             !  2 rows and 2 time step (now and before)
337
338   REAL(wp), DIMENSION(jpi,jpk) ::   &  !:
339      u_cysbnd, v_cysbnd            !: Meridional component of the phase speed ratio compu-
340      !                             !  ted with radiation of zonal and meridional velocity
341      !                             !  (repsectively) at the south OB (u_cynbnd = cx rdt )
342
343   REAL(wp), DIMENSION(jpi,jpk) ::   &  !:
344      usmsk, vsmsk, tsmsk           !: 2D mask for the South OB
345
346   ! Note that those arrays are optimized for mpp case
347   ! (hence the dimension jpj is the size of one processor subdomain)
348
349#else
350   !!----------------------------------------------------------------------
351   !!   Default option :                                       Empty module
352   !!----------------------------------------------------------------------
353#endif
354
355   !!======================================================================
356END MODULE obc_oce
Note: See TracBrowser for help on using the repository browser.