Changeset 11822 for NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/BDY/bdydyn3d.F90
- Timestamp:
- 2019-10-29T11:41:36+01:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/BDY/bdydyn3d.F90
r10957 r11822 44 44 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT( inout ) :: puu, pvv ! Ocean velocities (to be updated at open boundaries) 45 45 ! 46 INTEGER :: ib_bdy ! loop index 47 !!---------------------------------------------------------------------- 48 ! 49 DO ib_bdy=1, nb_bdy 46 INTEGER :: ib_bdy, ir ! BDY set index, rim index 47 LOGICAL :: llrim0 ! indicate if rim 0 is treated 48 LOGICAL, DIMENSION(4) :: llsend2, llrecv2, llsend3, llrecv3 ! indicate how communications are to be carried out 49 50 !!---------------------------------------------------------------------- 51 llsend2(:) = .false. ; llrecv2(:) = .false. 52 llsend3(:) = .false. ; llrecv3(:) = .false. 53 DO ir = 1, 0, -1 ! treat rim 1 before rim 0 54 IF( ir == 0 ) THEN ; llrim0 = .TRUE. 55 ELSE ; llrim0 = .FALSE. 56 END IF 57 DO ib_bdy=1, nb_bdy 58 ! 59 SELECT CASE( cn_dyn3d(ib_bdy) ) 60 CASE('none') ; CYCLE 61 CASE('frs' ) ! treat the whole boundary at once 62 IF( ir == 0) CALL bdy_dyn3d_frs( puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 63 CASE('specified') ! treat the whole rim at once 64 IF( ir == 0) CALL bdy_dyn3d_spe( puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 65 CASE('zero') ! treat the whole rim at once 66 IF( ir == 0) CALL bdy_dyn3d_zro( puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 67 CASE('orlanski' ) ; CALL bdy_dyn3d_orlanski( Kbb, puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, llrim0, ll_npo=.false. ) 68 CASE('orlanski_npo'); CALL bdy_dyn3d_orlanski( Kbb, puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, llrim0, ll_npo=.true. ) 69 CASE('zerograd') ; CALL bdy_dyn3d_zgrad( puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy, llrim0 ) 70 CASE('neumann') ; CALL bdy_dyn3d_nmn( puu, pvv, Kaa, idx_bdy(ib_bdy), ib_bdy, llrim0 ) 71 CASE DEFAULT ; CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' ) 72 END SELECT 73 END DO 50 74 ! 51 SELECT CASE( cn_dyn3d(ib_bdy) ) 52 CASE('none') ; CYCLE 53 CASE('frs' ) ; CALL bdy_dyn3d_frs( puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy ) 54 CASE('specified') ; CALL bdy_dyn3d_spe( puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy ) 55 CASE('zero') ; CALL bdy_dyn3d_zro( puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy ) 56 CASE('orlanski' ) ; CALL bdy_dyn3d_orlanski( Kbb, puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.false. ) 57 CASE('orlanski_npo'); CALL bdy_dyn3d_orlanski( Kbb, puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.true. ) 58 CASE('zerograd') ; CALL bdy_dyn3d_zgrad( puu, pvv, Kaa, idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy ) 59 CASE('neumann') ; CALL bdy_dyn3d_nmn( puu, pvv, Kaa, idx_bdy(ib_bdy), ib_bdy ) 60 CASE DEFAULT ; CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' ) 61 END SELECT 62 END DO 75 IF( nn_hls > 1 .AND. ir == 1 ) CYCLE ! at least 2 halos will be corrected -> no need to correct rim 1 before rim 0 76 IF( nn_hls == 1 ) THEN 77 llsend2(:) = .false. ; llrecv2(:) = .false. 78 llsend3(:) = .false. ; llrecv3(:) = .false. 79 END IF 80 DO ib_bdy=1, nb_bdy 81 SELECT CASE( cn_dyn3d(ib_bdy) ) 82 CASE('orlanski', 'orlanski_npo') 83 llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:,ir) ! possibly every direction, U points 84 llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:,ir) ! possibly every direction, U points 85 llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:,ir) ! possibly every direction, V points 86 llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:,ir) ! possibly every direction, V points 87 CASE('zerograd') 88 llsend2(3:4) = llsend2(3:4) .OR. lsend_bdyint(ib_bdy,2,3:4,ir) ! north/south, U points 89 llrecv2(3:4) = llrecv2(3:4) .OR. lrecv_bdyint(ib_bdy,2,3:4,ir) ! north/south, U points 90 llsend3(1:2) = llsend3(1:2) .OR. lsend_bdyint(ib_bdy,3,1:2,ir) ! west/east, V points 91 llrecv3(1:2) = llrecv3(1:2) .OR. lrecv_bdyint(ib_bdy,3,1:2,ir) ! west/east, V points 92 CASE('neumann') 93 llsend2(:) = llsend2(:) .OR. lsend_bdyint(ib_bdy,2,:,ir) ! possibly every direction, U points 94 llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(ib_bdy,2,:,ir) ! possibly every direction, U points 95 llsend3(:) = llsend3(:) .OR. lsend_bdyint(ib_bdy,3,:,ir) ! possibly every direction, V points 96 llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(ib_bdy,3,:,ir) ! possibly every direction, V points 97 END SELECT 98 END DO 99 ! 100 IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN ! if need to send/recv in at least one direction 101 CALL lbc_lnk( 'bdydyn2d', puu(:,:,:,Kaa), 'U', -1., kfillmode=jpfillnothing ,lsend=llsend2, lrecv=llrecv2 ) 102 END IF 103 IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN ! if need to send/recv in at least one direction 104 CALL lbc_lnk( 'bdydyn2d', pvv(:,:,:,Kaa), 'V', -1., kfillmode=jpfillnothing ,lsend=llsend3, lrecv=llrecv3 ) 105 END IF 106 END DO ! ir 63 107 ! 64 108 END SUBROUTINE bdy_dyn3d 65 109 66 110 67 SUBROUTINE bdy_dyn3d_spe( puu, pvv, Kaa, idx, dta, ib_bdy )111 SUBROUTINE bdy_dyn3d_spe( puu, pvv, Kaa, idx, dta, kt, ib_bdy ) 68 112 !!---------------------------------------------------------------------- 69 113 !! *** SUBROUTINE bdy_dyn3d_spe *** … … 77 121 TYPE(OBC_INDEX) , INTENT( in ) :: idx ! OBC indices 78 122 TYPE(OBC_DATA) , INTENT( in ) :: dta ! OBC external data 123 INTEGER , INTENT( in ) :: kt ! Time step 79 124 INTEGER , INTENT( in ) :: ib_bdy ! BDY set index 80 125 ! 81 126 INTEGER :: jb, jk ! dummy loop indices 82 127 INTEGER :: ii, ij, igrd ! local integers 83 REAL(wp) :: zwgt ! boundary weight84 128 !!---------------------------------------------------------------------- 85 129 ! … … 101 145 END DO 102 146 END DO 103 CALL lbc_bdy_lnk( 'bdydyn3d', puu(:,:,:,Kaa), 'U', -1., ib_bdy ) ! Boundary points should be updated104 CALL lbc_bdy_lnk( 'bdydyn3d', pvv(:,:,:,Kaa), 'V', -1., ib_bdy )105 147 ! 106 148 END SUBROUTINE bdy_dyn3d_spe 107 149 108 150 109 SUBROUTINE bdy_dyn3d_zgrad( puu, pvv, Kaa, idx, dta, ib_bdy)151 SUBROUTINE bdy_dyn3d_zgrad( puu, pvv, Kaa, idx, dta, kt, ib_bdy, llrim0 ) 110 152 !!---------------------------------------------------------------------- 111 153 !! *** SUBROUTINE bdy_dyn3d_zgrad *** … … 118 160 TYPE(OBC_INDEX) , INTENT( in ) :: idx ! OBC indices 119 161 TYPE(OBC_DATA) , INTENT( in ) :: dta ! OBC external data 120 INTEGER , INTENT( in ) :: ib_bdy ! BDY set index 162 INTEGER , INTENT( in ) :: kt 163 INTEGER , INTENT( in ) :: ib_bdy ! BDY set index 164 LOGICAL , INTENT( in ) :: llrim0 ! indicate if rim 0 is treated 121 165 !! 122 166 INTEGER :: jb, jk ! dummy loop indices 123 167 INTEGER :: ii, ij, igrd ! local integers 124 REAL(wp) :: zwgt ! boundary weight125 INTEGER :: fu, fv168 INTEGER :: flagu, flagv ! short cuts 169 INTEGER :: ibeg, iend ! length of rim to be treated (rim 0 or rim 1 or both) 126 170 !!---------------------------------------------------------------------- 127 171 ! 128 172 igrd = 2 ! Copying tangential velocity into bdy points 129 DO jb = 1, idx%nblenrim(igrd) 130 DO jk = 1, jpkm1 131 ii = idx%nbi(jb,igrd) 132 ij = idx%nbj(jb,igrd) 133 fu = ABS( ABS (NINT( idx%flagu(jb,igrd) ) ) - 1 ) 134 puu(ii,ij,jk,Kaa) = puu(ii,ij,jk,Kaa) * REAL( 1 - fu) + ( puu(ii,ij+fu,jk,Kaa) * umask(ii,ij+fu,jk) & 135 &+ puu(ii,ij-fu,jk,Kaa) * umask(ii,ij-fu,jk) ) * umask(ii,ij,jk) * REAL( fu ) 136 END DO 173 IF( llrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(igrd) 174 ELSE ; ibeg = idx%nblenrim0(igrd)+1 ; iend = idx%nblenrim(igrd) 175 ENDIF 176 DO jb = ibeg, iend 177 ii = idx%nbi(jb,igrd) 178 ij = idx%nbj(jb,igrd) 179 flagu = NINT(idx%flagu(jb,igrd)) 180 flagv = NINT(idx%flagv(jb,igrd)) 181 ! 182 IF( flagu == 0 ) THEN ! north/south bdy 183 IF( ij+flagv > jpj .OR. ij+flagv < 1 ) CYCLE 184 ! 185 DO jk = 1, jpkm1 186 puu(ii,ij,jk,Kaa) = puu(ii,ij+flagv,jk,Kaa) * umask(ii,ij+flagv,jk) 187 END DO 188 ! 189 END IF 137 190 END DO 138 191 ! 139 192 igrd = 3 ! Copying tangential velocity into bdy points 140 DO jb = 1, idx%nblenrim(igrd) 141 DO jk = 1, jpkm1 142 ii = idx%nbi(jb,igrd) 143 ij = idx%nbj(jb,igrd) 144 fv = ABS( ABS (NINT( idx%flagv(jb,igrd) ) ) - 1 ) 145 pvv(ii,ij,jk,Kaa) = pvv(ii,ij,jk,Kaa) * REAL( 1 - fv ) + ( pvv(ii+fv,ij,jk,Kaa) * vmask(ii+fv,ij,jk) & 146 &+ pvv(ii-fv,ij,jk,Kaa) * vmask(ii-fv,ij,jk) ) * vmask(ii,ij,jk) * REAL( fv ) 147 END DO 148 END DO 149 CALL lbc_bdy_lnk( 'bdydyn3d', puu(:,:,:,Kaa), 'U', -1., ib_bdy ) ! Boundary points should be updated 150 CALL lbc_bdy_lnk( 'bdydyn3d', pvv(:,:,:,Kaa), 'V', -1., ib_bdy ) 193 IF( llrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(igrd) 194 ELSE ; ibeg = idx%nblenrim0(igrd)+1 ; iend = idx%nblenrim(igrd) 195 ENDIF 196 DO jb = ibeg, iend 197 ii = idx%nbi(jb,igrd) 198 ij = idx%nbj(jb,igrd) 199 flagu = NINT(idx%flagu(jb,igrd)) 200 flagv = NINT(idx%flagv(jb,igrd)) 201 ! 202 IF( flagv == 0 ) THEN ! west/east bdy 203 IF( ii+flagu > jpi .OR. ii+flagu < 1 ) CYCLE 204 ! 205 DO jk = 1, jpkm1 206 pvv(ii,ij,jk,Kaa) = pvv(ii+flagu,ij,jk,Kaa) * vmask(ii+flagu,ij,jk) 207 END DO 208 ! 209 END IF 210 END DO 151 211 ! 152 212 END SUBROUTINE bdy_dyn3d_zgrad 153 213 154 214 155 SUBROUTINE bdy_dyn3d_zro( puu, pvv, Kaa, idx, dta, ib_bdy )215 SUBROUTINE bdy_dyn3d_zro( puu, pvv, Kaa, idx, dta, kt, ib_bdy ) 156 216 !!---------------------------------------------------------------------- 157 217 !! *** SUBROUTINE bdy_dyn3d_zro *** … … 160 220 !! 161 221 !!---------------------------------------------------------------------- 222 INTEGER , INTENT( in ) :: kt ! time step index 162 223 INTEGER , INTENT( in ) :: Kaa ! Time level index 163 224 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT( inout ) :: puu, pvv ! Ocean velocities (to be updated at open boundaries) … … 168 229 INTEGER :: ib, ik ! dummy loop indices 169 230 INTEGER :: ii, ij, igrd ! local integers 170 REAL(wp) :: zwgt ! boundary weight171 231 !!---------------------------------------------------------------------- 172 232 ! … … 179 239 END DO 180 240 END DO 181 241 ! 182 242 igrd = 3 ! Everything is at T-points here 183 243 DO ib = 1, idx%nblenrim(igrd) … … 189 249 END DO 190 250 ! 191 CALL lbc_bdy_lnk( 'bdydyn3d', puu(:,:,:,Kaa), 'U', -1., ib_bdy ) ; CALL lbc_bdy_lnk( 'bdydyn3d', pvv(:,:,:,Kaa), 'V', -1.,ib_bdy ) ! Boundary points should be updated192 !193 251 END SUBROUTINE bdy_dyn3d_zro 194 252 195 253 196 SUBROUTINE bdy_dyn3d_frs( puu, pvv, Kaa, idx, dta, ib_bdy )254 SUBROUTINE bdy_dyn3d_frs( puu, pvv, Kaa, idx, dta, kt, ib_bdy ) 197 255 !!---------------------------------------------------------------------- 198 256 !! *** SUBROUTINE bdy_dyn3d_frs *** … … 205 263 !! topography. Tellus, 365-382. 206 264 !!---------------------------------------------------------------------- 265 INTEGER , INTENT( in ) :: kt ! time step index 207 266 INTEGER , INTENT( in ) :: Kaa ! Time level index 208 267 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT( inout ) :: puu, pvv ! Ocean velocities (to be updated at open boundaries) … … 234 293 pvv(ii,ij,jk,Kaa) = ( pvv(ii,ij,jk,Kaa) + zwgt * ( dta%v3d(jb,jk) - pvv(ii,ij,jk,Kaa) ) ) * vmask(ii,ij,jk) 235 294 END DO 236 END DO 237 CALL lbc_bdy_lnk( 'bdydyn3d', puu(:,:,:,Kaa), 'U', -1., ib_bdy ) ! Boundary points should be updated 238 CALL lbc_bdy_lnk( 'bdydyn3d', pvv(:,:,:,Kaa), 'V', -1., ib_bdy ) 295 END DO 239 296 ! 240 297 END SUBROUTINE bdy_dyn3d_frs 241 298 242 299 243 SUBROUTINE bdy_dyn3d_orlanski( Kbb, puu, pvv, Kaa, idx, dta, ib_bdy, ll _npo )300 SUBROUTINE bdy_dyn3d_orlanski( Kbb, puu, pvv, Kaa, idx, dta, ib_bdy, llrim0, ll_npo ) 244 301 !!---------------------------------------------------------------------- 245 302 !! *** SUBROUTINE bdy_dyn3d_orlanski *** … … 253 310 INTEGER , INTENT( in ) :: Kbb, Kaa ! Time level indices 254 311 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT( inout ) :: puu, pvv ! Ocean velocities (to be updated at open boundaries) 255 TYPE(OBC_INDEX) , INTENT( in ) :: idx ! OBC indices 256 TYPE(OBC_DATA) , INTENT( in ) :: dta ! OBC external data 257 INTEGER , INTENT( in ) :: ib_bdy ! BDY set index 258 LOGICAL , INTENT( in ) :: ll_npo ! switch for NPO version 312 TYPE(OBC_INDEX) , INTENT( in ) :: idx ! OBC indices 313 TYPE(OBC_DATA) , INTENT( in ) :: dta ! OBC external data 314 INTEGER , INTENT( in ) :: ib_bdy ! BDY set index 315 LOGICAL , INTENT( in ) :: llrim0 ! indicate if rim 0 is treated 316 LOGICAL , INTENT( in ) :: ll_npo ! switch for NPO version 259 317 260 318 INTEGER :: jb, igrd ! dummy loop indices … … 265 323 igrd = 2 ! Orlanski bc on u-velocity; 266 324 ! 267 CALL bdy_orlanski_3d( idx, igrd, puu(:,:,:,Kbb), puu(:,:,:,Kaa), dta%u3d, ll_npo )325 CALL bdy_orlanski_3d( idx, igrd, puu(:,:,:,Kbb), puu(:,:,:,Kaa), dta%u3d, ll_npo, llrim0 ) 268 326 269 327 igrd = 3 ! Orlanski bc on v-velocity 270 328 ! 271 CALL bdy_orlanski_3d( idx, igrd, pvv(:,:,:,Kbb), pvv(:,:,:,Kaa), dta%v3d, ll_npo ) 272 ! 273 CALL lbc_bdy_lnk( 'bdydyn3d', puu(:,:,:,Kaa), 'U', -1., ib_bdy ) ! Boundary points should be updated 274 CALL lbc_bdy_lnk( 'bdydyn3d', pvv(:,:,:,Kaa), 'V', -1., ib_bdy ) 329 CALL bdy_orlanski_3d( idx, igrd, pvv(:,:,:,Kbb), pvv(:,:,:,Kaa), dta%v3d, ll_npo, llrim0 ) 275 330 ! 276 331 END SUBROUTINE bdy_dyn3d_orlanski … … 322 377 END DO 323 378 ! 324 CALL lbc_lnk_multi( 'bdydyn3d', puu(:,:,:,Krhs), 'U', -1., pvv(:,:,:,Krhs), 'V', -1. ) ! Boundary points should be updated325 !326 379 IF( ln_timing ) CALL timing_stop('bdy_dyn3d_dmp') 327 380 ! … … 329 382 330 383 331 SUBROUTINE bdy_dyn3d_nmn( puu, pvv, Kaa, idx, ib_bdy )384 SUBROUTINE bdy_dyn3d_nmn( puu, pvv, Kaa, idx, ib_bdy, llrim0 ) 332 385 !!---------------------------------------------------------------------- 333 386 !! *** SUBROUTINE bdy_dyn3d_nmn *** … … 342 395 TYPE(OBC_INDEX) , INTENT( in ) :: idx ! OBC indices 343 396 INTEGER , INTENT( in ) :: ib_bdy ! BDY set index 344 345 INTEGER :: jb, igrd ! dummy loop indices397 LOGICAL , INTENT( in ) :: llrim0 ! indicate if rim 0 is treated 398 INTEGER :: igrd ! dummy indice 346 399 !!---------------------------------------------------------------------- 347 400 ! … … 350 403 igrd = 2 ! Neumann bc on u-velocity; 351 404 ! 352 CALL bdy_nmn( idx, igrd, puu(:,:,:,Kaa) )405 CALL bdy_nmn( idx, igrd, puu(:,:,:,Kaa), llrim0 ) 353 406 354 407 igrd = 3 ! Neumann bc on v-velocity 355 408 ! 356 CALL bdy_nmn( idx, igrd, pvv(:,:,:,Kaa) ) 357 ! 358 CALL lbc_bdy_lnk( 'bdydyn3d', puu(:,:,:,Kaa), 'U', -1., ib_bdy ) ! Boundary points should be updated 359 CALL lbc_bdy_lnk( 'bdydyn3d', pvv(:,:,:,Kaa), 'V', -1., ib_bdy ) 409 CALL bdy_nmn( idx, igrd, pvv(:,:,:,Kaa), llrim0 ) 360 410 ! 361 411 END SUBROUTINE bdy_dyn3d_nmn
Note: See TracChangeset
for help on using the changeset viewer.