Changeset 11191 for NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn2d.F90
- Timestamp:
- 2019-06-27T10:14:39+02:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn2d.F90
r11149 r11191 50 50 REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pssh 51 51 !! 52 INTEGER :: ib_bdy ! Loop counter 52 INTEGER :: ib_bdy, ir ! BDY set index, rim index 53 LOGICAL :: llrim0 ! indicate if rim 0 is treated 53 54 LOGICAL, DIMENSION(4) :: llsend2, llrecv2, llsend3, llrecv3 ! indicate how communications are to be carried out 54 55 55 DO ib_bdy=1, nb_bdy 56 SELECT CASE( cn_dyn2d(ib_bdy) ) 57 CASE('none') 58 CYCLE 59 CASE('frs') 60 CALL bdy_dyn2d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d ) 61 CASE('flather') 62 CALL bdy_dyn2d_fla( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d, pssh, phur, phvr ) 63 CASE('orlanski') 64 CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, & 65 & pua2d, pva2d, pub2d, pvb2d, ll_npo=.false.) 66 CASE('orlanski_npo') 67 CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, & 68 & pua2d, pva2d, pub2d, pvb2d, ll_npo=.true. ) 69 CASE DEFAULT 70 CALL ctl_stop( 'bdy_dyn2d : unrecognised option for open boundaries for barotropic variables' ) 71 END SELECT 72 ENDDO 73 ! 74 llsend2(:) = .false. 75 llrecv2(:) = .false. 76 llsend3(:) = .false. 77 llrecv3(:) = .false. 78 DO ib_bdy=1, nb_bdy 79 SELECT CASE( cn_dyn2d(ib_bdy) ) 80 CASE('flather') 81 llsend2(1:2) = llsend2(1:2) .OR. lsend_bdyint(ib_bdy,2,1:2) ! west/east, U points 82 llsend2(1) = llsend2(1) .OR. lsend_bdyext(ib_bdy,2,1) ! neighbour might search point towards bdy on its east 83 llrecv2(1:2) = llrecv2(1:2) .OR. lrecv_bdyint(ib_bdy,2,1:2) ! west/east, U points 84 llrecv2(2) = llrecv2(2) .OR. lrecv_bdyext(ib_bdy,2,2) ! might search point towards bdy on the east 85 llsend3(3:4) = llsend3(3:4) .OR. lsend_bdyint(ib_bdy,3,3:4) ! north/south, V points 86 llsend3(3) = llsend3(3) .OR. lsend_bdyext(ib_bdy,3,3) ! neighbour might search point towards bdy on its north 87 llrecv3(3:4) = llrecv3(3:4) .OR. lrecv_bdyint(ib_bdy,3,3:4) ! north/south, V points 88 llrecv3(4) = llrecv3(4) .OR. lrecv_bdyext(ib_bdy,3,4) ! might search point towards bdy on the north 89 CASE('orlanski', 'orlanski_npo') 90 llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:) ! possibly every direction, U points 91 llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:) ! possibly every direction, U points 92 llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:) ! possibly every direction, V points 93 llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:) ! possibly every direction, V points 94 END SELECT 95 END DO 96 IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN ! if need to send/recv in at least one direction 97 CALL lbc_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, pua2d, 'U', -1. ) 98 END IF 99 IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN ! if need to send/recv in at least one direction 100 CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, pva2d, 'V', -1. ) 101 END IF 56 DO ir = 1, 0, -1 ! treat rim 1 before rim 0 57 IF( ir == 0 ) THEN ; llrim0 = .TRUE. 58 ELSE ; llrim0 = .FALSE. 59 END IF 60 DO ib_bdy=1, nb_bdy 61 SELECT CASE( cn_dyn2d(ib_bdy) ) 62 CASE('none') 63 CYCLE 64 CASE('frs') ! treat the whole boundary at once 65 IF( llrim0 ) CALL bdy_dyn2d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d ) 66 CASE('flather') 67 CALL bdy_dyn2d_fla( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d, pssh, phur, phvr, llrim0 ) 68 CASE('orlanski') 69 CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, & 70 & pua2d, pva2d, pub2d, pvb2d, llrim0, ll_npo=.false. ) 71 CASE('orlanski_npo') 72 CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, & 73 & pua2d, pva2d, pub2d, pvb2d, llrim0, ll_npo=.true. ) 74 CASE DEFAULT 75 CALL ctl_stop( 'bdy_dyn2d : unrecognised option for open boundaries for barotropic variables' ) 76 END SELECT 77 ENDDO 78 ! 79 llsend2(:) = .false. 80 llrecv2(:) = .false. 81 llsend3(:) = .false. 82 llrecv3(:) = .false. 83 DO ib_bdy=1, nb_bdy 84 SELECT CASE( cn_dyn2d(ib_bdy) ) 85 CASE('flather') 86 llsend2(1:2) = llsend2(1:2) .OR. lsend_bdyint(ib_bdy,2,1:2,ir) ! west/east, U points 87 llsend2(1) = llsend2(1) .OR. lsend_bdyext(ib_bdy,2,1,ir) ! neighbour might search point towards its east bdy 88 llrecv2(1:2) = llrecv2(1:2) .OR. lrecv_bdyint(ib_bdy,2,1:2,ir) ! west/east, U points 89 llrecv2(2) = llrecv2(2) .OR. lrecv_bdyext(ib_bdy,2,2,ir) ! might search point towards bdy on the east 90 llsend3(3:4) = llsend3(3:4) .OR. lsend_bdyint(ib_bdy,3,3:4,ir) ! north/south, V points 91 llsend3(3) = llsend3(3) .OR. lsend_bdyext(ib_bdy,3,3,ir) ! neighbour might search point towards its north bdy 92 llrecv3(3:4) = llrecv3(3:4) .OR. lrecv_bdyint(ib_bdy,3,3:4,ir) ! north/south, V points 93 llrecv3(4) = llrecv3(4) .OR. lrecv_bdyext(ib_bdy,3,4,ir) ! might search point towards bdy on the north 94 CASE('orlanski', 'orlanski_npo') 95 llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:,ir) ! possibly every direction, U points 96 llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:,ir) ! possibly every direction, U points 97 llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:,ir) ! possibly every direction, V points 98 llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:,ir) ! possibly every direction, V points 99 END SELECT 100 END DO 101 IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN ! if need to send/recv in at least one direction 102 CALL lbc_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, pua2d, 'U', -1. ) 103 END IF 104 IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN ! if need to send/recv in at least one direction 105 CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, pva2d, 'V', -1. ) 106 END IF 107 ! 108 END DO 102 109 ! 103 110 END SUBROUTINE bdy_dyn2d … … 119 126 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pua2d, pva2d 120 127 !! 121 INTEGER :: jb , jk! dummy loop indices128 INTEGER :: jb ! dummy loop indices 122 129 INTEGER :: ii, ij, igrd ! local integers 123 130 REAL(wp) :: zwgt ! boundary weight … … 143 150 144 151 145 SUBROUTINE bdy_dyn2d_fla( idx, dta, ib_bdy, pua2d, pva2d, pssh, phur, phvr )152 SUBROUTINE bdy_dyn2d_fla( idx, dta, ib_bdy, pua2d, pva2d, pssh, phur, phvr, llrim0 ) 146 153 !!---------------------------------------------------------------------- 147 154 !! *** SUBROUTINE bdy_dyn2d_fla *** … … 166 173 INTEGER, INTENT(in) :: ib_bdy ! BDY set index 167 174 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pua2d, pva2d 168 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pssh, phur, phvr 169 175 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pssh, phur, phvr 176 LOGICAL , INTENT(in) :: llrim0 ! indicate if rim 0 is treated 177 INTEGER :: ibeg, iend ! length of rim to be treated (rim 0 or rim 1) 170 178 INTEGER :: jb, igrd ! dummy loop indices 171 179 INTEGER :: ii, ij ! 2D addresses 172 180 INTEGER :: iiTrim, ijTrim ! T pts i/j-indice on the rim 173 181 INTEGER :: iiToce, ijToce, iiUoce, ijVoce ! T, U and V pts i/j-indice of the ocean next to the rim 174 REAL(wp) , POINTER :: flagu, flagv! short cuts182 REAL(wp) :: flagu, flagv ! short cuts 175 183 REAL(wp) :: zfla ! Flather correction 176 184 REAL(wp) :: z1_2 ! … … 185 193 ! Fill temporary array with ssh data (here we use spgu with the alias sshdta): 186 194 igrd = 1 195 IF( llrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(igrd) 196 ELSE ; ibeg = idx%nblenrim0(igrd)+1 ; iend = idx%nblenrim(igrd) 197 END IF 187 198 sshdta(:,:) = 0.0 188 DO jb = 1, idx%nblenrim(igrd)199 DO jb = ibeg, iend 189 200 ii = idx%nbi(jb,igrd) 190 201 ij = idx%nbj(jb,igrd) … … 196 207 END DO 197 208 ! 198 igrd = 2 ! Flather bc on u-velocity ;209 igrd = 2 ! Flather bc on u-velocity 199 210 ! ! remember that flagu=-1 if normal velocity direction is outward 200 211 ! ! I think we should rather use after ssh ? 201 DO jb = 1, idx%nblenrim(igrd) 202 ii = idx%nbi(jb,igrd) 203 ij = idx%nbj(jb,igrd) 204 flagu => idx%flagu(jb,igrd) 212 IF( llrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(igrd) 213 ELSE ; ibeg = idx%nblenrim0(igrd)+1 ; iend = idx%nblenrim(igrd) 214 END IF 215 DO jb = ibeg, iend 216 ii = idx%nbi(jb,igrd) 217 ij = idx%nbj(jb,igrd) 218 flagu = idx%flagu(jb,igrd) 205 219 IF( flagu == 0. ) THEN 206 220 pua2d(ii,ij) = dta%u2d(jb) … … 222 236 igrd = 3 ! Flather bc on v-velocity 223 237 ! ! remember that flagv=-1 if normal velocity direction is outward 224 DO jb = 1, idx%nblenrim(igrd) 225 ii = idx%nbi(jb,igrd) 226 ij = idx%nbj(jb,igrd) 227 flagv => idx%flagv(jb,igrd) 238 IF( llrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(igrd) 239 ELSE ; ibeg = idx%nblenrim0(igrd)+1 ; iend = idx%nblenrim(igrd) 240 END IF 241 DO jb = ibeg, iend 242 ii = idx%nbi(jb,igrd) 243 ij = idx%nbj(jb,igrd) 244 flagv = idx%flagv(jb,igrd) 228 245 IF( flagv == 0. ) THEN 229 246 pva2d(ii,ij) = dta%v2d(jb) … … 246 263 247 264 248 SUBROUTINE bdy_dyn2d_orlanski( idx, dta, ib_bdy, pua2d, pva2d, pub2d, pvb2d, ll _npo )265 SUBROUTINE bdy_dyn2d_orlanski( idx, dta, ib_bdy, pua2d, pva2d, pub2d, pvb2d, llrim0, ll_npo ) 249 266 !!---------------------------------------------------------------------- 250 267 !! *** SUBROUTINE bdy_dyn2d_orlanski *** … … 263 280 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pub2d, pvb2d 264 281 LOGICAL, INTENT(in) :: ll_npo ! flag for NPO version 265 282 LOGICAL, INTENT(in) :: llrim0 ! indicate if rim 0 is treated 266 283 INTEGER :: ib, igrd ! dummy loop indices 267 284 INTEGER :: ii, ij, iibm1, ijbm1 ! indices … … 270 287 igrd = 2 ! Orlanski bc on u-velocity; 271 288 ! 272 CALL bdy_orlanski_2d( idx, igrd, pub2d, pua2d, dta%u2d, ll _npo )289 CALL bdy_orlanski_2d( idx, igrd, pub2d, pua2d, dta%u2d, llrim0, ll_npo ) 273 290 274 291 igrd = 3 ! Orlanski bc on v-velocity 275 292 ! 276 CALL bdy_orlanski_2d( idx, igrd, pvb2d, pva2d, dta%v2d, ll _npo )293 CALL bdy_orlanski_2d( idx, igrd, pvb2d, pva2d, dta%v2d, llrim0, ll_npo ) 277 294 ! 278 295 END SUBROUTINE bdy_dyn2d_orlanski … … 288 305 REAL(wp), DIMENSION(jpi,jpj,1), INTENT(inout) :: zssh ! Sea level, need 3 dimensions to be used by bdy_nmn 289 306 !! 290 INTEGER :: ib_bdy ! bdy index 307 INTEGER :: ib_bdy, ir ! bdy index, rim index 308 INTEGER :: ibeg, iend ! length of rim to be treated (rim 0 or rim 1) 309 LOGICAL :: llrim0 ! indicate if rim 0 is treated 291 310 LOGICAL, DIMENSION(4) :: llsend1, llrecv1 ! indicate how communications are to be carried out 292 311 !!---------------------------------------------------------------------- 293 llsend1(:) = .false. 294 llrecv1(:) = .false. 295 DO ib_bdy = 1, nb_bdy 296 CALL bdy_nmn( idx_bdy(ib_bdy), 1, zssh ) ! zssh is masked 297 llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:) ! possibly every direction, T points 298 llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:) ! possibly every direction, T points 299 END DO 300 IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN ! if need to send/recv in at least one direction 301 CALL lbc_bdy_lnk( 'bdydyn2d', llsend1, llrecv1, zssh(:,:,1), 'T', 1. ) 302 END IF 312 DO ir = 1, 0, -1 ! treat rim 1 before rim 0 313 llsend1(:) = .false. 314 llrecv1(:) = .false. 315 IF( ir == 0 ) THEN ; llrim0 = .TRUE. 316 ELSE ; llrim0 = .FALSE. 317 END IF 318 DO ib_bdy = 1, nb_bdy 319 CALL bdy_nmn( idx_bdy(ib_bdy), 1, zssh, llrim0 ) ! zssh is masked 320 llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:,ir) ! possibly every direction, T points 321 llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:,ir) ! possibly every direction, T points 322 END DO 323 IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN ! if need to send/recv in at least one direction 324 CALL lbc_bdy_lnk( 'bdydyn2d', llsend1, llrecv1, zssh(:,:,1), 'T', 1. ) 325 END IF 326 END DO 303 327 ! 304 328 END SUBROUTINE bdy_ssh
Note: See TracChangeset
for help on using the changeset viewer.