- Timestamp:
- 2019-04-12T11:58:31+02:00 (5 years ago)
- Location:
- NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE/BDY
- Files:
-
- 2 edited
-
bdydta.F90 (modified) (4 diffs)
-
bdytides.F90 (modified) (2 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE/BDY/bdydta.F90
r10773 r10865 65 65 CONTAINS 66 66 67 SUBROUTINE bdy_dta( kt, jit,time_offset )67 SUBROUTINE bdy_dta( kt, time_offset ) 68 68 !!---------------------------------------------------------------------- 69 69 !! *** SUBROUTINE bdy_dta *** … … 75 75 !!---------------------------------------------------------------------- 76 76 INTEGER, INTENT(in) :: kt ! ocean time-step index 77 INTEGER, INTENT(in), OPTIONAL :: jit ! subcycle time-step index (for timesplitting option) 78 INTEGER, INTENT(in), OPTIONAL :: time_offset ! time offset in units of timesteps. NB. if jit 79 ! ! is present then units = subcycle timesteps. 77 INTEGER, INTENT(in), OPTIONAL :: time_offset ! time offset in units of timesteps. 80 78 ! ! time_offset = 0 => get data at "now" time level 81 79 ! ! time_offset = -1 => get data at "before" time level … … 94 92 ! Initialise data arrays once for all from initial conditions where required 95 93 !--------------------------------------------------------------------------- 96 IF( kt == nit000 .AND. .NOT.PRESENT(jit)) THEN94 IF( kt == nit000 ) THEN 97 95 98 96 ! Calculate depth-mean currents … … 228 226 IF( nn_dta(jbdy) == 1 ) THEN ! skip this bit if no external data required 229 227 230 IF( PRESENT(jit) ) THEN 231 ! Update barotropic boundary conditions only 232 ! jit is optional argument for fld_read and bdytide_update 233 IF( cn_dyn2d(jbdy) /= 'none' ) THEN 234 IF( nn_dyn2d_dta(jbdy) == 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays 235 IF( dta%ll_ssh ) dta%ssh(:) = 0._wp 236 IF( dta%ll_u2d ) dta%u2d(:) = 0._wp 237 IF( dta%ll_u3d ) dta%v2d(:) = 0._wp 238 ENDIF 239 IF (cn_tra(jbdy) /= 'runoff') THEN 240 IF( nn_dyn2d_dta(jbdy) == 1 .OR. nn_dyn2d_dta(jbdy) == 3 ) THEN 241 242 jend = jstart + dta%nread(2) - 1 243 IF( ln_full_vel_array(jbdy) ) THEN 244 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend), & 245 & kit=jit, kt_offset=time_offset , jpk_bdy=nb_jpk_bdy, & 246 & fvl=ln_full_vel_array(jbdy) ) 247 ELSE 248 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), map=nbmap_ptr(jstart:jend), & 249 & kit=jit, kt_offset=time_offset ) 250 ENDIF 251 252 ! If full velocities in boundary data then extract barotropic velocities from 3D fields 253 IF( ln_full_vel_array(jbdy) .AND. & 254 & ( nn_dyn2d_dta(jbdy) == 1 .OR. nn_dyn2d_dta(jbdy) == 3 .OR. & 255 & nn_dyn3d_dta(jbdy) == 1 ) )THEN 256 257 igrd = 2 ! zonal velocity 258 dta%u2d(:) = 0._wp 259 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 260 ii = idx_bdy(jbdy)%nbi(ib,igrd) 261 ij = idx_bdy(jbdy)%nbj(ib,igrd) 262 DO ik = 1, jpkm1 263 dta%u2d(ib) = dta%u2d(ib) & 264 & + e3u_n(ii,ij,ik) * umask(ii,ij,ik) * dta%u3d(ib,ik) 265 END DO 266 dta%u2d(ib) = dta%u2d(ib) * r1_hu_n(ii,ij) 267 END DO 268 igrd = 3 ! meridional velocity 269 dta%v2d(:) = 0._wp 270 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 271 ii = idx_bdy(jbdy)%nbi(ib,igrd) 272 ij = idx_bdy(jbdy)%nbj(ib,igrd) 273 DO ik = 1, jpkm1 274 dta%v2d(ib) = dta%v2d(ib) & 275 & + e3v_n(ii,ij,ik) * vmask(ii,ij,ik) * dta%v3d(ib,ik) 276 END DO 277 dta%v2d(ib) = dta%v2d(ib) * r1_hv_n(ii,ij) 278 END DO 279 ENDIF 280 ENDIF 281 IF( nn_dyn2d_dta(jbdy) .ge. 2 ) THEN ! update tidal harmonic forcing 282 CALL bdytide_update( kt=kt, idx=idx_bdy(jbdy), dta=dta, td=tides(jbdy), & 283 & jit=jit, time_offset=time_offset ) 284 ENDIF 285 ENDIF 286 ENDIF 228 IF (cn_tra(jbdy) == 'runoff') then ! runoff condition 229 jend = nb_bdy_fld(jbdy) 230 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), & 231 & map=nbmap_ptr(jstart:jend), kt_offset=time_offset ) 232 ! 233 igrd = 2 ! zonal velocity 234 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 235 ii = idx_bdy(jbdy)%nbi(ib,igrd) 236 ij = idx_bdy(jbdy)%nbj(ib,igrd) 237 dta%u2d(ib) = dta%u2d(ib) / ( e2u(ii,ij) * hu_0(ii,ij) ) 238 END DO 239 ! 240 igrd = 3 ! meridional velocity 241 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 242 ii = idx_bdy(jbdy)%nbi(ib,igrd) 243 ij = idx_bdy(jbdy)%nbj(ib,igrd) 244 dta%v2d(ib) = dta%v2d(ib) / ( e1v(ii,ij) * hv_0(ii,ij) ) 245 END DO 287 246 ELSE 288 IF (cn_tra(jbdy) == 'runoff') then ! runoff condition 289 jend = nb_bdy_fld(jbdy) 290 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), & 291 & map=nbmap_ptr(jstart:jend), kt_offset=time_offset ) 292 ! 247 IF( nn_dyn2d_dta(jbdy) == 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays 248 IF( dta%ll_ssh ) dta%ssh(:) = 0._wp 249 IF( dta%ll_u2d ) dta%u2d(:) = 0._wp 250 IF( dta%ll_v2d ) dta%v2d(:) = 0._wp 251 ENDIF 252 IF( dta%nread(1) .gt. 0 ) THEN ! update external data 253 jend = jstart + dta%nread(1) - 1 254 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), & 255 & map=nbmap_ptr(jstart:jend), kt_offset=time_offset, jpk_bdy=nb_jpk_bdy, & 256 & fvl=ln_full_vel_array(jbdy) ) 257 ENDIF 258 ! If full velocities in boundary data then split into barotropic and baroclinic data 259 IF( ln_full_vel_array(jbdy) .and. & 260 & ( nn_dyn2d_dta(jbdy) == 1 .OR. nn_dyn2d_dta(jbdy) == 3 .OR. & 261 & nn_dyn3d_dta(jbdy) == 1 ) ) THEN 293 262 igrd = 2 ! zonal velocity 263 dta%u2d(:) = 0._wp 294 264 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 295 265 ii = idx_bdy(jbdy)%nbi(ib,igrd) 296 266 ij = idx_bdy(jbdy)%nbj(ib,igrd) 297 dta%u2d(ib) = dta%u2d(ib) / ( e2u(ii,ij) * hu_0(ii,ij) ) 267 DO ik = 1, jpkm1 268 dta%u2d(ib) = dta%u2d(ib) & 269 & + e3u_n(ii,ij,ik) * umask(ii,ij,ik) * dta%u3d(ib,ik) 270 END DO 271 dta%u2d(ib) = dta%u2d(ib) * r1_hu_n(ii,ij) 272 DO ik = 1, jpkm1 273 dta%u3d(ib,ik) = dta%u3d(ib,ik) - dta%u2d(ib) 274 END DO 298 275 END DO 299 !300 276 igrd = 3 ! meridional velocity 277 dta%v2d(:) = 0._wp 301 278 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 302 279 ii = idx_bdy(jbdy)%nbi(ib,igrd) 303 280 ij = idx_bdy(jbdy)%nbj(ib,igrd) 304 dta%v2d(ib) = dta%v2d(ib) / ( e1v(ii,ij) * hv_0(ii,ij) ) 281 DO ik = 1, jpkm1 282 dta%v2d(ib) = dta%v2d(ib) & 283 & + e3v_n(ii,ij,ik) * vmask(ii,ij,ik) * dta%v3d(ib,ik) 284 END DO 285 dta%v2d(ib) = dta%v2d(ib) * r1_hv_n(ii,ij) 286 DO ik = 1, jpkm1 287 dta%v3d(ib,ik) = dta%v3d(ib,ik) - dta%v2d(ib) 288 END DO 305 289 END DO 306 ELSE 307 IF( nn_dyn2d_dta(jbdy) == 2 ) THEN ! tidal harmonic forcing ONLY: initialise arrays 308 IF( dta%ll_ssh ) dta%ssh(:) = 0._wp 309 IF( dta%ll_u2d ) dta%u2d(:) = 0._wp 310 IF( dta%ll_v2d ) dta%v2d(:) = 0._wp 311 ENDIF 312 IF( dta%nread(1) .gt. 0 ) THEN ! update external data 313 jend = jstart + dta%nread(1) - 1 314 CALL fld_read( kt=kt, kn_fsbc=1, sd=bf(jstart:jend), & 315 & map=nbmap_ptr(jstart:jend), kt_offset=time_offset, jpk_bdy=nb_jpk_bdy, & 316 & fvl=ln_full_vel_array(jbdy) ) 317 ENDIF 318 ! If full velocities in boundary data then split into barotropic and baroclinic data 319 IF( ln_full_vel_array(jbdy) .and. & 320 & ( nn_dyn2d_dta(jbdy) == 1 .OR. nn_dyn2d_dta(jbdy) == 3 .OR. & 321 & nn_dyn3d_dta(jbdy) == 1 ) ) THEN 322 igrd = 2 ! zonal velocity 323 dta%u2d(:) = 0._wp 324 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 325 ii = idx_bdy(jbdy)%nbi(ib,igrd) 326 ij = idx_bdy(jbdy)%nbj(ib,igrd) 327 DO ik = 1, jpkm1 328 dta%u2d(ib) = dta%u2d(ib) & 329 & + e3u_n(ii,ij,ik) * umask(ii,ij,ik) * dta%u3d(ib,ik) 330 END DO 331 dta%u2d(ib) = dta%u2d(ib) * r1_hu_n(ii,ij) 332 DO ik = 1, jpkm1 333 dta%u3d(ib,ik) = dta%u3d(ib,ik) - dta%u2d(ib) 334 END DO 335 END DO 336 igrd = 3 ! meridional velocity 337 dta%v2d(:) = 0._wp 338 DO ib = 1, idx_bdy(jbdy)%nblen(igrd) 339 ii = idx_bdy(jbdy)%nbi(ib,igrd) 340 ij = idx_bdy(jbdy)%nbj(ib,igrd) 341 DO ik = 1, jpkm1 342 dta%v2d(ib) = dta%v2d(ib) & 343 & + e3v_n(ii,ij,ik) * vmask(ii,ij,ik) * dta%v3d(ib,ik) 344 END DO 345 dta%v2d(ib) = dta%v2d(ib) * r1_hv_n(ii,ij) 346 DO ik = 1, jpkm1 347 dta%v3d(ib,ik) = dta%v3d(ib,ik) - dta%v2d(ib) 348 END DO 349 END DO 350 ENDIF 351 352 ENDIF 353 #if defined key_si3 354 ! convert N-cat fields (input) into jpl-cat (output) 355 IF( cn_ice(jbdy) /= 'none' .AND. nn_ice_dta(jbdy) == 1 ) THEN 356 jfld_hti = jfld_htit(jbdy) 357 jfld_hts = jfld_htst(jbdy) 358 jfld_ai = jfld_ait(jbdy) 359 IF ( jpl /= 1 .AND. nice_cat == 1 ) THEN ! case input cat = 1 360 CALL ice_var_itd ( bf(jfld_hti)%fnow(:,1,1), bf(jfld_hts)%fnow(:,1,1), bf(jfld_ai)%fnow(:,1,1), & 361 & dta_bdy(jbdy)%h_i , dta_bdy(jbdy)%h_s , dta_bdy(jbdy)%a_i ) 362 ELSEIF( jpl /= 1 .AND. nice_cat /= 1 .AND. nice_cat /= jpl ) THEN ! case input cat /=1 and /=jpl 363 CALL ice_var_itd2( bf(jfld_hti)%fnow(:,1,:), bf(jfld_hts)%fnow(:,1,:), bf(jfld_ai)%fnow(:,1,:), & 364 & dta_bdy(jbdy)%h_i , dta_bdy(jbdy)%h_s , dta_bdy(jbdy)%a_i ) 365 ENDIF 366 ENDIF 367 #endif 368 ENDIF 290 ENDIF 291 292 ENDIF 293 #if defined key_si3 294 ! convert N-cat fields (input) into jpl-cat (output) 295 IF( cn_ice(jbdy) /= 'none' .AND. nn_ice_dta(jbdy) == 1 ) THEN 296 jfld_hti = jfld_htit(jbdy) 297 jfld_hts = jfld_htst(jbdy) 298 jfld_ai = jfld_ait(jbdy) 299 IF ( jpl /= 1 .AND. nice_cat == 1 ) THEN ! case input cat = 1 300 CALL ice_var_itd ( bf(jfld_hti)%fnow(:,1,1), bf(jfld_hts)%fnow(:,1,1), bf(jfld_ai)%fnow(:,1,1), & 301 & dta_bdy(jbdy)%h_i , dta_bdy(jbdy)%h_s , dta_bdy(jbdy)%a_i ) 302 ELSEIF( jpl /= 1 .AND. nice_cat /= 1 .AND. nice_cat /= jpl ) THEN ! case input cat /=1 and /=jpl 303 CALL ice_var_itd2( bf(jfld_hti)%fnow(:,1,:), bf(jfld_hts)%fnow(:,1,:), bf(jfld_ai)%fnow(:,1,:), & 304 & dta_bdy(jbdy)%h_i , dta_bdy(jbdy)%h_s , dta_bdy(jbdy)%a_i ) 305 ENDIF 306 ENDIF 307 #endif 369 308 jstart = jstart + dta%nread(1) 370 309 ENDIF ! nn_dta(jbdy) = 1 -
NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE/BDY/bdytides.F90
r10861 r10865 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 … … 263 262 264 263 265 SUBROUTINE bdytide_update( kt, idx, dta, td, jit, time_offset )266 !!----------------------------------------------------------------------267 !! *** SUBROUTINE bdytide_update ***268 !!269 !! ** Purpose : - Add tidal forcing to ssh, u2d and v2d OBC data arrays.270 !!271 !!----------------------------------------------------------------------272 INTEGER , INTENT(in ) :: kt ! Main timestep counter273 TYPE(OBC_INDEX) , INTENT(in ) :: idx ! OBC indices274 TYPE(OBC_DATA) , INTENT(inout) :: dta ! OBC external data275 TYPE(TIDES_DATA) , INTENT(inout) :: td ! tidal harmonics data276 INTEGER, OPTIONAL, INTENT(in ) :: jit ! Barotropic timestep counter (for timesplitting option)277 INTEGER, OPTIONAL, INTENT(in ) :: time_offset ! time offset in units of timesteps. NB. if jit278 ! ! is present then units = subcycle timesteps.279 ! ! time_offset = 0 => get data at "now" time level280 ! ! time_offset = -1 => get data at "before" time level281 ! ! time_offset = +1 => get data at "after" time level282 ! ! etc.283 !284 INTEGER :: itide, igrd, ib ! dummy loop indices285 INTEGER :: time_add ! time offset in units of timesteps286 INTEGER, DIMENSION(3) :: ilen0 ! length of boundary data (from OBC arrays)287 REAL(wp) :: z_arg, z_sarg, zflag, zramp ! local scalars288 REAL(wp), DIMENSION(jpmax_harmo) :: z_sist, z_cost289 !!----------------------------------------------------------------------290 !291 ilen0(1) = SIZE(td%ssh(:,1,1))292 ilen0(2) = SIZE(td%u(:,1,1))293 ilen0(3) = SIZE(td%v(:,1,1))294 295 zflag=1296 IF ( PRESENT(jit) ) THEN297 IF ( jit /= 1 ) zflag=0298 ENDIF299 300 IF ( (nsec_day == NINT(0.5_wp * rdt) .OR. kt==nit000) .AND. zflag==1 ) THEN301 !302 kt_tide = kt - (nsec_day - 0.5_wp * rdt)/rdt303 !304 IF(lwp) THEN305 WRITE(numout,*)306 WRITE(numout,*) 'bdytide_update : (re)Initialization of the tidal bdy forcing at kt=',kt307 WRITE(numout,*) '~~~~~~~~~~~~~~ '308 ENDIF309 !310 CALL tide_init_elevation ( idx, td )311 CALL tide_init_velocities( idx, td )312 !313 ENDIF314 315 time_add = 0316 IF( PRESENT(time_offset) ) THEN317 time_add = time_offset318 ENDIF319 320 IF( PRESENT(jit) ) THEN321 z_arg = ((kt-kt_tide) * rdt + (jit+0.5_wp*(time_add-1)) * rdt / REAL(nn_baro,wp) )322 ELSE323 z_arg = ((kt-kt_tide)+time_add) * rdt324 ENDIF325 326 ! Linear ramp on tidal component at open boundaries327 zramp = 1._wp328 IF (ln_tide_ramp) zramp = MIN(MAX( (z_arg + (kt_tide-nit000)*rdt)/(rn_tide_ramp_dt*rday),0._wp),1._wp)329 330 DO itide = 1, nb_harmo331 z_sarg = z_arg * tide_harmonics(itide)%omega332 z_cost(itide) = COS( z_sarg )333 z_sist(itide) = SIN( z_sarg )334 END DO335 336 DO itide = 1, nb_harmo337 igrd=1 ! SSH on tracer grid338 DO ib = 1, ilen0(igrd)339 dta%ssh(ib) = dta%ssh(ib) + zramp*(td%ssh(ib,itide,1)*z_cost(itide) + td%ssh(ib,itide,2)*z_sist(itide))340 END DO341 igrd=2 ! U grid342 DO ib = 1, ilen0(igrd)343 dta%u2d(ib) = dta%u2d(ib) + zramp*(td%u (ib,itide,1)*z_cost(itide) + td%u (ib,itide,2)*z_sist(itide))344 END DO345 igrd=3 ! V grid346 DO ib = 1, ilen0(igrd)347 dta%v2d(ib) = dta%v2d(ib) + zramp*(td%v (ib,itide,1)*z_cost(itide) + td%v (ib,itide,2)*z_sist(itide))348 END DO349 END DO350 !351 END SUBROUTINE bdytide_update352 353 354 264 SUBROUTINE bdy_dta_tides( kt, kit, time_offset ) 355 265 !!----------------------------------------------------------------------
Note: See TracChangeset
for help on using the changeset viewer.
