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/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY – NEMO

source: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90 @ 8306

Last change on this file since 8306 was 8306, checked in by clem, 7 years ago

step1: remove LIM2 from the code

  • Property svn:keywords set to Id
File size: 37.3 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   !!            3.6  !  2012-01  (C. Rousset) add ice boundary conditions for lim3
14   !!----------------------------------------------------------------------
15   !!    bdy_dta        : read external data along open boundaries from file
16   !!    bdy_dta_init   : initialise arrays etc for reading of external data
17   !!----------------------------------------------------------------------
18   USE timing          ! Timing
19   USE oce             ! ocean dynamics and tracers
20   USE dom_oce         ! ocean space and time domain
21   USE phycst          ! physical constants
22   USE bdy_oce         ! ocean open boundary conditions 
23   USE bdytides        ! tidal forcing at boundaries
24   USE fldread         ! read input fields
25   USE iom             ! IOM library
26   USE in_out_manager  ! I/O logical units
27#if defined key_lim3
28   USE ice
29   USE limvar          ! redistribute ice input into categories
30#endif
31   USE sbcapr
32   USE sbctide         ! Tidal forcing or not
33
34   IMPLICIT NONE
35   PRIVATE
36
37   PUBLIC   bdy_dta          ! routine called by step.F90 and dynspg_ts.F90
38   PUBLIC   bdy_dta_init     ! routine called by nemogcm.F90
39
40   INTEGER, ALLOCATABLE, DIMENSION(:)   ::   nb_bdy_fld        ! Number of fields to update for each boundary set.
41   INTEGER                              ::   nb_bdy_fld_sum    ! Total number of fields to update for all boundary sets.
42   LOGICAL,           DIMENSION(jp_bdy) ::   ln_full_vel_array ! =T => full velocities in 3D boundary conditions
43                                                               ! =F => baroclinic velocities in 3D boundary conditions
44!$AGRIF_DO_NOT_TREAT
45   TYPE(FLD), PUBLIC, ALLOCATABLE, DIMENSION(:), TARGET ::   bf        ! structure of input fields (file informations, fields read)
46!$AGRIF_END_DO_NOT_TREAT
47   TYPE(MAP_POINTER), ALLOCATABLE, DIMENSION(:) :: nbmap_ptr   ! array of pointers to nbmap
48
49#if defined key_lim3
50   LOGICAL :: ll_bdylim3                  ! determine whether ice input is 1cat (F) or Xcat (T) type
51   INTEGER :: jfld_hti, jfld_hts, jfld_ai ! indices of ice thickness, snow thickness and concentration in bf structure
52#endif
53
54   !!----------------------------------------------------------------------
55   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
56   !! $Id$
57   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
58   !!----------------------------------------------------------------------
59CONTAINS
60
61      SUBROUTINE bdy_dta( kt, jit, time_offset )
62      !!----------------------------------------------------------------------
63      !!                   ***  SUBROUTINE bdy_dta  ***
64      !!                   
65      !! ** Purpose :   Update external data for open boundary conditions
66      !!
67      !! ** Method  :   Use fldread.F90
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, jl  ! local indices
80      INTEGER,          DIMENSION(jpbgrd) ::   ilen1 
81      INTEGER, POINTER, DIMENSION(:)      ::   nblen, nblenrim  ! short cuts
82      TYPE(OBC_DATA), POINTER             ::   dta              ! short cut
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 == nit000 .AND. .NOT.PRESENT(jit) ) THEN
90
91         ! Calculate depth-mean currents
92         !-----------------------------
93         
94         DO ib_bdy = 1, nb_bdy
95            !
96            nblen => idx_bdy(ib_bdy)%nblen
97            nblenrim => idx_bdy(ib_bdy)%nblenrim
98            dta => dta_bdy(ib_bdy)
99
100            IF( nn_dyn2d_dta(ib_bdy) == 0 ) THEN
101               ilen1(:) = nblen(:)
102               IF( dta%ll_ssh ) THEN
103                  igrd = 1
104                  DO ib = 1, ilen1(igrd)
105                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
106                     ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
107                     dta_bdy(ib_bdy)%ssh(ib) = sshn(ii,ij) * tmask(ii,ij,1)         
108                  END DO
109               END IF
110               IF( dta%ll_u2d ) THEN
111                  igrd = 2
112                  DO ib = 1, ilen1(igrd)
113                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
114                     ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
115                     dta_bdy(ib_bdy)%u2d(ib) = un_b(ii,ij) * umask(ii,ij,1)         
116                  END DO
117               END IF
118               IF( dta%ll_v2d ) THEN
119                  igrd = 3
120                  DO ib = 1, ilen1(igrd)
121                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
122                     ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
123                     dta_bdy(ib_bdy)%v2d(ib) = vn_b(ii,ij) * vmask(ii,ij,1)         
124                  END DO
125               END IF
126            ENDIF
127
128            IF( nn_dyn3d_dta(ib_bdy) == 0 ) THEN
129               ilen1(:) = nblen(:)
130               IF( dta%ll_u3d ) THEN
131                  igrd = 2 
132                  DO ib = 1, ilen1(igrd)
133                     DO ik = 1, jpkm1
134                        ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
135                        ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
136                        dta_bdy(ib_bdy)%u3d(ib,ik) =  ( un(ii,ij,ik) - un_b(ii,ij) ) * umask(ii,ij,ik)         
137                     END DO
138                  END DO
139               END IF
140               IF( dta%ll_v3d ) THEN
141                  igrd = 3 
142                  DO ib = 1, ilen1(igrd)
143                     DO ik = 1, jpkm1
144                        ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
145                        ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
146                        dta_bdy(ib_bdy)%v3d(ib,ik) =  ( vn(ii,ij,ik) - vn_b(ii,ij) ) * vmask(ii,ij,ik)         
147                        END DO
148                  END DO
149               END IF
150            ENDIF
151
152            IF( nn_tra_dta(ib_bdy) == 0 ) THEN
153               ilen1(:) = nblen(:)
154               IF( dta%ll_tem ) THEN
155                  igrd = 1 
156                  DO ib = 1, ilen1(igrd)
157                     DO ik = 1, jpkm1
158                        ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
159                        ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
160                        dta_bdy(ib_bdy)%tem(ib,ik) = tsn(ii,ij,ik,jp_tem) * tmask(ii,ij,ik)         
161                     END DO
162                  END DO
163               END IF
164               IF( dta%ll_sal ) THEN
165                  igrd = 1 
166                  DO ib = 1, ilen1(igrd)
167                     DO ik = 1, jpkm1
168                        ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
169                        ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
170                        dta_bdy(ib_bdy)%sal(ib,ik) = tsn(ii,ij,ik,jp_sal) * tmask(ii,ij,ik)         
171                     END DO
172                  END DO
173               END IF
174            ENDIF
175
176#if defined key_lim3
177            IF( nn_ice_lim_dta(ib_bdy) == 0 ) THEN
178               ilen1(:) = nblen(:)
179               IF( dta%ll_a_i ) THEN
180                  igrd = 1   
181                  DO jl = 1, jpl
182                     DO ib = 1, ilen1(igrd)
183                        ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
184                        ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
185                        dta_bdy(ib_bdy)%a_i (ib,jl) =  a_i(ii,ij,jl) * tmask(ii,ij,1) 
186                     END DO
187                  END DO
188               ENDIF
189               IF( dta%ll_ht_i ) THEN
190                  igrd = 1   
191                  DO jl = 1, jpl
192                     DO ib = 1, ilen1(igrd)
193                        ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
194                        ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
195                        dta_bdy(ib_bdy)%ht_i (ib,jl) =  ht_i(ii,ij,jl) * tmask(ii,ij,1) 
196                     END DO
197                  END DO
198               ENDIF
199               IF( dta%ll_ht_s ) THEN
200                  igrd = 1   
201                  DO jl = 1, jpl
202                     DO ib = 1, ilen1(igrd)
203                        ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
204                        ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
205                        dta_bdy(ib_bdy)%ht_s (ib,jl) =  ht_s(ii,ij,jl) * tmask(ii,ij,1) 
206                     END DO
207                  END DO
208               ENDIF
209            ENDIF
210#endif
211         END DO ! ib_bdy
212         !
213      ENDIF ! kt == nit000
214
215      ! update external data from files
216      !--------------------------------
217     
218      jstart = 1
219      DO ib_bdy = 1, nb_bdy   
220         dta => dta_bdy(ib_bdy)
221         IF( nn_dta(ib_bdy) == 1 ) THEN ! skip this bit if no external data required
222     
223            IF( PRESENT(jit) ) THEN
224               ! Update barotropic boundary conditions only
225               ! jit is optional argument for fld_read and bdytide_update
226               IF( cn_dyn2d(ib_bdy) /= 'none' ) THEN
227                  IF( nn_dyn2d_dta(ib_bdy) == 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays
228                     IF( dta%ll_ssh ) dta%ssh(:) = 0._wp
229                     IF( dta%ll_u2d ) dta%u2d(:) = 0._wp
230                     IF( dta%ll_u3d ) dta%v2d(:) = 0._wp
231                  ENDIF
232                  IF (cn_tra(ib_bdy) /= 'runoff') THEN
233                     IF( nn_dyn2d_dta(ib_bdy) == 1 .OR. nn_dyn2d_dta(ib_bdy) == 3 ) THEN
234
235                        jend = jstart + dta%nread(2) - 1
236                        IF( ln_full_vel_array(ib_bdy) ) THEN
237                           CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend),  &
238                                     & kit=jit, kt_offset=time_offset , jpk_bdy=nb_jpk_bdy, fvl=ln_full_vel_array(ib_bdy)  )
239                        ELSE
240                           CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend),  &
241                                     & kit=jit, kt_offset=time_offset  )
242                        ENDIF
243
244                        ! If full velocities in boundary data then extract barotropic velocities from 3D fields
245                        IF( ln_full_vel_array(ib_bdy) .AND.                                             &
246                          &    ( nn_dyn2d_dta(ib_bdy) == 1 .OR. nn_dyn2d_dta(ib_bdy) == 3 .OR.  &
247                          &      nn_dyn3d_dta(ib_bdy) == 1 ) )THEN
248
249                           igrd = 2                      ! zonal velocity
250                           dta%u2d(:) = 0._wp
251                           DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
252                              ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
253                              ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
254                              DO ik = 1, jpkm1
255                                 dta%u2d(ib) = dta%u2d(ib) &
256                       &                          + e3u_n(ii,ij,ik) * umask(ii,ij,ik) * dta%u3d(ib,ik)
257                              END DO
258                              dta%u2d(ib) =  dta%u2d(ib) * r1_hu_n(ii,ij)
259                           END DO
260                           igrd = 3                      ! meridional velocity
261                           dta%v2d(:) = 0._wp
262                           DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
263                              ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
264                              ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
265                              DO ik = 1, jpkm1
266                                 dta%v2d(ib) = dta%v2d(ib) &
267                       &                       + e3v_n(ii,ij,ik) * vmask(ii,ij,ik) * dta%v3d(ib,ik)
268                              END DO
269                              dta%v2d(ib) =  dta%v2d(ib) * r1_hv_n(ii,ij)
270                           END DO
271                        ENDIF                   
272                     ENDIF
273                     IF( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN ! update tidal harmonic forcing
274                        CALL bdytide_update( kt=kt, idx=idx_bdy(ib_bdy), dta=dta, td=tides(ib_bdy),   & 
275                          &                 jit=jit, time_offset=time_offset )
276                     ENDIF
277                  ENDIF
278               ENDIF
279            ELSE
280               IF (cn_tra(ib_bdy) == 'runoff') then      ! runoff condition
281                  jend = nb_bdy_fld(ib_bdy)
282                  CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend),  &
283                               & map=nbmap_ptr(jstart:jend), kt_offset=time_offset )
284                  !
285                  igrd = 2                      ! zonal velocity
286                  DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
287                     ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
288                     ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
289                     dta%u2d(ib) = dta%u2d(ib) / ( e2u(ii,ij) * hu_0(ii,ij) )
290                  END DO
291                  !
292                  igrd = 3                      ! meridional velocity
293                  DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
294                     ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
295                     ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
296                     dta%v2d(ib) = dta%v2d(ib) / ( e1v(ii,ij) * hv_0(ii,ij) )
297                  END DO
298               ELSE
299                  IF( nn_dyn2d_dta(ib_bdy) == 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays
300                     IF( dta%ll_ssh ) dta%ssh(:) = 0._wp
301                     IF( dta%ll_u2d ) dta%u2d(:) = 0._wp
302                     IF( dta%ll_v2d ) dta%v2d(:) = 0._wp
303                  ENDIF
304                  IF( dta%nread(1) .gt. 0 ) THEN ! update external data
305                     jend = jstart + dta%nread(1) - 1
306                     CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), &
307                                  & map=nbmap_ptr(jstart:jend), kt_offset=time_offset, jpk_bdy=nb_jpk_bdy, fvl=ln_full_vel_array(ib_bdy) )
308                  ENDIF
309                  ! If full velocities in boundary data then split into barotropic and baroclinic data
310                  IF( ln_full_vel_array(ib_bdy) .and.                                             &
311                    & ( nn_dyn2d_dta(ib_bdy) == 1 .OR. nn_dyn2d_dta(ib_bdy) == 3 .OR. &
312                    &   nn_dyn3d_dta(ib_bdy) == 1 ) ) THEN
313                     igrd = 2                      ! zonal velocity
314                     dta%u2d(:) = 0._wp
315                     DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
316                        ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
317                        ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
318                        DO ik = 1, jpkm1
319                           dta%u2d(ib) = dta%u2d(ib) &
320                 &                       + e3u_n(ii,ij,ik) * umask(ii,ij,ik) * dta%u3d(ib,ik)
321                        END DO
322                        dta%u2d(ib) =  dta%u2d(ib) * r1_hu_n(ii,ij)
323                        DO ik = 1, jpkm1
324                           dta%u3d(ib,ik) = dta%u3d(ib,ik) - dta%u2d(ib)
325                        END DO
326                     END DO
327                     igrd = 3                      ! meridional velocity
328                     dta%v2d(:) = 0._wp
329                     DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
330                        ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
331                        ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
332                        DO ik = 1, jpkm1
333                           dta%v2d(ib) = dta%v2d(ib) &
334                 &                       + e3v_n(ii,ij,ik) * vmask(ii,ij,ik) * dta%v3d(ib,ik)
335                        END DO
336                        dta%v2d(ib) =  dta%v2d(ib) * r1_hv_n(ii,ij)
337                        DO ik = 1, jpkm1
338                           dta%v3d(ib,ik) = dta%v3d(ib,ik) - dta%v2d(ib)
339                        END DO
340                     END DO
341                  ENDIF
342
343               ENDIF
344#if defined key_lim3
345               IF( .NOT. ll_bdylim3 .AND. cn_ice_lim(ib_bdy) /= 'none' .AND. nn_ice_lim_dta(ib_bdy) == 1 ) THEN ! bdy ice input (case input is 1cat)
346                CALL lim_var_itd ( bf(jfld_hti)%fnow(:,1,1), bf(jfld_hts)%fnow(:,1,1), bf(jfld_ai)%fnow(:,1,1), &
347                                  & dta_bdy(ib_bdy)%ht_i,     dta_bdy(ib_bdy)%ht_s,     dta_bdy(ib_bdy)%a_i     )
348               ENDIF
349#endif
350            ENDIF
351            jstart = jstart + dta%nread(1)
352         END IF ! nn_dta(ib_bdy) = 1
353      END DO  ! ib_bdy
354
355      IF ( ln_tide ) THEN
356         IF (ln_dynspg_ts) THEN      ! Fill temporary arrays with slow-varying bdy data                           
357            DO ib_bdy = 1, nb_bdy    ! Tidal component added in ts loop
358               IF ( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN
359                  nblen => idx_bdy(ib_bdy)%nblen
360                  nblenrim => idx_bdy(ib_bdy)%nblenrim
361                  IF( cn_dyn2d(ib_bdy) == 'frs' ) THEN; ilen1(:)=nblen(:) ; ELSE ; ilen1(:)=nblenrim(:) ; ENDIF
362                  IF ( dta_bdy(ib_bdy)%ll_ssh ) dta_bdy_s(ib_bdy)%ssh(1:ilen1(1)) = dta_bdy(ib_bdy)%ssh(1:ilen1(1))
363                  IF ( dta_bdy(ib_bdy)%ll_u2d ) dta_bdy_s(ib_bdy)%u2d(1:ilen1(2)) = dta_bdy(ib_bdy)%u2d(1:ilen1(2))
364                  IF ( dta_bdy(ib_bdy)%ll_v2d ) dta_bdy_s(ib_bdy)%v2d(1:ilen1(3)) = dta_bdy(ib_bdy)%v2d(1:ilen1(3))
365               ENDIF
366            END DO
367         ELSE ! Add tides if not split-explicit free surface else this is done in ts loop
368            !
369            CALL bdy_dta_tides( kt=kt, time_offset=time_offset )
370         ENDIF
371      ENDIF
372
373      IF ( ln_apr_obc ) THEN
374         DO ib_bdy = 1, nb_bdy
375            IF (cn_tra(ib_bdy) /= 'runoff')THEN
376               igrd = 1                      ! meridional velocity
377               DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
378                  ii   = idx_bdy(ib_bdy)%nbi(ib,igrd)
379                  ij   = idx_bdy(ib_bdy)%nbj(ib,igrd)
380                  dta_bdy(ib_bdy)%ssh(ib) = dta_bdy(ib_bdy)%ssh(ib) + ssh_ib(ii,ij)
381               END DO
382            ENDIF
383         END DO
384      ENDIF
385      !
386      IF( nn_timing == 1 ) CALL timing_stop('bdy_dta')
387      !
388   END SUBROUTINE bdy_dta
389
390
391   SUBROUTINE bdy_dta_init
392      !!----------------------------------------------------------------------
393      !!                   ***  SUBROUTINE bdy_dta_init  ***
394      !!                   
395      !! ** Purpose :   Initialise arrays for reading of external data
396      !!                for open boundary conditions
397      !!
398      !! ** Method  :   
399      !!               
400      !!----------------------------------------------------------------------
401      INTEGER ::   ib_bdy, jfld, jstart, jend, ierror, ios     ! Local integers
402      !
403      CHARACTER(len=100)                     ::   cn_dir        ! Root directory for location of data files
404      CHARACTER(len=100), DIMENSION(nb_bdy)  ::   cn_dir_array  ! Root directory for location of data files
405      CHARACTER(len = 256)::   clname                           ! temporary file name
406      LOGICAL                                ::   ln_full_vel   ! =T => full velocities in 3D boundary data
407                                                                ! =F => baroclinic velocities in 3D boundary data
408      INTEGER                                ::   ilen_global   ! Max length required for global bdy dta arrays
409      INTEGER, ALLOCATABLE, DIMENSION(:)     ::   ilen1, ilen3  ! size of 1st and 3rd dimensions of local arrays
410      INTEGER, ALLOCATABLE, DIMENSION(:)     ::   ibdy           ! bdy set for a particular jfld
411      INTEGER, ALLOCATABLE, DIMENSION(:)     ::   igrid         ! index for grid type (1,2,3 = T,U,V)
412      INTEGER, POINTER, DIMENSION(:)         ::   nblen, nblenrim  ! short cuts
413      TYPE(OBC_DATA), POINTER                ::   dta           ! short cut
414#if defined key_lim3
415      INTEGER               ::   zndims   ! number of dimensions in an array (i.e. 3 = wo ice cat; 4 = w ice cat)
416      INTEGER               ::   inum,id1 ! local integer
417#endif
418      TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) ::   blf_i         !  array of namelist information structures
419      TYPE(FLD_N) ::   bn_tem, bn_sal, bn_u3d, bn_v3d   !
420      TYPE(FLD_N) ::   bn_ssh, bn_u2d, bn_v2d           ! informations about the fields to be read
421#if defined key_lim3
422      TYPE(FLD_N) ::   bn_a_i, bn_ht_i, bn_ht_s     
423      NAMELIST/nambdy_dta/ cn_dir, bn_tem, bn_sal, bn_u3d, bn_v3d, bn_ssh, bn_u2d, bn_v2d 
424      NAMELIST/nambdy_dta/ bn_a_i, bn_ht_i, bn_ht_s
425#endif
426      NAMELIST/nambdy_dta/ ln_full_vel, nb_jpk_bdy
427      !!---------------------------------------------------------------------------
428      !
429      IF( nn_timing == 1 )   CALL timing_start('bdy_dta_init')
430      !
431      IF(lwp) WRITE(numout,*)
432      IF(lwp) WRITE(numout,*) 'bdy_dta_ini : initialization of data at the open boundaries'
433      IF(lwp) WRITE(numout,*) '~~~~~~~~~~'
434      IF(lwp) WRITE(numout,*) ''
435
436      ! Set nn_dta
437      DO ib_bdy = 1, nb_bdy
438         nn_dta(ib_bdy) = MAX(  nn_dyn2d_dta(ib_bdy)       &
439                               ,nn_dyn3d_dta(ib_bdy)       &
440                               ,nn_tra_dta(ib_bdy)         &
441#if defined key_lim3
442                              ,nn_ice_lim_dta(ib_bdy)    &
443#endif
444                              )
445         IF(nn_dta(ib_bdy) > 1) nn_dta(ib_bdy) = 1
446      END DO
447
448      ! Work out upper bound of how many fields there are to read in and allocate arrays
449      ! ---------------------------------------------------------------------------
450      ALLOCATE( nb_bdy_fld(nb_bdy) )
451      nb_bdy_fld(:) = 0
452      DO ib_bdy = 1, nb_bdy         
453         IF( cn_dyn2d(ib_bdy) /= 'none' .and. ( nn_dyn2d_dta(ib_bdy) == 1 .or. nn_dyn2d_dta(ib_bdy) == 3 ) ) THEN
454            nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 3
455         ENDIF
456         IF( cn_dyn3d(ib_bdy) /= 'none' .and. nn_dyn3d_dta(ib_bdy) == 1 ) THEN
457            nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 2
458         ENDIF
459         IF( cn_tra(ib_bdy) /= 'none' .and. nn_tra_dta(ib_bdy) == 1  ) THEN
460            nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 2
461         ENDIF
462#if defined key_lim3
463         IF( cn_ice_lim(ib_bdy) /= 'none' .and. nn_ice_lim_dta(ib_bdy) == 1  ) THEN
464            nb_bdy_fld(ib_bdy) = nb_bdy_fld(ib_bdy) + 3
465         ENDIF
466#endif               
467         IF(lwp) WRITE(numout,*) 'Maximum number of files to open =',nb_bdy_fld(ib_bdy)
468      END DO           
469
470      nb_bdy_fld_sum = SUM( nb_bdy_fld )
471
472      ALLOCATE( bf(nb_bdy_fld_sum), STAT=ierror )
473      IF( ierror > 0 ) THEN   
474         CALL ctl_stop( 'bdy_dta: unable to allocate bf structure' )   ;   RETURN 
475      ENDIF
476      ALLOCATE( blf_i(nb_bdy_fld_sum), STAT=ierror )
477      IF( ierror > 0 ) THEN   
478         CALL ctl_stop( 'bdy_dta: unable to allocate blf_i structure' )   ;   RETURN 
479      ENDIF
480      ALLOCATE( nbmap_ptr(nb_bdy_fld_sum), STAT=ierror )
481      IF( ierror > 0 ) THEN   
482         CALL ctl_stop( 'bdy_dta: unable to allocate nbmap_ptr structure' )   ;   RETURN 
483      ENDIF
484      ALLOCATE( ilen1(nb_bdy_fld_sum), ilen3(nb_bdy_fld_sum) ) 
485      ALLOCATE( ibdy(nb_bdy_fld_sum) ) 
486      ALLOCATE( igrid(nb_bdy_fld_sum) ) 
487
488      ! Read namelists
489      ! --------------
490      REWIND(numnam_ref)
491      REWIND(numnam_cfg)
492      jfld = 0 
493      DO ib_bdy = 1, nb_bdy         
494         IF( nn_dta(ib_bdy) == 1 ) THEN
495            READ  ( numnam_ref, nambdy_dta, IOSTAT = ios, ERR = 901)
496901         IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambdy_dta in reference namelist', lwp )
497
498            READ  ( numnam_cfg, nambdy_dta, IOSTAT = ios, ERR = 902 )
499902         IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambdy_dta in configuration namelist', lwp )
500            IF(lwm) WRITE( numond, nambdy_dta )
501
502            cn_dir_array(ib_bdy) = cn_dir
503            ln_full_vel_array(ib_bdy) = ln_full_vel
504
505            nblen => idx_bdy(ib_bdy)%nblen
506            nblenrim => idx_bdy(ib_bdy)%nblenrim
507            dta => dta_bdy(ib_bdy)
508            dta%nread(2) = 0
509
510            ! Only read in necessary fields for this set.
511            ! Important that barotropic variables come first.
512            IF( nn_dyn2d_dta(ib_bdy) == 1 .or. nn_dyn2d_dta(ib_bdy) == 3 ) THEN
513
514               IF( dta%ll_ssh ) THEN
515                  if(lwp) write(numout,*) '++++++ reading in ssh field'
516                  jfld = jfld + 1
517                  blf_i(jfld) = bn_ssh
518                  ibdy(jfld) = ib_bdy
519                  igrid(jfld) = 1
520                  ilen1(jfld) = nblen(igrid(jfld))
521                  ilen3(jfld) = 1
522                  dta%nread(2) = dta%nread(2) + 1
523               ENDIF
524
525               IF( dta%ll_u2d .and. .not. ln_full_vel_array(ib_bdy) ) THEN
526                  if(lwp) write(numout,*) '++++++ reading in u2d field'
527                  jfld = jfld + 1
528                  blf_i(jfld) = bn_u2d
529                  ibdy(jfld) = ib_bdy
530                  igrid(jfld) = 2
531                  ilen1(jfld) = nblen(igrid(jfld))
532                  ilen3(jfld) = 1
533                  dta%nread(2) = dta%nread(2) + 1
534               ENDIF
535
536               IF( dta%ll_v2d .and. .not. ln_full_vel_array(ib_bdy) ) THEN
537                  if(lwp) write(numout,*) '++++++ reading in v2d field'
538                  jfld = jfld + 1
539                  blf_i(jfld) = bn_v2d
540                  ibdy(jfld) = ib_bdy
541                  igrid(jfld) = 3
542                  ilen1(jfld) = nblen(igrid(jfld))
543                  ilen3(jfld) = 1
544                  dta%nread(2) = dta%nread(2) + 1
545               ENDIF
546
547            ENDIF
548
549            ! read 3D velocities if baroclinic velocities require OR if
550            ! barotropic velocities required and ln_full_vel set to .true.
551            IF( nn_dyn3d_dta(ib_bdy) == 1 .OR. &
552           &  ( ln_full_vel_array(ib_bdy) .AND. ( nn_dyn2d_dta(ib_bdy) == 1 .or. nn_dyn2d_dta(ib_bdy) == 3 ) ) ) THEN
553
554               IF( dta%ll_u3d .OR. ( ln_full_vel_array(ib_bdy) .and. dta%ll_u2d ) ) THEN
555                  if(lwp) write(numout,*) '++++++ reading in u3d field'
556                  jfld = jfld + 1
557                  blf_i(jfld) = bn_u3d
558                  ibdy(jfld) = ib_bdy
559                  igrid(jfld) = 2
560                  ilen1(jfld) = nblen(igrid(jfld))
561                  ilen3(jfld) = jpk
562                  IF( ln_full_vel_array(ib_bdy) .and. dta%ll_u2d ) dta%nread(2) = dta%nread(2) + 1
563               ENDIF
564
565               IF( dta%ll_v3d .OR. ( ln_full_vel_array(ib_bdy) .and. dta%ll_v2d ) ) THEN
566                  if(lwp) write(numout,*) '++++++ reading in v3d field'
567                  jfld = jfld + 1
568                  blf_i(jfld) = bn_v3d
569                  ibdy(jfld) = ib_bdy
570                  igrid(jfld) = 3
571                  ilen1(jfld) = nblen(igrid(jfld))
572                  ilen3(jfld) = jpk
573                  IF( ln_full_vel_array(ib_bdy) .and. dta%ll_v2d ) dta%nread(2) = dta%nread(2) + 1
574               ENDIF
575
576            ENDIF
577
578            ! temperature and salinity
579            IF( nn_tra_dta(ib_bdy) == 1 ) THEN
580
581               IF( dta%ll_tem ) THEN
582                  if(lwp) write(numout,*) '++++++ reading in tem field'
583                  jfld = jfld + 1
584                  blf_i(jfld) = bn_tem
585                  ibdy(jfld) = ib_bdy
586                  igrid(jfld) = 1
587                  ilen1(jfld) = nblen(igrid(jfld))
588                  ilen3(jfld) = jpk
589               ENDIF
590
591               IF( dta%ll_sal ) THEN
592                  if(lwp) write(numout,*) '++++++ reading in sal field'
593                  jfld = jfld + 1
594                  blf_i(jfld) = bn_sal
595                  ibdy(jfld) = ib_bdy
596                  igrid(jfld) = 1
597                  ilen1(jfld) = nblen(igrid(jfld))
598                  ilen3(jfld) = jpk
599               ENDIF
600
601            ENDIF
602
603#if defined key_lim3
604            ! sea ice
605            IF( nn_ice_lim_dta(ib_bdy) == 1 ) THEN
606               ! Test for types of ice input (1cat or Xcat)
607               ! Build file name to find dimensions
608               clname=TRIM( cn_dir )//TRIM(bn_a_i%clname)
609               IF( .NOT. bn_a_i%ln_clim ) THEN   
610                                                  WRITE(clname, '(a,"_y",i4.4)' ) TRIM( clname ), nyear    ! add year
611                  IF( bn_a_i%cltype /= 'yearly' ) WRITE(clname, '(a,"m" ,i2.2)' ) TRIM( clname ), nmonth   ! add month
612               ELSE
613                  IF( bn_a_i%cltype /= 'yearly' ) WRITE(clname, '(a,"_m",i2.2)' ) TRIM( clname ), nmonth   ! add month
614               ENDIF
615               IF( bn_a_i%cltype == 'daily' .OR. bn_a_i%cltype(1:4) == 'week' ) &
616               &                                  WRITE(clname, '(a,"d" ,i2.2)' ) TRIM( clname ), nday     ! add day
617               !
618               CALL iom_open  ( clname, inum )
619               id1 = iom_varid( inum, bn_a_i%clvar, kndims=zndims, ldstop = .FALSE. )
620               CALL iom_close ( inum )
621
622                IF ( zndims == 4 ) THEN
623                 ll_bdylim3 = .TRUE.   ! Xcat input
624               ELSE
625                 ll_bdylim3 = .FALSE.  ! 1cat input     
626               ENDIF
627               ! End test
628
629               IF( dta%ll_a_i ) THEN
630                  jfld = jfld + 1
631                  blf_i(jfld) = bn_a_i
632                  ibdy(jfld) = ib_bdy
633                  igrid(jfld) = 1
634                  ilen1(jfld) = nblen(igrid(jfld))
635                  IF ( ll_bdylim3 ) THEN ; ilen3(jfld)=jpl ; ELSE ; ilen3(jfld)=1 ; ENDIF
636               ENDIF
637
638               IF( dta%ll_ht_i ) THEN
639                  jfld = jfld + 1
640                  blf_i(jfld) = bn_ht_i
641                  ibdy(jfld) = ib_bdy
642                  igrid(jfld) = 1
643                  ilen1(jfld) = nblen(igrid(jfld))
644                  IF ( ll_bdylim3 ) THEN ; ilen3(jfld)=jpl ; ELSE ; ilen3(jfld)=1 ; ENDIF
645               ENDIF
646
647               IF( dta%ll_ht_s ) THEN
648                  jfld = jfld + 1
649                   blf_i(jfld) = bn_ht_s
650                  ibdy(jfld) = ib_bdy
651                  igrid(jfld) = 1
652                  ilen1(jfld) = nblen(igrid(jfld))
653                  IF ( ll_bdylim3 ) THEN ; ilen3(jfld)=jpl ; ELSE ; ilen3(jfld)=1 ; ENDIF
654               ENDIF
655
656            ENDIF
657#endif
658            ! Recalculate field counts
659            !-------------------------
660            IF( ib_bdy == 1 ) THEN
661               nb_bdy_fld_sum = 0
662               nb_bdy_fld(ib_bdy) = jfld
663               nb_bdy_fld_sum     = jfld             
664            ELSE
665               nb_bdy_fld(ib_bdy) = jfld - nb_bdy_fld_sum
666               nb_bdy_fld_sum = nb_bdy_fld_sum + nb_bdy_fld(ib_bdy)
667            ENDIF
668
669            dta%nread(1) = nb_bdy_fld(ib_bdy)
670
671         ENDIF ! nn_dta == 1
672      ENDDO ! ib_bdy
673
674      DO jfld = 1, nb_bdy_fld_sum
675         ALLOCATE( bf(jfld)%fnow(ilen1(jfld),1,ilen3(jfld)) )
676         IF( blf_i(jfld)%ln_tint ) ALLOCATE( bf(jfld)%fdta(ilen1(jfld),1,ilen3(jfld),2) )
677         nbmap_ptr(jfld)%ptr => idx_bdy(ibdy(jfld))%nbmap(:,igrid(jfld))
678         nbmap_ptr(jfld)%ll_unstruc = ln_coords_file(ibdy(jfld))
679      ENDDO
680
681      ! fill bf with blf_i and control print
682      !-------------------------------------
683      jstart = 1
684      DO ib_bdy = 1, nb_bdy
685         jend = jstart - 1 + nb_bdy_fld(ib_bdy) 
686         CALL fld_fill( bf(jstart:jend), blf_i(jstart:jend), cn_dir_array(ib_bdy), 'bdy_dta',   &
687         &              'open boundary conditions', 'nambdy_dta' )
688         jstart = jend + 1
689      ENDDO
690
691      ! Initialise local boundary data arrays
692      ! nn_xxx_dta=0 : allocate space - will be filled from initial conditions later
693      ! nn_xxx_dta=1 : point to "fnow" arrays
694      !-------------------------------------
695
696      jfld = 0
697      DO ib_bdy=1, nb_bdy
698
699         nblen => idx_bdy(ib_bdy)%nblen
700         dta => dta_bdy(ib_bdy)
701
702         if(lwp) then
703            write(numout,*) '++++++ dta%ll_ssh = ',dta%ll_ssh
704            write(numout,*) '++++++ dta%ll_u2d = ',dta%ll_u2d
705            write(numout,*) '++++++ dta%ll_v2d = ',dta%ll_v2d
706            write(numout,*) '++++++ dta%ll_u3d = ',dta%ll_u3d
707            write(numout,*) '++++++ dta%ll_v3d = ',dta%ll_v3d
708            write(numout,*) '++++++ dta%ll_tem = ',dta%ll_tem
709            write(numout,*) '++++++ dta%ll_sal = ',dta%ll_sal
710         endif
711
712         IF ( nn_dyn2d_dta(ib_bdy) == 0 .or. nn_dyn2d_dta(ib_bdy) == 2 ) THEN
713            if(lwp) write(numout,*) '++++++ dta%ssh/u2d/u3d allocated space'
714            IF( dta%ll_ssh ) ALLOCATE( dta%ssh(nblen(1)) )
715            IF( dta%ll_u2d ) ALLOCATE( dta%u2d(nblen(2)) )
716            IF( dta%ll_v2d ) ALLOCATE( dta%v2d(nblen(3)) )
717         ENDIF
718         IF ( nn_dyn2d_dta(ib_bdy) == 1 .or. nn_dyn2d_dta(ib_bdy) == 3 ) THEN
719            IF( dta%ll_ssh ) THEN
720               if(lwp) write(numout,*) '++++++ dta%ssh pointing to fnow'
721               jfld = jfld + 1
722               dta%ssh => bf(jfld)%fnow(:,1,1)
723            ENDIF
724            IF ( dta%ll_u2d ) THEN
725               IF ( ln_full_vel_array(ib_bdy) ) THEN
726                  if(lwp) write(numout,*) '++++++ dta%u2d allocated space'
727                  ALLOCATE( dta%u2d(nblen(2)) )
728               ELSE
729                  if(lwp) write(numout,*) '++++++ dta%u2d pointing to fnow'
730                  jfld = jfld + 1
731                  dta%u2d => bf(jfld)%fnow(:,1,1)
732               ENDIF
733            ENDIF
734            IF ( dta%ll_v2d ) THEN
735               IF ( ln_full_vel_array(ib_bdy) ) THEN
736                  if(lwp) write(numout,*) '++++++ dta%v2d allocated space'
737                  ALLOCATE( dta%v2d(nblen(3)) )
738               ELSE
739                  if(lwp) write(numout,*) '++++++ dta%v2d pointing to fnow'
740                  jfld = jfld + 1
741                  dta%v2d => bf(jfld)%fnow(:,1,1)
742               ENDIF
743            ENDIF
744         ENDIF
745
746         IF ( nn_dyn3d_dta(ib_bdy) == 0 ) THEN
747            if(lwp) write(numout,*) '++++++ dta%u3d/v3d allocated space'
748            IF( dta%ll_u3d ) ALLOCATE( dta_bdy(ib_bdy)%u3d(nblen(2),jpk) )
749            IF( dta%ll_v3d ) ALLOCATE( dta_bdy(ib_bdy)%v3d(nblen(3),jpk) )
750         ENDIF
751         IF ( nn_dyn3d_dta(ib_bdy) == 1 .or. &
752           &  ( ln_full_vel_array(ib_bdy) .and. ( nn_dyn2d_dta(ib_bdy) == 1 .or. nn_dyn2d_dta(ib_bdy) == 3 ) ) ) THEN
753            IF ( dta%ll_u3d .or. ( ln_full_vel_array(ib_bdy) .and. dta%ll_u2d ) ) THEN
754               if(lwp) write(numout,*) '++++++ dta%u3d pointing to fnow'
755               jfld = jfld + 1
756               dta_bdy(ib_bdy)%u3d => bf(jfld)%fnow(:,1,:)
757            ENDIF
758            IF ( dta%ll_v3d .or. ( ln_full_vel_array(ib_bdy) .and. dta%ll_v2d ) ) THEN
759               if(lwp) write(numout,*) '++++++ dta%v3d pointing to fnow'
760               jfld = jfld + 1
761               dta_bdy(ib_bdy)%v3d => bf(jfld)%fnow(:,1,:)
762            ENDIF
763         ENDIF
764
765         IF( nn_tra_dta(ib_bdy) == 0 ) THEN
766            if(lwp) write(numout,*) '++++++ dta%tem/sal allocated space'
767            IF( dta%ll_tem ) ALLOCATE( dta_bdy(ib_bdy)%tem(nblen(1),jpk) )
768            IF( dta%ll_sal ) ALLOCATE( dta_bdy(ib_bdy)%sal(nblen(1),jpk) )
769         ELSE
770            IF( dta%ll_tem ) THEN
771               if(lwp) write(numout,*) '++++++ dta%tem pointing to fnow'
772               jfld = jfld + 1
773               dta_bdy(ib_bdy)%tem => bf(jfld)%fnow(:,1,:)
774            ENDIF
775            IF( dta%ll_sal ) THEN
776               if(lwp) write(numout,*) '++++++ dta%sal pointing to fnow'
777               jfld = jfld + 1
778               dta_bdy(ib_bdy)%sal => bf(jfld)%fnow(:,1,:)
779            ENDIF
780         ENDIF
781
782#if defined key_lim3
783         IF (cn_ice_lim(ib_bdy) /= 'none') THEN
784            IF( nn_ice_lim_dta(ib_bdy) == 0 ) THEN
785               ALLOCATE( dta_bdy(ib_bdy)%a_i (nblen(1),jpl) )
786               ALLOCATE( dta_bdy(ib_bdy)%ht_i(nblen(1),jpl) )
787               ALLOCATE( dta_bdy(ib_bdy)%ht_s(nblen(1),jpl) )
788            ELSE
789               IF ( ll_bdylim3 ) THEN ! case input is Xcat
790                  jfld = jfld + 1
791                  dta_bdy(ib_bdy)%a_i  => bf(jfld)%fnow(:,1,:)
792                  jfld = jfld + 1
793                  dta_bdy(ib_bdy)%ht_i => bf(jfld)%fnow(:,1,:)
794                  jfld = jfld + 1
795                  dta_bdy(ib_bdy)%ht_s => bf(jfld)%fnow(:,1,:)
796               ELSE ! case input is 1cat
797                  jfld_ai  = jfld + 1
798                  jfld_hti = jfld + 2
799                  jfld_hts = jfld + 3
800                  jfld     = jfld + 3
801                  ALLOCATE( dta_bdy(ib_bdy)%a_i (nblen(1),jpl) )
802                  ALLOCATE( dta_bdy(ib_bdy)%ht_i(nblen(1),jpl) )
803                  ALLOCATE( dta_bdy(ib_bdy)%ht_s(nblen(1),jpl) )
804                  dta_bdy(ib_bdy)%a_i (:,:) = 0._wp
805                  dta_bdy(ib_bdy)%ht_i(:,:) = 0._wp
806                  dta_bdy(ib_bdy)%ht_s(:,:) = 0._wp
807               ENDIF
808
809            ENDIF
810         ENDIF
811#endif
812         !
813      END DO ! ib_bdy
814      !
815      IF( nn_timing == 1 ) CALL timing_stop('bdy_dta_init')
816      !
817   END SUBROUTINE bdy_dta_init
818
819   !!==============================================================================
820END MODULE bdydta
Note: See TracBrowser for help on using the repository browser.