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.
Changeset 12205 for NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/BDY/bdytides.F90 – NEMO

Ignore:
Timestamp:
2019-12-12T11:52:50+01:00 (4 years ago)
Author:
acc
Message:

2019/dev_r11943_MERGE_2019: Merge in dev_r11879_ENHANCE-05_SimonM-Harmonic_Analysis. SETTE tested

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/BDY/bdytides.F90

    r11960 r12205  
    1818   USE phycst         ! physical constants 
    1919   USE bdy_oce        ! ocean open boundary conditions 
    20    USE tideini        !  
     20   USE tide_mod       !  
    2121   USE daymod         ! calendar 
    2222   ! 
     
    3030 
    3131   PUBLIC   bdytide_init     ! routine called in bdy_init 
    32    PUBLIC   bdytide_update   ! routine called in bdy_dta 
    3332   PUBLIC   bdy_dta_tides    ! routine called in dyn_spg_ts 
    3433 
     
    4544   TYPE(OBC_DATA)  , PUBLIC, DIMENSION(jp_bdy) :: dta_bdy_s  !: bdy external data (slow component) 
    4645 
     46   INTEGER ::   kt_tide 
     47 
    4748   !!---------------------------------------------------------------------- 
    4849   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    6465      CHARACTER(len=80)                         ::   filtide             !: Filename root for tidal input files 
    6566      LOGICAL                                   ::   ln_bdytide_2ddta    !: If true, read 2d harmonic data 
    66       LOGICAL                                   ::   ln_bdytide_conj     !: If true, assume complex conjugate tidal data 
    6767      !! 
    6868      INTEGER                                   ::   ib_bdy, itide, ib   !: dummy loop indices 
     
    7979      TYPE(TIDES_DATA),  POINTER                ::   td                  !: local short cut    
    8080      !! 
    81       NAMELIST/nambdy_tide/filtide, ln_bdytide_2ddta, ln_bdytide_conj 
     81      NAMELIST/nambdy_tide/filtide, ln_bdytide_2ddta 
    8282      !!---------------------------------------------------------------------- 
    8383      ! 
     
    119119            IF(lwp) WRITE(numout,*) '          Namelist nambdy_tide : tidal harmonic forcing at open boundaries' 
    120120            IF(lwp) WRITE(numout,*) '             read tidal data in 2d files: ', ln_bdytide_2ddta 
    121             IF(lwp) WRITE(numout,*) '             assume complex conjugate   : ', ln_bdytide_conj 
    122121            IF(lwp) WRITE(numout,*) '             Number of tidal components to read: ', nb_harmo 
    123122            IF(lwp) THEN  
    124123                    WRITE(numout,*) '             Tidal components: '  
    125124               DO itide = 1, nb_harmo 
    126                   WRITE(numout,*)  '                 ', Wave(ntide(itide))%cname_tide  
     125                  WRITE(numout,*)  '                 ', tide_harmonics(itide)%cname_tide  
    127126               END DO 
    128127            ENDIF  
     
    165164               igrd = 1                       ! Everything is at T-points here 
    166165               DO itide = 1, nb_harmo 
    167                   CALL iom_get( inum, jpdom_autoglo, TRIM(Wave(ntide(itide))%cname_tide)//'_z1', ztr(:,:) ) 
    168                   CALL iom_get( inum, jpdom_autoglo, TRIM(Wave(ntide(itide))%cname_tide)//'_z2', zti(:,:) )  
     166                  CALL iom_get( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_z1', ztr(:,:) ) 
     167                  CALL iom_get( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_z2', zti(:,:) )  
    169168                  DO ib = 1, ilen0(igrd) 
    170169                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
     
    182181               igrd = 2                       ! Everything is at U-points here 
    183182               DO itide = 1, nb_harmo 
    184                   CALL iom_get  ( inum, jpdom_autoglo, TRIM(Wave(ntide(itide))%cname_tide)//'_u1', ztr(:,:) ) 
    185                   CALL iom_get  ( inum, jpdom_autoglo, TRIM(Wave(ntide(itide))%cname_tide)//'_u2', zti(:,:) ) 
     183                  CALL iom_get  ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_u1', ztr(:,:) ) 
     184                  CALL iom_get  ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_u2', zti(:,:) ) 
    186185                  DO ib = 1, ilen0(igrd) 
    187186                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
     
    199198               igrd = 3                       ! Everything is at V-points here 
    200199               DO itide = 1, nb_harmo 
    201                   CALL iom_get  ( inum, jpdom_autoglo, TRIM(Wave(ntide(itide))%cname_tide)//'_v1', ztr(:,:) ) 
    202                   CALL iom_get  ( inum, jpdom_autoglo, TRIM(Wave(ntide(itide))%cname_tide)//'_v2', zti(:,:) ) 
     200                  CALL iom_get  ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_v1', ztr(:,:) ) 
     201                  CALL iom_get  ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_v2', zti(:,:) ) 
    203202                  DO ib = 1, ilen0(igrd) 
    204203                     ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
     
    224223               DO itide = 1, nb_harmo 
    225224                  !                                                              ! SSH fields 
    226                   clfile = TRIM(filtide)//TRIM(Wave(ntide(itide))%cname_tide)//'_grid_T.nc' 
     225                  clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_T.nc' 
    227226                  CALL iom_open( clfile, inum ) 
    228227                  CALL fld_map( inum, 'z1' , dta_read(1:ilen0(1),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,1) ) 
     
    232231                  CALL iom_close( inum ) 
    233232                  !                                                              ! U fields 
    234                   clfile = TRIM(filtide)//TRIM(Wave(ntide(itide))%cname_tide)//'_grid_U.nc' 
     233                  clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_U.nc' 
    235234                  CALL iom_open( clfile, inum ) 
    236235                  CALL fld_map( inum, 'u1' , dta_read(1:ilen0(2),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,2) ) 
     
    240239                  CALL iom_close( inum ) 
    241240                  !                                                              ! V fields 
    242                   clfile = TRIM(filtide)//TRIM(Wave(ntide(itide))%cname_tide)//'_grid_V.nc' 
     241                  clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_V.nc' 
    243242                  CALL iom_open( clfile, inum ) 
    244243                  CALL fld_map( inum, 'v1' , dta_read(1:ilen0(3),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,3) ) 
     
    254253            ENDIF ! ln_bdytide_2ddta=.true. 
    255254            ! 
    256             IF( ln_bdytide_conj ) THEN    ! assume complex conjugate in data files 
    257                td%ssh0(:,:,2) = - td%ssh0(:,:,2) 
    258                td%u0  (:,:,2) = - td%u0  (:,:,2) 
    259                td%v0  (:,:,2) = - td%v0  (:,:,2) 
    260             ENDIF 
    261             ! 
    262255            ! Allocate slow varying data in the case of time splitting: 
    263256            ! Do it anyway because at this stage knowledge of free surface scheme is unknown 
     
    274267      ! 
    275268   END SUBROUTINE bdytide_init 
    276  
    277  
    278    SUBROUTINE bdytide_update( kt, idx, dta, td, kit, kt_offset ) 
    279       !!---------------------------------------------------------------------- 
    280       !!                 ***  SUBROUTINE bdytide_update  *** 
    281       !!                 
    282       !! ** Purpose : - Add tidal forcing to ssh, u2d and v2d OBC data arrays.  
    283       !!                 
    284       !!---------------------------------------------------------------------- 
    285       INTEGER          , INTENT(in   ) ::   kt          ! Main timestep counter 
    286       TYPE(OBC_INDEX)  , INTENT(in   ) ::   idx         ! OBC indices 
    287       TYPE(OBC_DATA)   , INTENT(inout) ::   dta         ! OBC external data 
    288       TYPE(TIDES_DATA) , INTENT(inout) ::   td          ! tidal harmonics data 
    289       INTEGER, OPTIONAL, INTENT(in   ) ::   kit         ! Barotropic timestep counter (for timesplitting option) 
    290       INTEGER, OPTIONAL, INTENT(in   ) ::   kt_offset   ! time offset in units of timesteps. NB. if kit 
    291       !                                                 ! is present then units = subcycle timesteps. 
    292       !                                                 ! kt_offset = 0  => get data at "now"    time level 
    293       !                                                 ! kt_offset = -1 => get data at "before" time level 
    294       !                                                 ! kt_offset = +1 => get data at "after"  time level 
    295       !                                                 ! etc. 
    296       ! 
    297       INTEGER  ::   itide, igrd, ib       ! dummy loop indices 
    298       INTEGER  ::   time_add              ! time offset in units of timesteps 
    299       INTEGER, DIMENSION(3) ::   ilen0    ! length of boundary data (from OBC arrays) 
    300       REAL(wp) ::   z_arg, z_sarg, zflag, zramp   ! local scalars     
    301       REAL(wp), DIMENSION(jpmax_harmo) :: z_sist, z_cost 
    302       !!---------------------------------------------------------------------- 
    303       ! 
    304       ilen0(1) =  SIZE(td%ssh(:,1,1)) 
    305       ilen0(2) =  SIZE(td%u(:,1,1)) 
    306       ilen0(3) =  SIZE(td%v(:,1,1)) 
    307  
    308       zflag=1 
    309       IF ( PRESENT(kit) ) THEN 
    310         IF ( kit /= 1 ) zflag=0 
    311       ENDIF 
    312  
    313       IF ( (nsec_day == NINT(0.5_wp * rdt) .OR. kt==nit000) .AND. zflag==1 ) THEN 
    314         ! 
    315         kt_tide = kt - (nsec_day - 0.5_wp * rdt)/rdt 
    316         ! 
    317         IF(lwp) THEN 
    318            WRITE(numout,*) 
    319            WRITE(numout,*) 'bdytide_update : (re)Initialization of the tidal bdy forcing at kt=',kt 
    320            WRITE(numout,*) '~~~~~~~~~~~~~~ ' 
    321         ENDIF 
    322         ! 
    323         CALL tide_init_elevation ( idx, td ) 
    324         CALL tide_init_velocities( idx, td ) 
    325         ! 
    326       ENDIF  
    327  
    328       time_add = 0 
    329       IF( PRESENT(kt_offset) ) THEN 
    330          time_add = kt_offset 
    331       ENDIF 
    332           
    333       IF( PRESENT(kit) ) THEN   
    334          z_arg = ((kt-kt_tide) * rdt + (kit+0.5_wp*(time_add-1)) * rdt / REAL(nn_baro,wp) ) 
    335       ELSE                               
    336          z_arg = ((kt-kt_tide)+time_add) * rdt 
    337       ENDIF 
    338  
    339       ! Linear ramp on tidal component at open boundaries  
    340       zramp = 1._wp 
    341       IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg + (kt_tide-nit000)*rdt)/(rdttideramp*rday),0._wp),1._wp) 
    342  
    343       DO itide = 1, nb_harmo 
    344          z_sarg = z_arg * omega_tide(itide) 
    345          z_cost(itide) = COS( z_sarg ) 
    346          z_sist(itide) = SIN( z_sarg ) 
    347       END DO 
    348  
    349       DO itide = 1, nb_harmo 
    350          igrd=1                              ! SSH on tracer grid 
    351          DO ib = 1, ilen0(igrd) 
    352             dta%ssh(ib) = dta%ssh(ib) + zramp*(td%ssh(ib,itide,1)*z_cost(itide) + td%ssh(ib,itide,2)*z_sist(itide)) 
    353          END DO 
    354          igrd=2                              ! U grid 
    355          DO ib = 1, ilen0(igrd) 
    356             dta%u2d(ib) = dta%u2d(ib) + zramp*(td%u  (ib,itide,1)*z_cost(itide) + td%u  (ib,itide,2)*z_sist(itide)) 
    357          END DO 
    358          igrd=3                              ! V grid 
    359          DO ib = 1, ilen0(igrd)  
    360             dta%v2d(ib) = dta%v2d(ib) + zramp*(td%v  (ib,itide,1)*z_cost(itide) + td%v  (ib,itide,2)*z_sist(itide)) 
    361          END DO 
    362       END DO 
    363       ! 
    364    END SUBROUTINE bdytide_update 
    365269 
    366270 
     
    406310      ! Linear ramp on tidal component at open boundaries  
    407311      zramp = 1. 
    408       IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg - nit000*rdt)/(rdttideramp*rday),0.),1.) 
     312      IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg - nit000*rdt)/(rn_tide_ramp_dt*rday),0.),1.) 
    409313 
    410314      DO ib_bdy = 1,nb_bdy 
     
    447351            DO itide = 1, nb_harmo 
    448352               ! 
    449                z_sarg = (z_arg + zoff) * omega_tide(itide) 
     353               z_sarg = (z_arg + zoff) * tide_harmonics(itide)%omega 
    450354               z_cost = zramp * COS( z_sarg ) 
    451355               z_sist = zramp * SIN( z_sarg ) 
     
    505409         END DO 
    506410         DO ib = 1 , ilen0(igrd) 
    507             mod_tide(ib)=mod_tide(ib)*ftide(itide) 
    508             phi_tide(ib)=phi_tide(ib)+v0tide(itide)+utide(itide) 
     411            mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f 
     412            phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0+tide_harmonics(itide)%u 
    509413         ENDDO 
    510414         DO ib = 1 , ilen0(igrd) 
     
    544448         END DO 
    545449         DO ib = 1, ilen0(igrd) 
    546             mod_tide(ib)=mod_tide(ib)*ftide(itide) 
    547             phi_tide(ib)=phi_tide(ib)+v0tide(itide)+utide(itide) 
     450            mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f 
     451            phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0 + tide_harmonics(itide)%u 
    548452         ENDDO 
    549453         DO ib = 1, ilen0(igrd) 
     
    565469         END DO 
    566470         DO ib = 1, ilen0(igrd) 
    567             mod_tide(ib)=mod_tide(ib)*ftide(itide) 
    568             phi_tide(ib)=phi_tide(ib)+v0tide(itide)+utide(itide) 
     471            mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f 
     472            phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0 + tide_harmonics(itide)%u 
    569473         ENDDO 
    570474         DO ib = 1, ilen0(igrd) 
Note: See TracChangeset for help on using the changeset viewer.