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 branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/OBC – NEMO

source: branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/OBC/obc_oce.F90 @ 2200

Last change on this file since 2200 was 2200, checked in by smasson, 13 years ago

merge trunk rev 2198 into DEV_r2106_LOCEAN2010

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 15.0 KB
Line 
1MODULE obc_oce
2   !!==============================================================================
3   !!                       ***  MODULE obc_oce   ***
4   !! Open Boundary Cond. :   define related variables
5   !!==============================================================================
6   !!----------------------------------------------------------------------
7   !!   'key_obc' :                                 Open Boundary Condition
8   !!----------------------------------------------------------------------
9   !! history :
10   !!  8.0   01/91   (CLIPPER)  Original code
11   !!  8.5   06/02   (C. Talandier)  modules
12   !!        06/04   (F. Durand) ORCA2_ZIND config
13   !!                 
14   !!----------------------------------------------------------------------
15   !! * Modules used
16   USE par_oce         ! ocean parameters
17   USE obc_par         ! open boundary condition parameters
18
19#if defined key_obc
20
21   IMPLICIT NONE
22   PUBLIC
23
24   !!----------------------------------------------------------------------
25   !! open boundary variables
26   !!----------------------------------------------------------------------
27   !
28   !                                  !!* Namelist namobc: open boundary condition *
29   INTEGER           ::   nn_obcdta   = 0        !:  = 0 use the initial state as obc data
30   !                                             !   = 1 read obc data in obcxxx.dta files
31   CHARACTER(len=20) ::   cn_obcdta   = 'annual' !: set to annual  if obc datafile hold 1 year of data
32   !                                             !  set to monthly if obc datafile hold 1 month of data
33   LOGICAL           ::   ln_obc_clim = .true.   !:  obc data files are climatological
34   LOGICAL           ::   ln_obc_fla  = .false.  !:  Flather open boundary condition not used
35   LOGICAL           ::   ln_vol_cst  = .true.   !:  Conservation of the whole volume
36   REAL(wp)          ::   rn_dpein    =  1.      !: damping time scale for inflow at East open boundary
37   REAL(wp)          ::   rn_dpwin    =  1.      !:    "                      "   at West open boundary
38   REAL(wp)          ::   rn_dpsin    =  1.      !:    "                      "   at South open boundary
39   REAL(wp)          ::   rn_dpnin    =  1.      !:    "                      "   at North open boundary
40   REAL(wp)          ::   rn_dpeob    = 15.      !: damping time scale for the climatology at East open boundary
41   REAL(wp)          ::   rn_dpwob    = 15.      !:    "                           "       at West open boundary
42   REAL(wp)          ::   rn_dpsob    = 15.      !:    "                           "       at South open boundary
43   REAL(wp)          ::   rn_dpnob    = 15.      !:    "                           "       at North open boundary
44   REAL(wp)          ::   rn_volemp   =  1.      !: = 0 the total volume will have the variability of the
45   !                                             !      surface Flux E-P else (volemp = 1) the volume will be constant
46   !                                             !  = 1 the volume will be constant during all the integration.
47
48   !                                  !!! OLD non-DOCTOR name of namelist variables
49   INTEGER  ::   nbobc                 !: number of open boundaries ( 1=< nbobc =< 4 )
50   INTEGER  ::   nobc_dta              !:  = 0 use the initial state as obc data
51   REAL(wp) ::   rdpein                !: damping time scale for inflow at East open boundary
52   REAL(wp) ::   rdpwin                !:    "                      "   at West open boundary
53   REAL(wp) ::   rdpsin                !:    "                      "   at South open boundary
54   REAL(wp) ::   rdpnin                !:    "                      "   at North open boundary
55   REAL(wp) ::   rdpeob                !: damping time scale for the climatology at East open boundary
56   REAL(wp) ::   rdpwob                !:    "                           "       at West open boundary
57   REAL(wp) ::   rdpsob                !:    "                           "       at South open boundary
58   REAL(wp) ::   rdpnob                !:    "                           "       at North open boundary
59   REAL(wp) ::   volemp                !: = 0 the total volume will have the variability of the
60   CHARACTER(len=20) :: cffile
61
62
63   !!General variables for open boundaries:
64   !!--------------------------------------
65   LOGICAL ::              &  !:
66      lfbceast, lfbcwest,  &  !: logical flag for a fixed East and West open boundaries       
67      lfbcnorth, lfbcsouth    !: logical flag for a fixed North and South open boundaries       
68      !                       !  These logical flags are set to 'true' if damping time
69      !                       !  scale are set to 0 in the namelist, for both inflow and outflow).
70
71   REAL(wp), PUBLIC ::   obcsurftot       !: Total lateral surface of open boundaries
72   
73   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: &  !:
74      obctmsk,            &  !: mask array identical to tmask, execpt along OBC where it is set to 0
75      !                      !  it used to calculate the cumulate flux E-P in the obcvol.F90 routine
76      obcumask, obcvmask     !: u-, v- Force filtering mask for the open
77      !                      !  boundary condition on grad D
78
79   !!--------------------
80   !! East open boundary:
81   !!--------------------
82   INTEGER ::   nie0  , nie1      !: do loop index in mpp case for jpieob
83   INTEGER ::   nie0p1, nie1p1    !: do loop index in mpp case for jpieob+1
84   INTEGER ::   nie0m1, nie1m1    !: do loop index in mpp case for jpieob-1
85   INTEGER ::   nje0  , nje1      !: do loop index in mpp case for jpjed, jpjef
86   INTEGER ::   nje0p1, nje1m1    !: do loop index in mpp case for jpjedp1,jpjefm1
87   INTEGER ::   nje1m2, nje0m1    !: do loop index in mpp case for jpjefm1-1,jpjed
88
89   REAL(wp), DIMENSION(jpj) ::   &  !:
90      sshfoe,           & !: now climatology of the east boundary sea surface height
91      ubtfoe,vbtfoe       !: now climatology of the east boundary barotropic transport
92     
93   REAL(wp), DIMENSION(jpj,jpk) ::   &  !:
94      ufoe, vfoe,       & !: now climatology of the east boundary velocities
95      tfoe, sfoe,       & !: now climatology of the east boundary temperature and salinity
96      uclie               !: baroclinic componant of the zonal velocity after radiation
97      !                   ! in the obcdyn.F90 routine
98
99   REAL(wp), DIMENSION(jpi,jpj) ::   sshfoe_b      !: east boundary ssh correction averaged over the barotropic loop
100      !                                            !  (if Flather's algoritm applied at open boundary)
101
102   !!-------------------------------
103   !! Arrays for radiative East OBC:
104   !!-------------------------------
105   REAL(wp), DIMENSION(jpj,jpk,3,3) ::   uebnd, vebnd      !: baroclinic u & v component of the velocity over 3 rows
106      !                                                    !  and 3 time step (now, before, and before before)
107   REAL(wp), DIMENSION(jpj,jpk,2,2) ::   tebnd, sebnd      !: East boundary temperature and salinity over 2 rows
108      !                                                    !  and 2 time step (now and before)
109   REAL(wp), DIMENSION(jpj,jpk) ::   u_cxebnd, v_cxebnd    !: Zonal component of the phase speed ratio computed with
110      !                                                    !  radiation of u and v velocity (respectively) at the
111      !                                                    !  east open boundary (u_cxebnd = cx rdt )
112   REAL(wp), DIMENSION(jpj,jpk) ::   uemsk, vemsk, temsk   !: 2D mask for the East OB
113
114   ! Note that those arrays are optimized for mpp case
115   ! (hence the dimension jpj is the size of one processor subdomain)
116
117   !!--------------------
118   !! West open boundary
119   !!--------------------
120   INTEGER ::   niw0  , niw1       !: do loop index in mpp case for jpiwob
121   INTEGER ::   niw0p1, niw1p1     !: do loop index in mpp case for jpiwob+1
122   INTEGER ::   njw0  , njw1       !: do loop index in mpp case for jpjwd, jpjwf
123   INTEGER ::   njw0p1, njw1m1     !: do loop index in mpp case for jpjwdp1,jpjwfm1
124   INTEGER ::   njw1m2, njw0m1     !: do loop index in mpp case for jpjwfm2,jpjwd
125
126   REAL(wp), DIMENSION(jpj) ::   &  !:
127      sshfow,           & !: now climatology of the west boundary sea surface height
128      ubtfow,vbtfow       !: now climatology of the west boundary barotropic transport
129
130   REAL(wp), DIMENSION(jpj,jpk) ::   &  !:
131      ufow, vfow,       & !: now climatology of the west velocities
132      tfow, sfow,       & !: now climatology of the west temperature and salinity
133      ucliw               !: baroclinic componant of the zonal velocity after the radiation
134      !                   !  in the obcdyn.F90 routine
135
136   REAL(wp), DIMENSION(jpi,jpj) ::   sshfow_b    !: west boundary ssh correction averaged over the barotropic loop
137      !                                          !  (if Flather's algoritm applied at open boundary)
138
139   !!-------------------------------
140   !! Arrays for radiative West OBC
141   !!-------------------------------
142   REAL(wp), DIMENSION(jpj,jpk,3,3) ::   uwbnd, vwbnd     !: baroclinic u & v components of the velocity over 3 rows
143      !                                                   !  and 3 time step (now, before, and before before)
144   REAL(wp), DIMENSION(jpj,jpk,2,2) ::   twbnd, swbnd     !: west boundary temperature and salinity over 2 rows and
145      !                                                   !  2 time step (now and before)
146   REAL(wp), DIMENSION(jpj,jpk) ::   u_cxwbnd, v_cxwbnd   !: Zonal component of the phase speed ratio computed with
147      !                                                   !  radiation of zonal and meridional velocity (respectively)
148      !                                                   !  at the west open boundary (u_cxwbnd = cx rdt )
149   REAL(wp), DIMENSION(jpj,jpk) ::   uwmsk, vwmsk, twmsk  !: 2D mask for the West OB
150
151   ! Note that those arrays are optimized for mpp case
152   ! (hence the dimension jpj is the size of one processor subdomain)
153
154   !!---------------------
155   !! North open boundary
156   !!---------------------
157   INTEGER ::   nin0  , nin1       !: do loop index in mpp case for jpind, jpinf
158   INTEGER ::   nin0p1, nin1m1     !: do loop index in mpp case for jpindp1, jpinfm1
159   INTEGER ::   nin1m2, nin0m1     !: do loop index in mpp case for jpinfm1-1,jpind
160   INTEGER ::   njn0  , njn1       !: do loop index in mpp case for jpnob
161   INTEGER ::   njn0p1, njn1p1     !: do loop index in mpp case for jpnob+1
162   INTEGER ::   njn0m1, njn1m1     !: do loop index in mpp case for jpnob-1
163
164   REAL(wp), DIMENSION(jpi) ::   &  !:
165      sshfon,           & !: now climatology of the north boundary sea surface height
166      ubtfon,vbtfon       !: now climatology of the north boundary barotropic transport
167
168   REAL(wp), DIMENSION(jpi,jpk) ::   &    !:
169      ufon, vfon,       & !: now climatology of the north boundary velocities
170      tfon, sfon,       & !: now climatology of the north boundary temperature and salinity
171      vclin               !: baroclinic componant of the meridian velocity after the radiation
172      !                   !  in yhe obcdyn.F90 routine
173
174   REAL(wp), DIMENSION(jpi,jpj) ::   sshfon_b      !: north boundary ssh correction averaged over the barotropic loop
175      !                                            !  (if Flather's algoritm applied at open boundary)
176
177   !!--------------------------------
178   !! Arrays for radiative North OBC
179   !!--------------------------------
180   REAL(wp), DIMENSION(jpi,jpk,3,3) ::   unbnd, vnbnd      !: baroclinic u & v components of the velocity over 3
181      !                                                    !  rows and 3 time step (now, before, and before before)
182   REAL(wp), DIMENSION(jpi,jpk,2,2) ::   tnbnd, snbnd      !: north boundary temperature and salinity over
183      !                                                    !  2 rows and 2 time step (now and before)
184   REAL(wp), DIMENSION(jpi,jpk) ::   u_cynbnd, v_cynbnd    !: Meridional component of the phase speed ratio compu-
185      !                                                    !  ted with radiation of zonal and meridional velocity
186      !                                                    !  (respectively) at the north OB (u_cynbnd = cx rdt )
187   REAL(wp), DIMENSION(jpi,jpk) ::   unmsk, vnmsk, tnmsk   !: 2D mask for the North OB
188
189   ! Note that those arrays are optimized for mpp case
190   ! (hence the dimension jpj is the size of one processor subdomain)
191   
192   !!---------------------
193   !! South open boundary
194   !!---------------------
195   INTEGER ::   nis0  , nis1       !: do loop index in mpp case for jpisd, jpisf
196   INTEGER ::   nis0p1, nis1m1     !: do loop index in mpp case for jpisdp1, jpisfm1
197   INTEGER ::   nis1m2, nis0m1     !: do loop index in mpp case for jpisfm1-1,jpisd
198   INTEGER ::   njs0  , njs1       !: do loop index in mpp case for jpsob
199   INTEGER ::   njs0p1, njs1p1     !: do loop index in mpp case for jpsob+1
200
201   REAL(wp), DIMENSION(jpi) ::    &   !:
202      sshfos,           & !: now climatology of the south boundary sea surface height
203      ubtfos,vbtfos       !: now climatology of the south boundary barotropic transport
204
205   REAL(wp), DIMENSION(jpi,jpk) ::    &   !:
206      ufos, vfos,       & !: now climatology of the south boundary velocities
207      tfos, sfos,       & !: now climatology of the south boundary temperature and salinity
208      vclis               !: baroclinic componant of the meridian velocity after the radiation
209      !                   !  in the obcdyn.F90 routine
210
211   REAL(wp), DIMENSION(jpi,jpj) ::   sshfos_b     !: south boundary ssh correction averaged over the barotropic loop
212      !                                           !  (if Flather's algoritm applied at open boundary)
213
214   !!--------------------------------
215   !! Arrays for radiative South OBC   (computed by the forward time step in dynspg)
216   !!--------------------------------
217   REAL(wp), DIMENSION(jpi,jpk,3,3) ::   usbnd, vsbnd     !: baroclinic u & v components of the velocity over 3
218      !                                                   !  rows and 3 time step (now, before, and before before)
219   REAL(wp), DIMENSION(jpi,jpk,2,2) ::   tsbnd, ssbnd     !: south boundary temperature and salinity over
220      !                                                   !  2 rows and 2 time step (now and before)
221   REAL(wp), DIMENSION(jpi,jpk) ::   u_cysbnd, v_cysbnd   !: Meridional component of the phase speed ratio
222      !                                                   !  computed with radiation of zonal and meridional velocity
223      !                                                   !  (repsectively) at the south OB (u_cynbnd = cx rdt )
224   REAL(wp), DIMENSION(jpi,jpk) ::   usmsk, vsmsk, tsmsk  !: 2D mask for the South OB
225
226#else
227   !!----------------------------------------------------------------------
228   !!   Default option :                                       Empty module
229   !!----------------------------------------------------------------------
230#endif
231
232   !!----------------------------------------------------------------------
233   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)
234   !! $Id$
235   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
236   !!======================================================================
237END MODULE obc_oce
Note: See TracBrowser for help on using the repository browser.