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 3970 for branches/2013/dev_r3867_MERCATOR1_DYN/NEMOGCM/NEMO/OPA_SRC/BDY/bdytides.F90 – NEMO

Ignore:
Timestamp:
2013-07-11T15:59:14+02:00 (11 years ago)
Author:
cbricaud
Message:

Time splitting update, see ticket #1079

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2013/dev_r3867_MERCATOR1_DYN/NEMOGCM/NEMO/OPA_SRC/BDY/bdytides.F90

    r3651 r3970  
    99   !!            3.3  !  2010-09  (D.Storkey and E.O'Dea)  bug fixes 
    1010   !!            3.4  !  2012-09  (G. Reffray and J. Chanut) New inputs + mods 
     11   !!            3.5  !  2013-07  (J. Chanut) Compliant with time splitting changes 
    1112   !!---------------------------------------------------------------------- 
    1213#if defined key_bdy 
     
    3233!   USE tide_mod       ! Useless ?? 
    3334   USE fldread, ONLY: fld_map 
     35   USE dynspg_oce, ONLY: lk_dynspg_ts 
    3436 
    3537   IMPLICIT NONE 
     
    3840   PUBLIC   bdytide_init     ! routine called in bdy_init 
    3941   PUBLIC   bdytide_update   ! routine called in bdy_dta 
     42   PUBLIC   bdy_dta_tides    ! routine called in dyn_spg_ts 
    4043 
    4144   TYPE, PUBLIC ::   TIDES_DATA     !: Storage for external tidal harmonics data 
     
    4952 
    5053   TYPE(TIDES_DATA), PUBLIC, DIMENSION(jp_bdy), TARGET :: tides  !: External tidal harmonics data 
     54   TYPE(OBC_DATA)  , PRIVATE, DIMENSION(jp_bdy) :: dta_bdy_s  !: bdy external data (slow component) 
    5155 
    5256   !!---------------------------------------------------------------------- 
     
    252256            ENDIF 
    253257            ! 
     258            IF ( lk_dynspg_ts ) THEN ! Allocate arrays to save slowly varying boundary data during 
     259                                     ! time splitting integration 
     260               ALLOCATE( dta_bdy_s(ib_bdy)%ssh ( ilen0(1) ) ) 
     261               ALLOCATE( dta_bdy_s(ib_bdy)%u2d ( ilen0(2) ) ) 
     262               ALLOCATE( dta_bdy_s(ib_bdy)%v2d ( ilen0(3) ) ) 
     263               dta_bdy_s(ib_bdy)%ssh(:) = 0.e0 
     264               dta_bdy_s(ib_bdy)%u2d(:) = 0.e0 
     265               dta_bdy_s(ib_bdy)%v2d(:) = 0.e0 
     266            ENDIF 
     267            ! 
    254268         ENDIF ! nn_dyn2d_dta(ib_bdy) .ge. 2 
    255269         ! 
     
    318332          
    319333      IF( PRESENT(jit) ) THEN   
    320          z_arg = ( ((kt-kt_tide)-1) * rdt + (jit+time_add) * rdt / REAL(nn_baro,wp) ) 
     334         z_arg = ((kt-kt_tide) * rdt + (jit+0.5*(time_add-1)) * rdt / REAL(nn_baro,wp) ) 
    321335      ELSE                               
    322336         z_arg = ((kt-kt_tide)+time_add) * rdt 
     
    324338 
    325339      ! Linear ramp on tidal component at open boundaries  
    326       zramp = 1. 
    327       IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg + (kt_tide-nit000)*rdt)/(rdttideramp*rday),0.),1.) 
     340      zramp = 1._wp 
     341      IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg + (kt_tide-nit000)*rdt)/(rdttideramp*rday),0._wp),1._wp) 
    328342 
    329343      DO itide = 1, nb_harmo 
     
    351365      ! 
    352366   END SUBROUTINE bdytide_update 
     367 
     368   SUBROUTINE bdy_dta_tides( kt, kit, time_offset ) 
     369      !!---------------------------------------------------------------------- 
     370      !!                 ***  SUBROUTINE bdy_dta_tides  *** 
     371      !!                 
     372      !! ** Purpose : - Add tidal forcing to ssh, u2d and v2d OBC data arrays.  
     373      !!                 
     374      !!---------------------------------------------------------------------- 
     375      INTEGER, INTENT( in )            ::   kt          ! Main timestep counter 
     376      INTEGER, INTENT( in ),OPTIONAL   ::   kit         ! Barotropic timestep counter (for timesplitting option) 
     377      INTEGER, INTENT( in ),OPTIONAL   ::   time_offset ! time offset in units of timesteps. NB. if kit 
     378                                                        ! is present then units = subcycle timesteps. 
     379                                                        ! time_offset = 0  => get data at "now"    time level 
     380                                                        ! time_offset = -1 => get data at "before" time level 
     381                                                        ! time_offset = +1 => get data at "after"  time level 
     382                                                        ! etc. 
     383      !! 
     384      LOGICAL  :: lk_first_btstp  ! =.TRUE. if time splitting and first barotropic step 
     385      INTEGER,          DIMENSION(jpbgrd) :: ilen0  
     386      INTEGER, DIMENSION(1:jpbgrd) :: nblen, nblenrim  ! short cuts 
     387      INTEGER  :: itide, ib_bdy, ib, igrd                     ! loop indices 
     388      INTEGER  :: time_add                                    ! time offset in units of timesteps 
     389      REAL(wp) :: z_arg, z_sarg, zramp, zoff, z_cost, z_sist       
     390      !!---------------------------------------------------------------------- 
     391 
     392      IF( nn_timing == 1 ) CALL timing_start('bdy_dta_tides') 
     393 
     394      lk_first_btstp=.TRUE. 
     395      IF ( PRESENT(kit).AND.( kit /= 1 ) ) THEN ; lk_first_btstp=.FALSE. ; ENDIF 
     396 
     397      time_add = 0 
     398      IF( PRESENT(time_offset) ) THEN 
     399         time_add = time_offset 
     400      ENDIF 
     401       
     402      ! Absolute time from model initialization:    
     403      IF( PRESENT(kit) ) THEN   
     404         z_arg = ( kt + (kit+0.5*(time_add-1)) / REAL(nn_baro,wp) ) * rdt 
     405      ELSE                               
     406         z_arg = ( kt + time_add ) * rdt 
     407      ENDIF 
     408 
     409      ! Linear ramp on tidal component at open boundaries  
     410      zramp = 1. 
     411      IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg - nit000*rdt)/(rdttideramp*rday),0.),1.) 
     412 
     413      DO ib_bdy = 1,nb_bdy 
     414 
     415         ! line below should be simplified (runoff case) 
     416         IF (( nn_dyn2d_dta(ib_bdy) .ge. 2 ).AND.(nn_tra(ib_bdy).NE.4)) THEN 
     417 
     418            nblen(1:jpbgrd) = idx_bdy(ib_bdy)%nblen(1:jpbgrd) 
     419            nblenrim(1:jpbgrd) = idx_bdy(ib_bdy)%nblenrim(1:jpbgrd) 
     420 
     421            IF( nn_dyn2d(ib_bdy) .eq. jp_frs ) THEN 
     422               ilen0(:)=nblen(:) 
     423            ELSE 
     424               ilen0(:)=nblenrim(:) 
     425            ENDIF      
     426 
     427            ! We refresh nodal factors every day below 
     428            ! This should be done somewhere else 
     429            IF ( nsec_day == NINT(0.5 * rdttra(1)) .AND. lk_first_btstp ) THEN 
     430               ! 
     431               kt_tide = kt                
     432               ! 
     433               IF(lwp) THEN 
     434               WRITE(numout,*) 
     435               WRITE(numout,*) 'bdy_tide_dta : Refresh nodal factors for tidal open bdy data at kt=',kt 
     436               WRITE(numout,*) '~~~~~~~~~~~~~~ ' 
     437               ENDIF 
     438               ! 
     439               CALL tide_init_elevation ( idx=idx_bdy(ib_bdy), td=tides(ib_bdy) ) 
     440               CALL tide_init_velocities( idx=idx_bdy(ib_bdy), td=tides(ib_bdy) ) 
     441               ! 
     442            ENDIF 
     443            zoff = -kt_tide * rdt ! time offset relative to nodal factor computation time 
     444            ! 
     445            ! If time splitting, save data at first barotropic iteration 
     446            IF ( PRESENT(kit) ) THEN 
     447               IF ( lk_first_btstp ) THEN ! Save slow varying open boundary data: 
     448                  dta_bdy_s(ib_bdy)%ssh(1:ilen0(1)) = dta_bdy(ib_bdy)%ssh(1:ilen0(1)) 
     449                  dta_bdy_s(ib_bdy)%u2d(1:ilen0(2)) = dta_bdy(ib_bdy)%u2d(1:ilen0(2)) 
     450                  dta_bdy_s(ib_bdy)%v2d(1:ilen0(3)) = dta_bdy(ib_bdy)%v2d(1:ilen0(3)) 
     451 
     452               ELSE ! Initialize arrays from slow varying open boundary data:             
     453                  dta_bdy(ib_bdy)%ssh(1:ilen0(1)) = dta_bdy_s(ib_bdy)%ssh(1:ilen0(1)) 
     454                  dta_bdy(ib_bdy)%u2d(1:ilen0(2)) = dta_bdy_s(ib_bdy)%u2d(1:ilen0(2)) 
     455                  dta_bdy(ib_bdy)%v2d(1:ilen0(3)) = dta_bdy_s(ib_bdy)%v2d(1:ilen0(3)) 
     456               ENDIF 
     457            ENDIF 
     458            ! 
     459            ! Update open boundary data arrays: 
     460            DO itide = 1, nb_harmo 
     461               ! 
     462               z_sarg = (z_arg + zoff) * omega_tide(itide) 
     463               z_cost = zramp * COS( z_sarg ) 
     464               z_sist = zramp * SIN( z_sarg ) 
     465               ! 
     466               igrd=1                              ! SSH on tracer grid 
     467               DO ib = 1, ilen0(igrd) 
     468                  dta_bdy(ib_bdy)%ssh(ib) = dta_bdy(ib_bdy)%ssh(ib) + & 
     469                     &                      ( tides(ib_bdy)%ssh(ib,itide,1)*z_cost + & 
     470                     &                        tides(ib_bdy)%ssh(ib,itide,2)*z_sist ) 
     471               END DO 
     472               ! 
     473               igrd=2                              ! U grid 
     474               DO ib = 1, ilen0(igrd) 
     475                  dta_bdy(ib_bdy)%u2d(ib) = dta_bdy(ib_bdy)%u2d(ib) + & 
     476                     &                      ( tides(ib_bdy)%u(ib,itide,1)*z_cost + & 
     477                     &                        tides(ib_bdy)%u(ib,itide,2)*z_sist ) 
     478               END DO 
     479               ! 
     480               igrd=3                              ! V grid 
     481               DO ib = 1, ilen0(igrd)  
     482                  dta_bdy(ib_bdy)%v2d(ib) = dta_bdy(ib_bdy)%v2d(ib) + & 
     483                     &                      ( tides(ib_bdy)%v(ib,itide,1)*z_cost + & 
     484                     &                        tides(ib_bdy)%v(ib,itide,2)*z_sist ) 
     485               END DO 
     486            END DO 
     487         END IF 
     488      END DO 
     489      ! 
     490      IF( nn_timing == 1 ) CALL timing_stop('bdy_dta_tides') 
     491      ! 
     492   END SUBROUTINE bdy_dta_tides 
    353493 
    354494   SUBROUTINE tide_init_elevation( idx, td ) 
     
    457597      WRITE(*,*) 'bdytide_update: You should not have seen this print! error?', kt, jit 
    458598   END SUBROUTINE bdytide_update 
     599   SUBROUTINE bdy_dta_tides( kt, jit )   ! Empty routine 
     600      WRITE(*,*) 'bdy_dta_tides: You should not have seen this print! error?', kt, jit 
     601   END SUBROUTINE bdy_dta_tides 
    459602#endif 
    460603 
    461604   !!====================================================================== 
    462605END MODULE bdytides 
     606 
Note: See TracChangeset for help on using the changeset viewer.