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.
bdytides.F90 in NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE/BDY – NEMO

source: NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE/BDY/bdytides.F90 @ 11574

Last change on this file since 11574 was 10865, checked in by smueller, 5 years ago

Removal of redundant source code in modules bdydta and bdytides (ticket #2194)
  • Property svn:keywords set to Id
File size: 23.1 KB
Line 
1MODULE bdytides
2   !!======================================================================
3   !!                       ***  MODULE  bdytides  ***
4   !! Ocean dynamics:   Tidal forcing at open boundaries
5   !!======================================================================
6   !! History :  2.0  !  2007-01  (D.Storkey)  Original code
7   !!            2.3  !  2008-01  (J.Holt)  Add date correction. Origins POLCOMS v6.3 2007
8   !!            3.0  !  2008-04  (NEMO team)  add in the reference version
9   !!            3.3  !  2010-09  (D.Storkey and E.O'Dea)  bug fixes
10   !!            3.4  !  2012-09  (G. Reffray and J. Chanut) New inputs + mods
11   !!            3.5  !  2013-07  (J. Chanut) Compliant with time splitting changes
12   !!----------------------------------------------------------------------
13   !!   bdytide_init  : read of namelist and initialisation of tidal harmonics data
14   !!   tide_update   : calculation of tidal forcing at each timestep
15   !!----------------------------------------------------------------------
16   USE oce            ! ocean dynamics and tracers
17   USE dom_oce        ! ocean space and time domain
18   USE phycst         ! physical constants
19   USE bdy_oce        ! ocean open boundary conditions
20   USE tide_mod       !
21   USE daymod         ! calendar
22   !
23   USE in_out_manager ! I/O units
24   USE iom            ! xIO server
25   USE fldread        !
26   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
27
28   IMPLICIT NONE
29   PRIVATE
30
31   PUBLIC   bdytide_init     ! routine called in bdy_init
32   PUBLIC   bdy_dta_tides    ! routine called in dyn_spg_ts
33
34   TYPE, PUBLIC ::   TIDES_DATA     !: Storage for external tidal harmonics data
35      REAL(wp), POINTER, DIMENSION(:,:,:)    ::   ssh0     !: Tidal constituents : SSH0   (read in file)
36      REAL(wp), POINTER, DIMENSION(:,:,:)    ::   u0, v0   !: Tidal constituents : U0, V0 (read in file)
37      REAL(wp), POINTER, DIMENSION(:,:,:)    ::   ssh      !: Tidal constituents : SSH    (after nodal cor.)
38      REAL(wp), POINTER, DIMENSION(:,:,:)    ::   u , v    !: Tidal constituents : U , V  (after nodal cor.)
39   END TYPE TIDES_DATA
40
41!$AGRIF_DO_NOT_TREAT
42   TYPE(TIDES_DATA), PUBLIC, DIMENSION(jp_bdy), TARGET :: tides  !: External tidal harmonics data
43!$AGRIF_END_DO_NOT_TREAT
44   TYPE(OBC_DATA)  , PUBLIC, DIMENSION(jp_bdy) :: dta_bdy_s  !: bdy external data (slow component)
45
46   INTEGER ::   kt_tide
47
48   !!----------------------------------------------------------------------
49   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
50   !! $Id$
51   !! Software governed by the CeCILL license (see ./LICENSE)
52   !!----------------------------------------------------------------------
53CONTAINS
54
55   SUBROUTINE bdytide_init
56      !!----------------------------------------------------------------------
57      !!                    ***  SUBROUTINE bdytide_init  ***
58      !!                     
59      !! ** Purpose : - Read in namelist for tides and initialise external
60      !!                tidal harmonics data
61      !!
62      !!----------------------------------------------------------------------
63      !! namelist variables
64      !!-------------------
65      CHARACTER(len=80)                         ::   filtide             !: Filename root for tidal input files
66      LOGICAL                                   ::   ln_bdytide_2ddta    !: If true, read 2d harmonic data
67      !!
68      INTEGER                                   ::   ib_bdy, itide, ib   !: dummy loop indices
69      INTEGER                                   ::   ii, ij              !: dummy loop indices
70      INTEGER                                   ::   inum, igrd
71      INTEGER, DIMENSION(3)                     ::   ilen0       !: length of boundary data (from OBC arrays)
72      INTEGER, POINTER, DIMENSION(:)            ::   nblen, nblenrim     ! short cuts
73      INTEGER                                   ::   ios                 ! Local integer output status for namelist read
74      CHARACTER(len=80)                         ::   clfile              !: full file name for tidal input file
75      REAL(wp),ALLOCATABLE, DIMENSION(:,:,:)    ::   dta_read            !: work space to read in tidal harmonics data
76      REAL(wp),ALLOCATABLE, DIMENSION(:,:)      ::   ztr, zti            !:  "     "    "   "   "   "        "      "
77      !!
78      TYPE(TIDES_DATA),  POINTER                ::   td                  !: local short cut   
79      TYPE(MAP_POINTER), DIMENSION(jpbgrd)      ::   ibmap_ptr           !: array of pointers to nbmap
80      !!
81      NAMELIST/nambdy_tide/filtide, ln_bdytide_2ddta
82      !!----------------------------------------------------------------------
83      !
84      IF (nb_bdy>0) THEN
85         IF(lwp) WRITE(numout,*)
86         IF(lwp) WRITE(numout,*) 'bdytide_init : initialization of tidal harmonic forcing at open boundaries'
87         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
88      ENDIF
89
90      REWIND(numnam_cfg)
91
92      DO ib_bdy = 1, nb_bdy
93         IF( nn_dyn2d_dta(ib_bdy) >= 2 ) THEN
94            !
95            td => tides(ib_bdy)
96            nblen => idx_bdy(ib_bdy)%nblen
97            nblenrim => idx_bdy(ib_bdy)%nblenrim
98
99            ! Namelist nambdy_tide : tidal harmonic forcing at open boundaries
100            filtide(:) = ''
101
102            ! Don't REWIND here - may need to read more than one of these namelists.
103            READ  ( numnam_ref, nambdy_tide, IOSTAT = ios, ERR = 901)
104901         IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambdy_tide in reference namelist', lwp )
105            READ  ( numnam_cfg, nambdy_tide, IOSTAT = ios, ERR = 902 )
106902         IF( ios >  0 )   CALL ctl_nam ( ios , 'nambdy_tide in configuration namelist', lwp )
107            IF(lwm) WRITE ( numond, nambdy_tide )
108            !                                               ! Parameter control and print
109            IF(lwp) WRITE(numout,*) '  '
110            IF(lwp) WRITE(numout,*) '          Namelist nambdy_tide : tidal harmonic forcing at open boundaries'
111            IF(lwp) WRITE(numout,*) '             read tidal data in 2d files: ', ln_bdytide_2ddta
112            IF(lwp) WRITE(numout,*) '             Number of tidal components to read: ', nb_harmo
113            IF(lwp) THEN
114                    WRITE(numout,*) '             Tidal components: ' 
115               DO itide = 1, nb_harmo
116                  WRITE(numout,*)  '                 ', tide_harmonics(itide)%cname_tide 
117               END DO
118            ENDIF
119            IF(lwp) WRITE(numout,*) ' '
120
121            ! Allocate space for tidal harmonics data - get size from OBC data arrays
122            ! -----------------------------------------------------------------------
123
124            ! JC: If FRS scheme is used, we assume that tidal is needed over the whole
125            ! relaxation area     
126            IF( cn_dyn2d(ib_bdy) == 'frs' ) THEN   ;   ilen0(:) = nblen   (:)
127            ELSE                                   ;   ilen0(:) = nblenrim(:)
128            ENDIF
129
130            ALLOCATE( td%ssh0( ilen0(1), nb_harmo, 2 ) )
131            ALLOCATE( td%ssh ( ilen0(1), nb_harmo, 2 ) )
132
133            ALLOCATE( td%u0( ilen0(2), nb_harmo, 2 ) )
134            ALLOCATE( td%u ( ilen0(2), nb_harmo, 2 ) )
135
136            ALLOCATE( td%v0( ilen0(3), nb_harmo, 2 ) )
137            ALLOCATE( td%v ( ilen0(3), nb_harmo, 2 ) )
138
139            td%ssh0(:,:,:) = 0._wp
140            td%ssh (:,:,:) = 0._wp
141            td%u0  (:,:,:) = 0._wp
142            td%u   (:,:,:) = 0._wp
143            td%v0  (:,:,:) = 0._wp
144            td%v   (:,:,:) = 0._wp
145
146            IF( ln_bdytide_2ddta ) THEN
147               ! It is assumed that each data file contains all complex harmonic amplitudes
148               ! given on the global domain (ie global, jpiglo x jpjglo)
149               !
150               ALLOCATE( zti(jpi,jpj), ztr(jpi,jpj) )
151               !
152               ! SSH fields
153               clfile = TRIM(filtide)//'_grid_T.nc'
154               CALL iom_open( clfile , inum ) 
155               igrd = 1                       ! Everything is at T-points here
156               DO itide = 1, nb_harmo
157                  CALL iom_get( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_z1', ztr(:,:) )
158                  CALL iom_get( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_z2', zti(:,:) ) 
159                  DO ib = 1, ilen0(igrd)
160                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
161                     ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
162                     td%ssh0(ib,itide,1) = ztr(ii,ij)
163                     td%ssh0(ib,itide,2) = zti(ii,ij)
164                  END DO
165               END DO
166               CALL iom_close( inum )
167               !
168               ! U fields
169               clfile = TRIM(filtide)//'_grid_U.nc'
170               CALL iom_open( clfile , inum ) 
171               igrd = 2                       ! Everything is at U-points here
172               DO itide = 1, nb_harmo
173                  CALL iom_get  ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_u1', ztr(:,:) )
174                  CALL iom_get  ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_u2', zti(:,:) )
175                  DO ib = 1, ilen0(igrd)
176                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
177                     ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
178                     td%u0(ib,itide,1) = ztr(ii,ij)
179                     td%u0(ib,itide,2) = zti(ii,ij)
180                  END DO
181               END DO
182               CALL iom_close( inum )
183               !
184               ! V fields
185               clfile = TRIM(filtide)//'_grid_V.nc'
186               CALL iom_open( clfile , inum ) 
187               igrd = 3                       ! Everything is at V-points here
188               DO itide = 1, nb_harmo
189                  CALL iom_get  ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_v1', ztr(:,:) )
190                  CALL iom_get  ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_v2', zti(:,:) )
191                  DO ib = 1, ilen0(igrd)
192                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
193                     ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
194                     td%v0(ib,itide,1) = ztr(ii,ij)
195                     td%v0(ib,itide,2) = zti(ii,ij)
196                  END DO
197               END DO 
198               CALL iom_close( inum )
199               !
200               DEALLOCATE( ztr, zti ) 
201               !
202            ELSE           
203               !
204               ! Read tidal data only on bdy segments
205               !
206               ALLOCATE( dta_read( MAXVAL(ilen0(1:3)), 1, 1 ) )
207               !
208               ! Set map structure
209               ibmap_ptr(1)%ptr => idx_bdy(ib_bdy)%nbmap(:,1)   ;   ibmap_ptr(1)%ll_unstruc = ln_coords_file(ib_bdy)
210               ibmap_ptr(2)%ptr => idx_bdy(ib_bdy)%nbmap(:,2)   ;   ibmap_ptr(2)%ll_unstruc = ln_coords_file(ib_bdy)
211               ibmap_ptr(3)%ptr => idx_bdy(ib_bdy)%nbmap(:,3)   ;   ibmap_ptr(3)%ll_unstruc = ln_coords_file(ib_bdy)
212
213               ! Open files and read in tidal forcing data
214               ! -----------------------------------------
215
216               DO itide = 1, nb_harmo
217                  !                                                              ! SSH fields
218                  clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_T.nc'
219                  CALL iom_open( clfile, inum )
220                  CALL fld_map( inum, 'z1' , dta_read(1:ilen0(1),1:1,1:1) , 1,  ibmap_ptr(1) )
221                  td%ssh0(:,itide,1) = dta_read(1:ilen0(1),1,1)
222                  CALL fld_map( inum, 'z2' , dta_read(1:ilen0(1),1:1,1:1) , 1,  ibmap_ptr(1) )
223                  td%ssh0(:,itide,2) = dta_read(1:ilen0(1),1,1)
224                  CALL iom_close( inum )
225                  !                                                              ! U fields
226                  clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_U.nc'
227                  CALL iom_open( clfile, inum )
228                  CALL fld_map( inum, 'u1' , dta_read(1:ilen0(2),1:1,1:1) , 1, ibmap_ptr(2) )
229                  td%u0(:,itide,1) = dta_read(1:ilen0(2),1,1)
230                  CALL fld_map( inum, 'u2' , dta_read(1:ilen0(2),1:1,1:1) , 1, ibmap_ptr(2) )
231                  td%u0(:,itide,2) = dta_read(1:ilen0(2),1,1)
232                  CALL iom_close( inum )
233                  !                                                              ! V fields
234                  clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_V.nc'
235                  CALL iom_open( clfile, inum )
236                  CALL fld_map( inum, 'v1' , dta_read(1:ilen0(3),1:1,1:1) , 1, ibmap_ptr(3) )
237                  td%v0(:,itide,1) = dta_read(1:ilen0(3),1,1)
238                  CALL fld_map( inum, 'v2' , dta_read(1:ilen0(3),1:1,1:1) , 1, ibmap_ptr(3) )
239                  td%v0(:,itide,2) = dta_read(1:ilen0(3),1,1)
240                  CALL iom_close( inum )
241                  !
242               END DO ! end loop on tidal components
243               !
244               DEALLOCATE( dta_read )
245               !
246            ENDIF ! ln_bdytide_2ddta=.true.
247            !
248            ! Allocate slow varying data in the case of time splitting:
249            ! Do it anyway because at this stage knowledge of free surface scheme is unknown
250            ALLOCATE( dta_bdy_s(ib_bdy)%ssh ( ilen0(1) ) )
251            ALLOCATE( dta_bdy_s(ib_bdy)%u2d ( ilen0(2) ) )
252            ALLOCATE( dta_bdy_s(ib_bdy)%v2d ( ilen0(3) ) )
253            dta_bdy_s(ib_bdy)%ssh(:) = 0._wp
254            dta_bdy_s(ib_bdy)%u2d(:) = 0._wp
255            dta_bdy_s(ib_bdy)%v2d(:) = 0._wp
256            !
257         ENDIF ! nn_dyn2d_dta(ib_bdy) >= 2
258         !
259      END DO ! loop on ib_bdy
260      !
261   END SUBROUTINE bdytide_init
262
263
264   SUBROUTINE bdy_dta_tides( kt, kit, time_offset )
265      !!----------------------------------------------------------------------
266      !!                 ***  SUBROUTINE bdy_dta_tides  ***
267      !!               
268      !! ** Purpose : - Add tidal forcing to ssh, u2d and v2d OBC data arrays.
269      !!               
270      !!----------------------------------------------------------------------
271      INTEGER,           INTENT(in) ::   kt          ! Main timestep counter
272      INTEGER, OPTIONAL, INTENT(in) ::   kit         ! Barotropic timestep counter (for timesplitting option)
273      INTEGER, OPTIONAL, INTENT(in) ::   time_offset ! time offset in units of timesteps. NB. if kit
274      !                                              ! is present then units = subcycle timesteps.
275      !                                              ! time_offset = 0  => get data at "now"    time level
276      !                                              ! time_offset = -1 => get data at "before" time level
277      !                                              ! time_offset = +1 => get data at "after"  time level
278      !                                              ! etc.
279      !
280      LOGICAL  ::   lk_first_btstp            ! =.TRUE. if time splitting and first barotropic step
281      INTEGER  ::   itide, ib_bdy, ib, igrd   ! loop indices
282      INTEGER  ::   time_add                  ! time offset in units of timesteps
283      INTEGER, DIMENSION(jpbgrd)   ::   ilen0 
284      INTEGER, DIMENSION(1:jpbgrd) ::   nblen, nblenrim  ! short cuts
285      REAL(wp) ::   z_arg, z_sarg, zramp, zoff, z_cost, z_sist     
286      !!----------------------------------------------------------------------
287      !
288      lk_first_btstp=.TRUE.
289      IF ( PRESENT(kit).AND.( kit /= 1 ) ) THEN ; lk_first_btstp=.FALSE. ; ENDIF
290
291      time_add = 0
292      IF( PRESENT(time_offset) ) THEN
293         time_add = time_offset
294      ENDIF
295     
296      ! Absolute time from model initialization:   
297      IF( PRESENT(kit) ) THEN 
298         z_arg = ( kt + (kit+time_add-1) / REAL(nn_baro,wp) ) * rdt
299      ELSE                             
300         z_arg = ( kt + time_add ) * rdt
301      ENDIF
302
303      ! Linear ramp on tidal component at open boundaries
304      zramp = 1.
305      IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg - nit000*rdt)/(rn_tide_ramp_dt*rday),0.),1.)
306
307      DO ib_bdy = 1,nb_bdy
308         !
309         IF( nn_dyn2d_dta(ib_bdy) >= 2 ) THEN
310            !
311            nblen(1:jpbgrd) = idx_bdy(ib_bdy)%nblen(1:jpbgrd)
312            nblenrim(1:jpbgrd) = idx_bdy(ib_bdy)%nblenrim(1:jpbgrd)
313            !
314            IF( cn_dyn2d(ib_bdy) == 'frs' ) THEN   ;   ilen0(:) = nblen   (:)
315            ELSE                                   ;   ilen0(:) = nblenrim(:)
316            ENDIF     
317            !
318            ! We refresh nodal factors every day below
319            ! This should be done somewhere else
320            IF ( ( nsec_day == NINT(0.5_wp * rdt) .OR. kt==nit000 ) .AND. lk_first_btstp ) THEN
321               !
322               kt_tide = kt - (nsec_day - 0.5_wp * rdt)/rdt
323               !
324               IF(lwp) THEN
325               WRITE(numout,*)
326               WRITE(numout,*) 'bdy_tide_dta : Refresh nodal factors for tidal open bdy data at kt=',kt
327               WRITE(numout,*) '~~~~~~~~~~~~~~ '
328               ENDIF
329               !
330               CALL tide_init_elevation ( idx=idx_bdy(ib_bdy), td=tides(ib_bdy) )
331               CALL tide_init_velocities( idx=idx_bdy(ib_bdy), td=tides(ib_bdy) )
332               !
333            ENDIF
334            zoff = -kt_tide * rdt ! time offset relative to nodal factor computation time
335            !
336            ! If time splitting, initialize arrays from slow varying open boundary data:
337            IF ( PRESENT(kit) ) THEN           
338               IF ( dta_bdy(ib_bdy)%ll_ssh ) dta_bdy(ib_bdy)%ssh(1:ilen0(1)) = dta_bdy_s(ib_bdy)%ssh(1:ilen0(1))
339               IF ( dta_bdy(ib_bdy)%ll_u2d ) dta_bdy(ib_bdy)%u2d(1:ilen0(2)) = dta_bdy_s(ib_bdy)%u2d(1:ilen0(2))
340               IF ( dta_bdy(ib_bdy)%ll_v2d ) dta_bdy(ib_bdy)%v2d(1:ilen0(3)) = dta_bdy_s(ib_bdy)%v2d(1:ilen0(3))
341            ENDIF
342            !
343            ! Update open boundary data arrays:
344            DO itide = 1, nb_harmo
345               !
346               z_sarg = (z_arg + zoff) * tide_harmonics(itide)%omega
347               z_cost = zramp * COS( z_sarg )
348               z_sist = zramp * SIN( z_sarg )
349               !
350               IF ( dta_bdy(ib_bdy)%ll_ssh ) THEN
351                  igrd=1                              ! SSH on tracer grid
352                  DO ib = 1, ilen0(igrd)
353                     dta_bdy(ib_bdy)%ssh(ib) = dta_bdy(ib_bdy)%ssh(ib) + &
354                        &                      ( tides(ib_bdy)%ssh(ib,itide,1)*z_cost + &
355                        &                        tides(ib_bdy)%ssh(ib,itide,2)*z_sist )
356                  END DO
357               ENDIF
358               !
359               IF ( dta_bdy(ib_bdy)%ll_u2d ) THEN
360                  igrd=2                              ! U grid
361                  DO ib = 1, ilen0(igrd)
362                     dta_bdy(ib_bdy)%u2d(ib) = dta_bdy(ib_bdy)%u2d(ib) + &
363                        &                      ( tides(ib_bdy)%u(ib,itide,1)*z_cost + &
364                        &                        tides(ib_bdy)%u(ib,itide,2)*z_sist )
365                  END DO
366               ENDIF
367               !
368               IF ( dta_bdy(ib_bdy)%ll_v2d ) THEN
369                  igrd=3                              ! V grid
370                  DO ib = 1, ilen0(igrd) 
371                     dta_bdy(ib_bdy)%v2d(ib) = dta_bdy(ib_bdy)%v2d(ib) + &
372                        &                      ( tides(ib_bdy)%v(ib,itide,1)*z_cost + &
373                        &                        tides(ib_bdy)%v(ib,itide,2)*z_sist )
374                  END DO
375               ENDIF
376            END DO             
377         END IF
378      END DO
379      !
380   END SUBROUTINE bdy_dta_tides
381
382
383   SUBROUTINE tide_init_elevation( idx, td )
384      !!----------------------------------------------------------------------
385      !!                 ***  ROUTINE tide_init_elevation  ***
386      !!----------------------------------------------------------------------
387      TYPE(OBC_INDEX) , INTENT(in   ) ::   idx   ! OBC indices
388      TYPE(TIDES_DATA), INTENT(inout) ::   td    ! tidal harmonics data
389      !
390      INTEGER ::   itide, igrd, ib       ! dummy loop indices
391      INTEGER, DIMENSION(1) ::   ilen0   ! length of boundary data (from OBC arrays)
392      REAL(wp),ALLOCATABLE, DIMENSION(:) ::   mod_tide, phi_tide
393      !!----------------------------------------------------------------------
394      !
395      igrd=1   
396                              ! SSH on tracer grid.
397      ilen0(1) =  SIZE(td%ssh0(:,1,1))
398      !
399      ALLOCATE( mod_tide(ilen0(igrd)), phi_tide(ilen0(igrd)) )
400      !
401      DO itide = 1, nb_harmo
402         DO ib = 1, ilen0(igrd)
403            mod_tide(ib)=SQRT(td%ssh0(ib,itide,1)**2.+td%ssh0(ib,itide,2)**2.)
404            phi_tide(ib)=ATAN2(-td%ssh0(ib,itide,2),td%ssh0(ib,itide,1))
405         END DO
406         DO ib = 1 , ilen0(igrd)
407            mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f
408            phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0+tide_harmonics(itide)%u
409         ENDDO
410         DO ib = 1 , ilen0(igrd)
411            td%ssh(ib,itide,1)= mod_tide(ib)*COS(phi_tide(ib))
412            td%ssh(ib,itide,2)=-mod_tide(ib)*SIN(phi_tide(ib))
413         ENDDO
414      END DO
415      !
416      DEALLOCATE( mod_tide, phi_tide )
417      !
418   END SUBROUTINE tide_init_elevation
419
420
421   SUBROUTINE tide_init_velocities( idx, td )
422      !!----------------------------------------------------------------------
423      !!                 ***  ROUTINE tide_init_elevation  ***
424      !!----------------------------------------------------------------------
425      TYPE(OBC_INDEX) , INTENT(in   ) ::   idx   ! OBC indices
426      TYPE(TIDES_DATA), INTENT(inout) ::   td    ! tidal harmonics data
427      !
428      INTEGER ::   itide, igrd, ib       ! dummy loop indices
429      INTEGER, DIMENSION(3) ::   ilen0   ! length of boundary data (from OBC arrays)
430      REAL(wp),ALLOCATABLE, DIMENSION(:) ::   mod_tide, phi_tide
431      !!----------------------------------------------------------------------
432      !
433      ilen0(2) =  SIZE(td%u0(:,1,1))
434      ilen0(3) =  SIZE(td%v0(:,1,1))
435      !
436      igrd=2                                 ! U grid.
437      !
438      ALLOCATE( mod_tide(ilen0(igrd)) , phi_tide(ilen0(igrd)) )
439      !
440      DO itide = 1, nb_harmo
441         DO ib = 1, ilen0(igrd)
442            mod_tide(ib)=SQRT(td%u0(ib,itide,1)**2.+td%u0(ib,itide,2)**2.)
443            phi_tide(ib)=ATAN2(-td%u0(ib,itide,2),td%u0(ib,itide,1))
444         END DO
445         DO ib = 1, ilen0(igrd)
446            mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f
447            phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0 + tide_harmonics(itide)%u
448         ENDDO
449         DO ib = 1, ilen0(igrd)
450            td%u(ib,itide,1)= mod_tide(ib)*COS(phi_tide(ib))
451            td%u(ib,itide,2)=-mod_tide(ib)*SIN(phi_tide(ib))
452         ENDDO
453      END DO
454      !
455      DEALLOCATE( mod_tide , phi_tide )
456      !
457      igrd=3                                 ! V grid.
458      !
459      ALLOCATE( mod_tide(ilen0(igrd)) , phi_tide(ilen0(igrd)) )
460
461      DO itide = 1, nb_harmo
462         DO ib = 1, ilen0(igrd)
463            mod_tide(ib)=SQRT(td%v0(ib,itide,1)**2.+td%v0(ib,itide,2)**2.)
464            phi_tide(ib)=ATAN2(-td%v0(ib,itide,2),td%v0(ib,itide,1))
465         END DO
466         DO ib = 1, ilen0(igrd)
467            mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f
468            phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0 + tide_harmonics(itide)%u
469         ENDDO
470         DO ib = 1, ilen0(igrd)
471            td%v(ib,itide,1)= mod_tide(ib)*COS(phi_tide(ib))
472            td%v(ib,itide,2)=-mod_tide(ib)*SIN(phi_tide(ib))
473         ENDDO
474      END DO
475      !
476      DEALLOCATE( mod_tide, phi_tide )
477      !
478  END SUBROUTINE tide_init_velocities
479
480   !!======================================================================
481END MODULE bdytides
482
Note: See TracBrowser for help on using the repository browser.