Changeset 12117 for NEMO/branches/2019/dev_r11879_ENHANCE-05_SimonM-Harmonic_Analysis/src/OCE/BDY/bdytides.F90
- Timestamp:
- 2019-12-09T10:46:17+01:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11879_ENHANCE-05_SimonM-Harmonic_Analysis/src/OCE/BDY/bdytides.F90
r11536 r12117 18 18 USE phycst ! physical constants 19 19 USE bdy_oce ! ocean open boundary conditions 20 USE tide ini!20 USE tide_mod ! 21 21 USE daymod ! calendar 22 22 ! … … 30 30 31 31 PUBLIC bdytide_init ! routine called in bdy_init 32 PUBLIC bdytide_update ! routine called in bdy_dta33 32 PUBLIC bdy_dta_tides ! routine called in dyn_spg_ts 34 33 … … 45 44 TYPE(OBC_DATA) , PUBLIC, DIMENSION(jp_bdy) :: dta_bdy_s !: bdy external data (slow component) 46 45 46 INTEGER :: kt_tide 47 47 48 !!---------------------------------------------------------------------- 48 49 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 64 65 CHARACTER(len=80) :: filtide !: Filename root for tidal input files 65 66 LOGICAL :: ln_bdytide_2ddta !: If true, read 2d harmonic data 66 LOGICAL :: ln_bdytide_conj !: If true, assume complex conjugate tidal data67 67 !! 68 68 INTEGER :: ib_bdy, itide, ib !: dummy loop indices … … 77 77 TYPE(TIDES_DATA), POINTER :: td !: local short cut 78 78 !! 79 NAMELIST/nambdy_tide/filtide, ln_bdytide_2ddta , ln_bdytide_conj79 NAMELIST/nambdy_tide/filtide, ln_bdytide_2ddta 80 80 !!---------------------------------------------------------------------- 81 81 ! … … 105 105 IF(lwp) WRITE(numout,*) ' Namelist nambdy_tide : tidal harmonic forcing at open boundaries' 106 106 IF(lwp) WRITE(numout,*) ' read tidal data in 2d files: ', ln_bdytide_2ddta 107 IF(lwp) WRITE(numout,*) ' assume complex conjugate : ', ln_bdytide_conj108 107 IF(lwp) WRITE(numout,*) ' Number of tidal components to read: ', nb_harmo 109 108 IF(lwp) THEN 110 109 WRITE(numout,*) ' Tidal components: ' 111 110 DO itide = 1, nb_harmo 112 WRITE(numout,*) ' ', Wave(ntide(itide))%cname_tide111 WRITE(numout,*) ' ', tide_harmonics(itide)%cname_tide 113 112 END DO 114 113 ENDIF … … 151 150 igrd = 1 ! Everything is at T-points here 152 151 DO itide = 1, nb_harmo 153 CALL iom_get( inum, jpdom_autoglo, TRIM( Wave(ntide(itide))%cname_tide)//'_z1', ztr(:,:) )154 CALL iom_get( inum, jpdom_autoglo, TRIM( Wave(ntide(itide))%cname_tide)//'_z2', zti(:,:) )152 CALL iom_get( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_z1', ztr(:,:) ) 153 CALL iom_get( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_z2', zti(:,:) ) 155 154 DO ib = 1, ilen0(igrd) 156 155 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) … … 168 167 igrd = 2 ! Everything is at U-points here 169 168 DO itide = 1, nb_harmo 170 CALL iom_get ( inum, jpdom_autoglo, TRIM( Wave(ntide(itide))%cname_tide)//'_u1', ztr(:,:) )171 CALL iom_get ( inum, jpdom_autoglo, TRIM( Wave(ntide(itide))%cname_tide)//'_u2', zti(:,:) )169 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_u1', ztr(:,:) ) 170 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_u2', zti(:,:) ) 172 171 DO ib = 1, ilen0(igrd) 173 172 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) … … 185 184 igrd = 3 ! Everything is at V-points here 186 185 DO itide = 1, nb_harmo 187 CALL iom_get ( inum, jpdom_autoglo, TRIM( Wave(ntide(itide))%cname_tide)//'_v1', ztr(:,:) )188 CALL iom_get ( inum, jpdom_autoglo, TRIM( Wave(ntide(itide))%cname_tide)//'_v2', zti(:,:) )186 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_v1', ztr(:,:) ) 187 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_v2', zti(:,:) ) 189 188 DO ib = 1, ilen0(igrd) 190 189 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) … … 210 209 DO itide = 1, nb_harmo 211 210 ! ! SSH fields 212 clfile = TRIM(filtide)//TRIM( Wave(ntide(itide))%cname_tide)//'_grid_T.nc'211 clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_T.nc' 213 212 CALL iom_open( clfile, inum ) 214 213 CALL fld_map( inum, 'z1' , dta_read(1:ilen0(1),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,1) ) … … 218 217 CALL iom_close( inum ) 219 218 ! ! U fields 220 clfile = TRIM(filtide)//TRIM( Wave(ntide(itide))%cname_tide)//'_grid_U.nc'219 clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_U.nc' 221 220 CALL iom_open( clfile, inum ) 222 221 CALL fld_map( inum, 'u1' , dta_read(1:ilen0(2),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,2) ) … … 226 225 CALL iom_close( inum ) 227 226 ! ! V fields 228 clfile = TRIM(filtide)//TRIM( Wave(ntide(itide))%cname_tide)//'_grid_V.nc'227 clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_V.nc' 229 228 CALL iom_open( clfile, inum ) 230 229 CALL fld_map( inum, 'v1' , dta_read(1:ilen0(3),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,3) ) … … 240 239 ENDIF ! ln_bdytide_2ddta=.true. 241 240 ! 242 IF( ln_bdytide_conj ) THEN ! assume complex conjugate in data files243 td%ssh0(:,:,2) = - td%ssh0(:,:,2)244 td%u0 (:,:,2) = - td%u0 (:,:,2)245 td%v0 (:,:,2) = - td%v0 (:,:,2)246 ENDIF247 !248 241 ! Allocate slow varying data in the case of time splitting: 249 242 ! Do it anyway because at this stage knowledge of free surface scheme is unknown … … 260 253 ! 261 254 END SUBROUTINE bdytide_init 262 263 264 SUBROUTINE bdytide_update( kt, idx, dta, td, kit, kt_offset )265 !!----------------------------------------------------------------------266 !! *** SUBROUTINE bdytide_update ***267 !!268 !! ** Purpose : - Add tidal forcing to ssh, u2d and v2d OBC data arrays.269 !!270 !!----------------------------------------------------------------------271 INTEGER , INTENT(in ) :: kt ! Main timestep counter272 TYPE(OBC_INDEX) , INTENT(in ) :: idx ! OBC indices273 TYPE(OBC_DATA) , INTENT(inout) :: dta ! OBC external data274 TYPE(TIDES_DATA) , INTENT(inout) :: td ! tidal harmonics data275 INTEGER, OPTIONAL, INTENT(in ) :: kit ! Barotropic timestep counter (for timesplitting option)276 INTEGER, OPTIONAL, INTENT(in ) :: kt_offset ! time offset in units of timesteps. NB. if kit277 ! ! is present then units = subcycle timesteps.278 ! ! kt_offset = 0 => get data at "now" time level279 ! ! kt_offset = -1 => get data at "before" time level280 ! ! kt_offset = +1 => get data at "after" time level281 ! ! etc.282 !283 INTEGER :: itide, igrd, ib ! dummy loop indices284 INTEGER :: time_add ! time offset in units of timesteps285 INTEGER, DIMENSION(3) :: ilen0 ! length of boundary data (from OBC arrays)286 REAL(wp) :: z_arg, z_sarg, zflag, zramp ! local scalars287 REAL(wp), DIMENSION(jpmax_harmo) :: z_sist, z_cost288 !!----------------------------------------------------------------------289 !290 ilen0(1) = SIZE(td%ssh(:,1,1))291 ilen0(2) = SIZE(td%u(:,1,1))292 ilen0(3) = SIZE(td%v(:,1,1))293 294 zflag=1295 IF ( PRESENT(kit) ) THEN296 IF ( kit /= 1 ) zflag=0297 ENDIF298 299 IF ( (nsec_day == NINT(0.5_wp * rdt) .OR. kt==nit000) .AND. zflag==1 ) THEN300 !301 kt_tide = kt - (nsec_day - 0.5_wp * rdt)/rdt302 !303 IF(lwp) THEN304 WRITE(numout,*)305 WRITE(numout,*) 'bdytide_update : (re)Initialization of the tidal bdy forcing at kt=',kt306 WRITE(numout,*) '~~~~~~~~~~~~~~ '307 ENDIF308 !309 CALL tide_init_elevation ( idx, td )310 CALL tide_init_velocities( idx, td )311 !312 ENDIF313 314 time_add = 0315 IF( PRESENT(kt_offset) ) THEN316 time_add = kt_offset317 ENDIF318 319 IF( PRESENT(kit) ) THEN320 z_arg = ((kt-kt_tide) * rdt + (kit+0.5_wp*(time_add-1)) * rdt / REAL(nn_baro,wp) )321 ELSE322 z_arg = ((kt-kt_tide)+time_add) * rdt323 ENDIF324 325 ! Linear ramp on tidal component at open boundaries326 zramp = 1._wp327 IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg + (kt_tide-nit000)*rdt)/(rdttideramp*rday),0._wp),1._wp)328 329 DO itide = 1, nb_harmo330 z_sarg = z_arg * omega_tide(itide)331 z_cost(itide) = COS( z_sarg )332 z_sist(itide) = SIN( z_sarg )333 END DO334 335 DO itide = 1, nb_harmo336 igrd=1 ! SSH on tracer grid337 DO ib = 1, ilen0(igrd)338 dta%ssh(ib) = dta%ssh(ib) + zramp*(td%ssh(ib,itide,1)*z_cost(itide) + td%ssh(ib,itide,2)*z_sist(itide))339 END DO340 igrd=2 ! U grid341 DO ib = 1, ilen0(igrd)342 dta%u2d(ib) = dta%u2d(ib) + zramp*(td%u (ib,itide,1)*z_cost(itide) + td%u (ib,itide,2)*z_sist(itide))343 END DO344 igrd=3 ! V grid345 DO ib = 1, ilen0(igrd)346 dta%v2d(ib) = dta%v2d(ib) + zramp*(td%v (ib,itide,1)*z_cost(itide) + td%v (ib,itide,2)*z_sist(itide))347 END DO348 END DO349 !350 END SUBROUTINE bdytide_update351 255 352 256 … … 392 296 ! Linear ramp on tidal component at open boundaries 393 297 zramp = 1. 394 IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg - nit000*rdt)/(r dttideramp*rday),0.),1.)298 IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg - nit000*rdt)/(rn_tide_ramp_dt*rday),0.),1.) 395 299 396 300 DO ib_bdy = 1,nb_bdy … … 433 337 DO itide = 1, nb_harmo 434 338 ! 435 z_sarg = (z_arg + zoff) * omega_tide(itide)339 z_sarg = (z_arg + zoff) * tide_harmonics(itide)%omega 436 340 z_cost = zramp * COS( z_sarg ) 437 341 z_sist = zramp * SIN( z_sarg ) … … 491 395 END DO 492 396 DO ib = 1 , ilen0(igrd) 493 mod_tide(ib)=mod_tide(ib)* ftide(itide)494 phi_tide(ib)=phi_tide(ib)+ v0tide(itide)+utide(itide)397 mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f 398 phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0+tide_harmonics(itide)%u 495 399 ENDDO 496 400 DO ib = 1 , ilen0(igrd) … … 530 434 END DO 531 435 DO ib = 1, ilen0(igrd) 532 mod_tide(ib)=mod_tide(ib)* ftide(itide)533 phi_tide(ib)=phi_tide(ib)+ v0tide(itide)+utide(itide)436 mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f 437 phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0 + tide_harmonics(itide)%u 534 438 ENDDO 535 439 DO ib = 1, ilen0(igrd) … … 551 455 END DO 552 456 DO ib = 1, ilen0(igrd) 553 mod_tide(ib)=mod_tide(ib)* ftide(itide)554 phi_tide(ib)=phi_tide(ib)+ v0tide(itide)+utide(itide)457 mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f 458 phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0 + tide_harmonics(itide)%u 555 459 ENDDO 556 460 DO ib = 1, ilen0(igrd)
Note: See TracChangeset
for help on using the changeset viewer.
