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.
bdydta.F90 in branches/2013/dev_LOCEAN_2013/NEMOGCM/NEMO/OPA_SRC/BDY – NEMO

source: branches/2013/dev_LOCEAN_2013/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90 @ 4148

Last change on this file since 4148 was 4148, checked in by cetlod, 10 years ago

merge in trunk changes between r3853 and r3940 and commit the changes, see ticket #1169
  • Property svn:keywords set to Id
File size: 31.2 KB
Line 
1MODULE bdydta
2   !!======================================================================
3   !!                       ***  MODULE bdydta  ***
4   !! Open boundary data : read the data for the unstructured open boundaries.
5   !!======================================================================
6   !! History :  1.0  !  2005-01  (J. Chanut, A. Sellar)  Original code
7   !!             -   !  2007-01  (D. Storkey) Update to use IOM module
8   !!             -   !  2007-07  (D. Storkey) add bdy_dta_fla
9   !!            3.0  !  2008-04  (NEMO team)  add in the reference version
10   !!            3.3  !  2010-09  (E.O'Dea) modifications for Shelf configurations
11   !!            3.3  !  2010-09  (D.Storkey) add ice boundary conditions
12   !!            3.4  !  2011     (D. Storkey) rewrite in preparation for OBC-BDY merge
13   !!----------------------------------------------------------------------
14#if defined key_bdy
15   !!----------------------------------------------------------------------
16   !!   'key_bdy'                     Open Boundary Conditions
17   !!----------------------------------------------------------------------
18   !!    bdy_dta        : read external data along open boundaries from file
19   !!    bdy_dta_init   : initialise arrays etc for reading of external data
20   !!----------------------------------------------------------------------
21   USE wrk_nemo        ! Memory Allocation
22   USE timing          ! Timing
23   USE oce             ! ocean dynamics and tracers
24   USE dom_oce         ! ocean space and time domain
25   USE phycst          ! physical constants
26   USE bdy_oce         ! ocean open boundary conditions 
27   USE bdytides        ! tidal forcing at boundaries
28   USE fldread         ! read input fields
29   USE iom             ! IOM library
30   USE in_out_manager  ! I/O logical units
31#if defined key_lim2
32   USE ice_2
33#endif
34   USE sbcapr
35
36   IMPLICIT NONE
37   PRIVATE
38
39   PUBLIC   bdy_dta          ! routine called by step.F90 and dynspg_ts.F90
40   PUBLIC   bdy_dta_init     ! routine called by nemogcm.F90
41
42   INTEGER, ALLOCATABLE, DIMENSION(:)   ::   nb_bdy_fld        ! Number of fields to update for each boundary set.
43   INTEGER                              ::   nb_bdy_fld_sum    ! Total number of fields to update for all boundary sets.
44
45   LOGICAL,           DIMENSION(jp_bdy) ::   ln_full_vel_array ! =T => full velocities in 3D boundary conditions
46                                                               ! =F => baroclinic velocities in 3D boundary conditions
47
48   TYPE(FLD), PUBLIC, ALLOCATABLE, DIMENSION(:), TARGET ::   bf        ! structure of input fields (file informations, fields read)
49
50   TYPE(MAP_POINTER), ALLOCATABLE, DIMENSION(:) :: nbmap_ptr   ! array of pointers to nbmap
51
52#  include "domzgr_substitute.h90"
53   !!----------------------------------------------------------------------
54   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
55   !! $Id$
56   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
57   !!----------------------------------------------------------------------
58CONTAINS
59
60      SUBROUTINE bdy_dta( kt, jit, time_offset )
61      !!----------------------------------------------------------------------
62      !!                   ***  SUBROUTINE bdy_dta  ***
63      !!                   
64      !! ** Purpose :   Update external data for open boundary conditions
65      !!
66      !! ** Method  :   Use fldread.F90
67      !!               
68      !!----------------------------------------------------------------------
69      !!
70      INTEGER, INTENT( in )           ::   kt    ! ocean time-step index
71      INTEGER, INTENT( in ), OPTIONAL ::   jit   ! subcycle time-step index (for timesplitting option)
72      INTEGER, INTENT( in ), OPTIONAL ::   time_offset  ! time offset in units of timesteps. NB. if jit
73                                                        ! is present then units = subcycle timesteps.
74                                                        ! time_offset = 0 => get data at "now" time level
75                                                        ! time_offset = -1 => get data at "before" time level
76                                                        ! time_offset = +1 => get data at "after" time level
77                                                        ! etc.
78      !!
79      INTEGER     ::  ib_bdy, jfld, jstart, jend, ib, ii, ij, ik, igrd  ! local indices
80      INTEGER,          DIMENSION(jpbgrd) ::   ilen1 
81      INTEGER, POINTER, DIMENSION(:)      ::   nblen, nblenrim  ! short cuts
82      !!
83      !!---------------------------------------------------------------------------
84      !!
85      IF( nn_timing == 1 ) CALL timing_start('bdy_dta')
86
87      ! Initialise data arrays once for all from initial conditions where required
88      !---------------------------------------------------------------------------
89      IF( kt .eq. nit000 .and. .not. PRESENT(jit) ) THEN
90
91         ! Calculate depth-mean currents
92         !-----------------------------
93         CALL wrk_alloc(jpi,jpj,pu2d,pv2d) 
94
95         pu2d(:,:) = 0.e0
96         pv2d(:,:) = 0.e0
97
98         DO ik = 1, jpkm1   !! Vertically integrated momentum trends
99             pu2d(:,:) = pu2d(:,:) + fse3u(:,:,ik) * umask(:,:,ik) * un(:,:,ik)
100             pv2d(:,:) = pv2d(:,:) + fse3v(:,:,ik) * vmask(:,:,ik) * vn(:,:,ik)
101         END DO
102         pu2d(:,:) = pu2d(:,:) * hur(:,:)
103         pv2d(:,:) = pv2d(:,:) * hvr(:,:)
104         
105         DO ib_bdy = 1, nb_bdy
106
107            nblen => idx_bdy(ib_bdy)%nblen
108            nblenrim => idx_bdy(ib_bdy)%nblenrim
109
110            IF( nn_dyn2d(ib_bdy) .gt. 0 .and. nn_dyn2d_dta(ib_bdy) .eq. 0 ) THEN
111               ilen1(:) = nblen(:)
112               igrd = 1
113               DO ib = 1, ilen1(igrd)
114                  ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
115                  ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
116                  dta_bdy(ib_bdy)%ssh(ib) = sshn(ii,ij) * tmask(ii,ij,1)         
117               END DO
118               igrd = 2
119               DO ib = 1, ilen1(igrd)
120                  ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
121                  ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
122                  dta_bdy(ib_bdy)%u2d(ib) = pu2d(ii,ij) * umask(ii,ij,1)         
123               END DO
124               igrd = 3
125               DO ib = 1, ilen1(igrd)
126                  ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
127                  ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
128                  dta_bdy(ib_bdy)%v2d(ib) = pv2d(ii,ij) * vmask(ii,ij,1)         
129               END DO
130            ENDIF
131
132            IF( nn_dyn3d(ib_bdy) .gt. 0 .and. nn_dyn3d_dta(ib_bdy) .eq. 0 ) THEN
133               ilen1(:) = nblen(:)
134               igrd = 2 
135               DO ib = 1, ilen1(igrd)
136                  DO ik = 1, jpkm1
137                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
138                     ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
139                     dta_bdy(ib_bdy)%u3d(ib,ik) =  ( un(ii,ij,ik) - pu2d(ii,ij) ) * umask(ii,ij,ik)         
140                  END DO
141               END DO
142               igrd = 3 
143               DO ib = 1, ilen1(igrd)
144                  DO ik = 1, jpkm1
145                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
146                     ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
147                     dta_bdy(ib_bdy)%v3d(ib,ik) =  ( vn(ii,ij,ik) - pv2d(ii,ij) ) * vmask(ii,ij,ik)         
148                     END DO
149               END DO
150            ENDIF
151
152            IF( nn_tra(ib_bdy) .gt. 0 .and. nn_tra_dta(ib_bdy) .eq. 0 ) THEN
153               ilen1(:) = nblen(:)
154               igrd = 1                       ! Everything is at T-points here
155               DO ib = 1, ilen1(igrd)
156                  DO ik = 1, jpkm1
157                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
158                     ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
159                     dta_bdy(ib_bdy)%tem(ib,ik) = tsn(ii,ij,ik,jp_tem) * tmask(ii,ij,ik)         
160                     dta_bdy(ib_bdy)%sal(ib,ik) = tsn(ii,ij,ik,jp_sal) * tmask(ii,ij,ik)         
161                  END DO
162               END DO
163            ENDIF
164
165#if defined key_lim2
166            IF( nn_ice_lim2(ib_bdy) .gt. 0 .and. nn_ice_lim2_dta(ib_bdy) .eq. 0 ) THEN
167               ilen1(:) = nblen(:)
168               igrd = 1                       ! Everything is at T-points here
169               DO ib = 1, ilen1(igrd)
170                  ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
171                  ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
172                  dta_bdy(ib_bdy)%frld(ib) = frld(ii,ij) * tmask(ii,ij,1)         
173                  dta_bdy(ib_bdy)%hicif(ib) = hicif(ii,ij) * tmask(ii,ij,1)         
174                  dta_bdy(ib_bdy)%hsnif(ib) = hsnif(ii,ij) * tmask(ii,ij,1)         
175               END DO
176            ENDIF
177#endif
178
179         ENDDO ! ib_bdy
180
181         CALL wrk_dealloc(jpi,jpj,pu2d,pv2d) 
182
183      ENDIF ! kt .eq. nit000
184
185      ! update external data from files
186      !--------------------------------
187     
188      jstart = 1
189      DO ib_bdy = 1, nb_bdy   
190         IF( nn_dta(ib_bdy) .eq. 1 ) THEN ! skip this bit if no external data required
191     
192            IF( PRESENT(jit) ) THEN
193               ! Update barotropic boundary conditions only
194               ! jit is optional argument for fld_read and bdytide_update
195               IF( nn_dyn2d(ib_bdy) .gt. 0 ) THEN
196                  IF( nn_dyn2d_dta(ib_bdy) .eq. 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays
197                     dta_bdy(ib_bdy)%ssh(:) = 0.0
198                     dta_bdy(ib_bdy)%u2d(:) = 0.0
199                     dta_bdy(ib_bdy)%v2d(:) = 0.0
200                  ENDIF
201                  IF (nn_tra(ib_bdy).ne.4) THEN
202                     IF( nn_dyn2d_dta(ib_bdy) .EQ. 1 .OR. nn_dyn2d_dta(ib_bdy) .EQ. 3 .OR.  &
203                       & (ln_full_vel_array(ib_bdy) .AND. nn_dyn3d_dta(ib_bdy).eq.1) )THEN
204
205                        ! For the runoff case, no need to update the forcing (already done in the baroclinic part)
206                        jend = nb_bdy_fld(ib_bdy)
207                        IF ( nn_tra(ib_bdy) .GT. 0 .AND. nn_tra_dta(ib_bdy) .GE. 1 ) jend = jend - 2
208                        CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend),  &
209                                     & kit=jit, kt_offset=time_offset )
210                        IF ( nn_tra(ib_bdy) .GT. 0 .AND. nn_tra_dta(ib_bdy) .GE. 1 ) jend = jend + 2
211
212                        ! If full velocities in boundary data then split into barotropic and baroclinic data
213                        IF( ln_full_vel_array(ib_bdy) .AND.                                             &
214                          &    ( nn_dyn2d_dta(ib_bdy) .EQ. 1 .OR. nn_dyn2d_dta(ib_bdy) .EQ. 3 .OR.  &
215                          &      nn_dyn3d_dta(ib_bdy) .EQ. 1 ) )THEN
216
217                           igrd = 2                      ! zonal velocity
218                           dta_bdy(ib_bdy)%u2d(:) = 0.0
219                           DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
220                              ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
221                              ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
222                              DO ik = 1, jpkm1
223                                 dta_bdy(ib_bdy)%u2d(ib) = dta_bdy(ib_bdy)%u2d(ib) &
224                       &                          + fse3u(ii,ij,ik) * umask(ii,ij,ik) * dta_bdy(ib_bdy)%u3d(ib,ik)
225                              END DO
226                              dta_bdy(ib_bdy)%u2d(ib) =  dta_bdy(ib_bdy)%u2d(ib) * hur(ii,ij)
227                              DO ik = 1, jpkm1
228                                 dta_bdy(ib_bdy)%u3d(ib,ik) = dta_bdy(ib_bdy)%u3d(ib,ik) - dta_bdy(ib_bdy)%u2d(ib)
229                              END DO
230                           END DO
231                           igrd = 3                      ! meridional velocity
232                           dta_bdy(ib_bdy)%v2d(:) = 0.0
233                           DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
234                              ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
235                              ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
236                              DO ik = 1, jpkm1
237                                 dta_bdy(ib_bdy)%v2d(ib) = dta_bdy(ib_bdy)%v2d(ib) &
238                       &                       + fse3v(ii,ij,ik) * vmask(ii,ij,ik) * dta_bdy(ib_bdy)%v3d(ib,ik)
239                              END DO
240                              dta_bdy(ib_bdy)%v2d(ib) =  dta_bdy(ib_bdy)%v2d(ib) * hvr(ii,ij)
241                              DO ik = 1, jpkm1
242                                 dta_bdy(ib_bdy)%v3d(ib,ik) = dta_bdy(ib_bdy)%v3d(ib,ik) - dta_bdy(ib_bdy)%v2d(ib)
243                              END DO
244                           END DO
245                        ENDIF                   
246                     ENDIF
247                     IF( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN ! update tidal harmonic forcing
248                        CALL bdytide_update( kt=kt, idx=idx_bdy(ib_bdy), dta=dta_bdy(ib_bdy), td=tides(ib_bdy),   & 
249                          &                 jit=jit, time_offset=time_offset )
250                     ENDIF
251                  ENDIF
252               ENDIF
253            ELSE
254               IF (nn_tra(ib_bdy).eq.4) then      ! runoff condition
255                  jend = nb_bdy_fld(ib_bdy)
256                  CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend),  &
257                               & map=nbmap_ptr(jstart:jend), kt_offset=time_offset )
258                  !
259                  igrd = 2                      ! zonal velocity
260                  DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
261                     ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
262                     ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
263                     dta_bdy(ib_bdy)%u2d(ib) = dta_bdy(ib_bdy)%u2d(ib) / ( e2u(ii,ij) * hu_0(ii,ij) )
264                  END DO
265                  !
266                  igrd = 3                      ! meridional velocity
267                  DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
268                     ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
269                     ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
270                     dta_bdy(ib_bdy)%v2d(ib) = dta_bdy(ib_bdy)%v2d(ib) / ( e1v(ii,ij) * hv_0(ii,ij) )
271                  END DO
272               ELSE
273                  IF( nn_dyn2d(ib_bdy) .gt. 0 .and. nn_dyn2d_dta(ib_bdy) .eq. 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays
274                     dta_bdy(ib_bdy)%ssh(:) = 0.0
275                     dta_bdy(ib_bdy)%u2d(:) = 0.0
276                     dta_bdy(ib_bdy)%v2d(:) = 0.0
277                  ENDIF
278                  IF( nb_bdy_fld(ib_bdy) .gt. 0 ) THEN ! update external data
279                     jend = nb_bdy_fld(ib_bdy)
280                     CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), &
281                                  & map=nbmap_ptr(jstart:jend), kt_offset=time_offset )
282                  ENDIF
283                  ! If full velocities in boundary data then split into barotropic and baroclinic data
284                  IF( ln_full_vel_array(ib_bdy) .and.                                             &
285                    & ( nn_dyn2d_dta(ib_bdy) .EQ. 1 .OR. nn_dyn2d_dta(ib_bdy) .EQ. 3 .OR. &
286                    &   nn_dyn3d_dta(ib_bdy) .EQ. 1 ) ) THEN
287                     igrd = 2                      ! zonal velocity
288                     dta_bdy(ib_bdy)%u2d(:) = 0.0
289                     DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
290                        ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
291                        ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
292                        DO ik = 1, jpkm1
293                           dta_bdy(ib_bdy)%u2d(ib) = dta_bdy(ib_bdy)%u2d(ib) &
294                 &                       + fse3u(ii,ij,ik) * umask(ii,ij,ik) * dta_bdy(ib_bdy)%u3d(ib,ik)
295                        END DO
296                        dta_bdy(ib_bdy)%u2d(ib) =  dta_bdy(ib_bdy)%u2d(ib) * hur(ii,ij)
297                        DO ik = 1, jpkm1
298                           dta_bdy(ib_bdy)%u3d(ib,ik) = dta_bdy(ib_bdy)%u3d(ib,ik) - dta_bdy(ib_bdy)%u2d(ib)
299                        END DO
300                     END DO
301                     igrd = 3                      ! meridional velocity
302                     dta_bdy(ib_bdy)%v2d(:) = 0.0
303                     DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
304                        ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
305                        ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
306                        DO ik = 1, jpkm1
307                           dta_bdy(ib_bdy)%v2d(ib) = dta_bdy(ib_bdy)%v2d(ib) &
308                 &                       + fse3v(ii,ij,ik) * vmask(ii,ij,ik) * dta_bdy(ib_bdy)%v3d(ib,ik)
309                        END DO
310                        dta_bdy(ib_bdy)%v2d(ib) =  dta_bdy(ib_bdy)%v2d(ib) * hvr(ii,ij)
311                        DO ik = 1, jpkm1
312                           dta_bdy(ib_bdy)%v3d(ib,ik) = dta_bdy(ib_bdy)%v3d(ib,ik) - dta_bdy(ib_bdy)%v2d(ib)
313                        END DO
314                     END DO
315                  ENDIF
316                  IF( nn_dyn2d(ib_bdy) .gt. 0 .and. nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN ! update tidal harmonic forcing
317                     CALL bdytide_update( kt=kt, idx=idx_bdy(ib_bdy), dta=dta_bdy(ib_bdy),  &
318                                        & td=tides(ib_bdy), time_offset=time_offset )
319                  ENDIF
320               ENDIF
321            ENDIF
322            jstart = jend+1
323         END IF ! nn_dta(ib_bdy) = 1
324      END DO  ! ib_bdy
325
326      IF ( ln_apr_obc ) THEN
327         DO ib_bdy = 1, nb_bdy
328            IF (nn_tra(ib_bdy).NE.4)THEN
329               igrd = 1                      ! meridional velocity
330               DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
331                  ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
332                  ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
333                  dta_bdy(ib_bdy)%ssh(ib) = dta_bdy(ib_bdy)%ssh(ib) + ssh_ib(ii,ij)
334               ENDDO
335            ENDIF
336         ENDDO
337      ENDIF
338
339      IF( nn_timing == 1 ) CALL timing_stop('bdy_dta')
340
341      END SUBROUTINE bdy_dta
342
343
344      SUBROUTINE bdy_dta_init
345      !!----------------------------------------------------------------------
346      !!                   ***  SUBROUTINE bdy_dta_init  ***
347      !!                   
348      !! ** Purpose :   Initialise arrays for reading of external data
349      !!                for open boundary conditions
350      !!
351      !! ** Method  :   Use fldread.F90
352      !!               
353      !!----------------------------------------------------------------------
354      USE dynspg_oce, ONLY: lk_dynspg_ts
355      !!
356      INTEGER     ::  ib_bdy, jfld, jstart, jend, ierror  ! local indices
357      INTEGER      ::   ios                               ! Local integer output status for namelist read
358      !!
359      CHARACTER(len=100)                     ::   cn_dir        ! Root directory for location of data files
360      CHARACTER(len=100), DIMENSION(nb_bdy)  ::   cn_dir_array  ! Root directory for location of data files
361      LOGICAL                                ::   ln_full_vel   ! =T => full velocities in 3D boundary data
362                                                                ! =F => baroclinic velocities in 3D boundary data
363      INTEGER                                ::   ilen_global   ! Max length required for global bdy dta arrays
364      INTEGER,              DIMENSION(jpbgrd) ::  ilen0         ! size of local arrays
365      INTEGER, ALLOCATABLE, DIMENSION(:)     ::   ilen1, ilen3  ! size of 1st and 3rd dimensions of local arrays
366      INTEGER, ALLOCATABLE, DIMENSION(:)     ::   ibdy           ! bdy set for a particular jfld
367      INTEGER, ALLOCATABLE, DIMENSION(:)     ::   igrid         ! index for grid type (1,2,3 = T,U,V)
368      INTEGER, POINTER, DIMENSION(:)         ::   nblen, nblenrim  ! short cuts
369      TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) ::   blf_i         !  array of namelist information structures
370      TYPE(FLD_N) ::   bn_tem, bn_sal, bn_u3d, bn_v3d   !
371      TYPE(FLD_N) ::   bn_ssh, bn_u2d, bn_v2d           ! informations about the fields to be read
372#if defined key_lim2
373      TYPE(FLD_N) ::   bn_frld, bn_hicif, bn_hsnif      !
374#endif
375      NAMELIST/nambdy_dta/ cn_dir, bn_tem, bn_sal, bn_u3d, bn_v3d, bn_ssh, bn_u2d, bn_v2d 
376#if defined key_lim2
377      NAMELIST/nambdy_dta/ bn_frld, bn_hicif, bn_hsnif
378#endif
379      NAMELIST/nambdy_dta/ ln_full_vel
380      !!---------------------------------------------------------------------------
381
382      IF( nn_timing == 1 ) CALL timing_start('bdy_dta_init')
383
384      IF(lwp) WRITE(numout,*)
385      IF(lwp) WRITE(numout,*) 'bdy_dta_ini : initialization of data at the open boundaries'
386      IF(lwp) WRITE(numout,*) '~~~~~~~~~~'
387      IF(lwp) WRITE(numout,*) ''
388
389      ! Set nn_dta
390      DO ib_bdy = 1, nb_bdy
391         nn_dta(ib_bdy) = MAX(  nn_dyn2d_dta(ib_bdy)       &
392                               ,nn_dyn3d_dta(ib_bdy)       &
393                               ,nn_tra_dta(ib_bdy)         &
394#if defined key_lim2
395                               ,nn_ice_lim2_dta(ib_bdy)    &
396#endif
397                              )
398         IF(nn_dta(ib_bdy) .gt. 1) nn_dta(ib_bdy) = 1
399      END DO
400
401      ! Work out upper bound of how many fields there are to read in and allocate arrays
402      ! ---------------------------------------------------------------------------
403      ALLOCATE( nb_bdy_fld(nb_bdy) )
404      nb_bdy_fld(:) = 0
405      DO ib_bdy = 1, nb_bdy         
406         IF( nn_dyn2d(ib_bdy) .gt. 0 .and. ( nn_dyn2d_dta(ib_bdy) .eq. 1 .or. nn_dyn2d_dta(ib_bdy) .eq. 3 ) ) THEN
407            nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 3
408         ENDIF
409         IF( nn_dyn3d(ib_bdy) .gt. 0 .and. nn_dyn3d_dta(ib_bdy) .eq. 1 ) THEN
410            nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 2
411         ENDIF
412         IF( nn_tra(ib_bdy) .gt. 0 .and. nn_tra_dta(ib_bdy) .eq. 1  ) THEN
413            nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 2
414         ENDIF
415#if defined key_lim2
416         IF( nn_ice_lim2(ib_bdy) .gt. 0 .and. nn_ice_lim2_dta(ib_bdy) .eq. 1  ) THEN
417            nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 3
418         ENDIF
419#endif               
420         IF(lwp) WRITE(numout,*) 'Maximum number of files to open =',nb_bdy_fld(ib_bdy)
421      ENDDO           
422
423      nb_bdy_fld_sum = SUM( nb_bdy_fld )
424
425      ALLOCATE( bf(nb_bdy_fld_sum), STAT=ierror )
426      IF( ierror > 0 ) THEN   
427         CALL ctl_stop( 'bdy_dta: unable to allocate bf structure' )   ;   RETURN 
428      ENDIF
429      ALLOCATE( blf_i(nb_bdy_fld_sum), STAT=ierror )
430      IF( ierror > 0 ) THEN   
431         CALL ctl_stop( 'bdy_dta: unable to allocate blf_i structure' )   ;   RETURN 
432      ENDIF
433      ALLOCATE( nbmap_ptr(nb_bdy_fld_sum), STAT=ierror )
434      IF( ierror > 0 ) THEN   
435         CALL ctl_stop( 'bdy_dta: unable to allocate nbmap_ptr structure' )   ;   RETURN 
436      ENDIF
437      ALLOCATE( ilen1(nb_bdy_fld_sum), ilen3(nb_bdy_fld_sum) ) 
438      ALLOCATE( ibdy(nb_bdy_fld_sum) ) 
439      ALLOCATE( igrid(nb_bdy_fld_sum) ) 
440
441      ! Read namelists
442      ! --------------
443      REWIND(numnam_ref)
444      REWIND(numnam_cfg)
445      jfld = 0 
446      DO ib_bdy = 1, nb_bdy         
447         IF( nn_dta(ib_bdy) .eq. 1 ) THEN
448
449            READ  ( numnam_ref, nambdy_dta, IOSTAT = ios, ERR = 901)
450901         IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy_dta in reference namelist', lwp )
451
452            READ  ( numnam_cfg, nambdy_dta, IOSTAT = ios, ERR = 902 )
453902         IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy_dta in configuration namelist', lwp )
454            WRITE ( numond, nambdy_dta )
455
456            cn_dir_array(ib_bdy) = cn_dir
457            ln_full_vel_array(ib_bdy) = ln_full_vel
458
459            nblen => idx_bdy(ib_bdy)%nblen
460            nblenrim => idx_bdy(ib_bdy)%nblenrim
461
462            ! Only read in necessary fields for this set.
463            ! Important that barotropic variables come first.
464            IF( nn_dyn2d(ib_bdy) .gt. 0 .and. ( nn_dyn2d_dta(ib_bdy) .eq. 1 .or. nn_dyn2d_dta(ib_bdy) .eq. 3 ) ) THEN
465
466               IF( nn_dyn2d(ib_bdy) .ne. jp_frs .and. nn_tra(ib_bdy) .ne. 4 ) THEN ! runoff condition : no ssh reading
467                  jfld = jfld + 1
468                  blf_i(jfld) = bn_ssh
469                  ibdy(jfld) = ib_bdy
470                  igrid(jfld) = 1
471                  ilen1(jfld) = nblen(igrid(jfld))
472                  ilen3(jfld) = 1
473               ENDIF
474
475               IF( .not. ln_full_vel_array(ib_bdy) ) THEN
476                  jfld = jfld + 1
477                  blf_i(jfld) = bn_u2d
478                  ibdy(jfld) = ib_bdy
479                  igrid(jfld) = 2
480                  ilen1(jfld) = nblen(igrid(jfld))
481                  ilen3(jfld) = 1
482
483                  jfld = jfld + 1
484                  blf_i(jfld) = bn_v2d
485                  ibdy(jfld) = ib_bdy
486                  igrid(jfld) = 3
487                  ilen1(jfld) = nblen(igrid(jfld))
488                  ilen3(jfld) = 1
489               ENDIF
490
491            ENDIF
492
493            ! baroclinic velocities
494            IF( ( nn_dyn3d(ib_bdy) .gt. 0 .and. nn_dyn3d_dta(ib_bdy) .eq. 1 ) .or. &
495           &      ( ln_full_vel_array(ib_bdy) .and. nn_dyn2d(ib_bdy) .gt. 0 .and.  &
496           &        ( nn_dyn2d_dta(ib_bdy) .eq. 1 .or. nn_dyn2d_dta(ib_bdy) .eq. 3 ) ) ) THEN
497
498               jfld = jfld + 1
499               blf_i(jfld) = bn_u3d
500               ibdy(jfld) = ib_bdy
501               igrid(jfld) = 2
502               ilen1(jfld) = nblen(igrid(jfld))
503               ilen3(jfld) = jpk
504
505               jfld = jfld + 1
506               blf_i(jfld) = bn_v3d
507               ibdy(jfld) = ib_bdy
508               igrid(jfld) = 3
509               ilen1(jfld) = nblen(igrid(jfld))
510               ilen3(jfld) = jpk
511
512            ENDIF
513
514            ! temperature and salinity
515            IF( nn_tra(ib_bdy) .gt. 0 .and. nn_tra_dta(ib_bdy) .eq. 1 ) THEN
516
517               jfld = jfld + 1
518               blf_i(jfld) = bn_tem
519               ibdy(jfld) = ib_bdy
520               igrid(jfld) = 1
521               ilen1(jfld) = nblen(igrid(jfld))
522               ilen3(jfld) = jpk
523
524               jfld = jfld + 1
525               blf_i(jfld) = bn_sal
526               ibdy(jfld) = ib_bdy
527               igrid(jfld) = 1
528               ilen1(jfld) = nblen(igrid(jfld))
529               ilen3(jfld) = jpk
530
531            ENDIF
532
533#if defined key_lim2
534            ! sea ice
535            IF( nn_ice_lim2(ib_bdy) .gt. 0 .and. nn_ice_lim2_dta(ib_bdy) .eq. 1 ) THEN
536
537               jfld = jfld + 1
538               blf_i(jfld) = bn_frld
539               ibdy(jfld) = ib_bdy
540               igrid(jfld) = 1
541               ilen1(jfld) = nblen(igrid(jfld))
542               ilen3(jfld) = 1
543
544               jfld = jfld + 1
545               blf_i(jfld) = bn_hicif
546               ibdy(jfld) = ib_bdy
547               igrid(jfld) = 1
548               ilen1(jfld) = nblen(igrid(jfld))
549               ilen3(jfld) = 1
550
551               jfld = jfld + 1
552               blf_i(jfld) = bn_hsnif
553               ibdy(jfld) = ib_bdy
554               igrid(jfld) = 1
555               ilen1(jfld) = nblen(igrid(jfld))
556               ilen3(jfld) = 1
557
558            ENDIF
559#endif
560            ! Recalculate field counts
561            !-------------------------
562            IF( ib_bdy .eq. 1 ) THEN
563               nb_bdy_fld_sum = 0
564               nb_bdy_fld(ib_bdy) = jfld
565               nb_bdy_fld_sum     = jfld             
566            ELSE
567               nb_bdy_fld(ib_bdy) = jfld - nb_bdy_fld_sum
568               nb_bdy_fld_sum = nb_bdy_fld_sum + nb_bdy_fld(ib_bdy)
569            ENDIF
570
571         ENDIF ! nn_dta .eq. 1
572      ENDDO ! ib_bdy
573
574      DO jfld = 1, nb_bdy_fld_sum
575         ALLOCATE( bf(jfld)%fnow(ilen1(jfld),1,ilen3(jfld)) )
576         IF( blf_i(jfld)%ln_tint ) ALLOCATE( bf(jfld)%fdta(ilen1(jfld),1,ilen3(jfld),2) )
577         nbmap_ptr(jfld)%ptr => idx_bdy(ibdy(jfld))%nbmap(:,igrid(jfld))
578      ENDDO
579
580      ! fill bf with blf_i and control print
581      !-------------------------------------
582      jstart = 1
583      DO ib_bdy = 1, nb_bdy
584         jend = nb_bdy_fld(ib_bdy) 
585         CALL fld_fill( bf(jstart:jend), blf_i(jstart:jend), cn_dir_array(ib_bdy), 'bdy_dta',   &
586         &              'open boundary conditions', 'nambdy_dta' )
587         jstart = jend + 1
588      ENDDO
589
590      ! Initialise local boundary data arrays
591      ! nn_xxx_dta=0 : allocate space - will be filled from initial conditions later
592      ! nn_xxx_dta=1 : point to "fnow" arrays
593      !-------------------------------------
594
595      jfld = 0
596      DO ib_bdy=1, nb_bdy
597
598         nblen => idx_bdy(ib_bdy)%nblen
599         nblenrim => idx_bdy(ib_bdy)%nblenrim
600
601         IF (nn_dyn2d(ib_bdy) .gt. 0) THEN
602            IF( nn_dyn2d_dta(ib_bdy) .eq. 0 .or. nn_dyn2d_dta(ib_bdy) .eq. 2 .or. ln_full_vel_array(ib_bdy) ) THEN
603               ilen0(1:3) = nblen(1:3)
604               ALLOCATE( dta_bdy(ib_bdy)%u2d(ilen0(2)) )
605               ALLOCATE( dta_bdy(ib_bdy)%v2d(ilen0(3)) )
606               IF ( nn_dyn2d(ib_bdy) .ne. jp_frs .and. (nn_dyn2d_dta(ib_bdy).eq.1.or.nn_dyn2d_dta(ib_bdy).eq.3) ) THEN
607                  jfld = jfld + 1
608                  dta_bdy(ib_bdy)%ssh => bf(jfld)%fnow(:,1,1)
609               ELSE
610                  ALLOCATE( dta_bdy(ib_bdy)%ssh(nblen(1)) )
611               ENDIF
612            ELSE
613               IF( nn_dyn2d(ib_bdy) .ne. jp_frs ) THEN
614                  jfld = jfld + 1
615                  dta_bdy(ib_bdy)%ssh => bf(jfld)%fnow(:,1,1)
616               ENDIF
617               jfld = jfld + 1
618               dta_bdy(ib_bdy)%u2d => bf(jfld)%fnow(:,1,1)
619               jfld = jfld + 1
620               dta_bdy(ib_bdy)%v2d => bf(jfld)%fnow(:,1,1)
621            ENDIF
622         ENDIF
623
624         IF ( nn_dyn3d(ib_bdy) .gt. 0 .and. nn_dyn3d_dta(ib_bdy) .eq. 0 ) THEN
625            ilen0(1:3) = nblen(1:3)
626            ALLOCATE( dta_bdy(ib_bdy)%u3d(ilen0(2),jpk) )
627            ALLOCATE( dta_bdy(ib_bdy)%v3d(ilen0(3),jpk) )
628         ENDIF
629         IF ( ( nn_dyn3d(ib_bdy) .gt. 0 .and. nn_dyn3d_dta(ib_bdy) .eq. 1 ).or. &
630           &  ( ln_full_vel_array(ib_bdy) .and. nn_dyn2d(ib_bdy) .gt. 0 .and.   &
631           &    ( nn_dyn2d_dta(ib_bdy) .eq. 1 .or. nn_dyn2d_dta(ib_bdy) .eq. 3 ) ) ) THEN
632            jfld = jfld + 1
633            dta_bdy(ib_bdy)%u3d => bf(jfld)%fnow(:,1,:)
634            jfld = jfld + 1
635            dta_bdy(ib_bdy)%v3d => bf(jfld)%fnow(:,1,:)
636         ENDIF
637
638         IF (nn_tra(ib_bdy) .gt. 0) THEN
639            IF( nn_tra_dta(ib_bdy) .eq. 0 ) THEN
640               ilen0(1:3) = nblen(1:3)
641               ALLOCATE( dta_bdy(ib_bdy)%tem(ilen0(1),jpk) )
642               ALLOCATE( dta_bdy(ib_bdy)%sal(ilen0(1),jpk) )
643            ELSE
644               jfld = jfld + 1
645               dta_bdy(ib_bdy)%tem => bf(jfld)%fnow(:,1,:)
646               jfld = jfld + 1
647               dta_bdy(ib_bdy)%sal => bf(jfld)%fnow(:,1,:)
648            ENDIF
649         ENDIF
650
651#if defined key_lim2
652         IF (nn_ice_lim2(ib_bdy) .gt. 0) THEN
653            IF( nn_ice_lim2_dta(ib_bdy) .eq. 0 ) THEN
654               ilen0(1:3) = nblen(1:3)
655               ALLOCATE( dta_bdy(ib_bdy)%frld(ilen0(1)) )
656               ALLOCATE( dta_bdy(ib_bdy)%hicif(ilen0(1)) )
657               ALLOCATE( dta_bdy(ib_bdy)%hsnif(ilen0(1)) )
658            ELSE
659               jfld = jfld + 1
660               dta_bdy(ib_bdy)%frld  => bf(jfld)%fnow(:,1,1)
661               jfld = jfld + 1
662               dta_bdy(ib_bdy)%hicif => bf(jfld)%fnow(:,1,1)
663               jfld = jfld + 1
664               dta_bdy(ib_bdy)%hsnif => bf(jfld)%fnow(:,1,1)
665            ENDIF
666         ENDIF
667#endif
668
669      ENDDO ! ib_bdy
670
671      IF( nn_timing == 1 ) CALL timing_stop('bdy_dta_init')
672
673      END SUBROUTINE bdy_dta_init
674
675#else
676   !!----------------------------------------------------------------------
677   !!   Dummy module                   NO Open Boundary Conditions
678   !!----------------------------------------------------------------------
679CONTAINS
680   SUBROUTINE bdy_dta( kt, jit, time_offset ) ! Empty routine
681      INTEGER, INTENT( in )           ::   kt   
682      INTEGER, INTENT( in ), OPTIONAL ::   jit   
683      INTEGER, INTENT( in ), OPTIONAL ::   time_offset
684      WRITE(*,*) 'bdy_dta: You should not have seen this print! error?', kt
685   END SUBROUTINE bdy_dta
686   SUBROUTINE bdy_dta_init()                  ! Empty routine
687      WRITE(*,*) 'bdy_dta_init: You should not have seen this print! error?'
688   END SUBROUTINE bdy_dta_init
689#endif
690
691   !!==============================================================================
692END MODULE bdydta
Note: See TracBrowser for help on using the repository browser.