Changeset 11191
- Timestamp:
- 2019-06-27T10:14:39+02:00 (5 years ago)
- Location:
- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE
- Files:
-
- 8 edited
-
BDY/bdy_oce.F90 (modified) (2 diffs)
-
BDY/bdydyn2d.F90 (modified) (11 diffs)
-
BDY/bdydyn3d.F90 (modified) (10 diffs)
-
BDY/bdyice.F90 (modified) (8 diffs)
-
BDY/bdyini.F90 (modified) (23 diffs)
-
BDY/bdylib.F90 (modified) (14 diffs)
-
BDY/bdytra.F90 (modified) (3 diffs)
-
DYN/dynkeg.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdy_oce.F90
r11071 r11191 22 22 INTEGER , DIMENSION(jpbgrd) :: nblen 23 23 INTEGER , DIMENSION(jpbgrd) :: nblenrim 24 INTEGER , DIMENSION(jpbgrd) :: nblenrim0 24 25 INTEGER , POINTER, DIMENSION(:,:) :: nbi 25 26 INTEGER , POINTER, DIMENSION(:,:) :: nbj … … 139 140 TYPE(OBC_DATA) , DIMENSION(jp_bdy), TARGET :: dta_bdy !: bdy external data (local process) 140 141 !$AGRIF_END_DO_NOT_TREAT 141 LOGICAL, ALLOCATABLE, DIMENSION(:,:,: ) :: lsend_bdy !: mark needed communication for given boundary, grid and neighbour142 LOGICAL, ALLOCATABLE, DIMENSION(:,:,: ) :: lrecv_bdy !: when searching in any direction143 LOGICAL, ALLOCATABLE, DIMENSION(:,:,: ) :: lsend_bdyint !: mark needed communication for given boundary, grid and neighbour144 LOGICAL, ALLOCATABLE, DIMENSION(:,:,: ) :: lrecv_bdyint !: when searching towards the interior of the computational domain145 LOGICAL, ALLOCATABLE, DIMENSION(:,:,: ) :: lsend_bdyext !: mark needed communication for given boundary, grid and neighbour146 LOGICAL, ALLOCATABLE, DIMENSION(:,:,: ) :: lrecv_bdyext !: when searching towards the exterior of the computational domain142 LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) :: lsend_bdy !: mark needed communication for given boundary, grid and neighbour 143 LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) :: lrecv_bdy !: when searching in any direction 144 LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) :: lsend_bdyint !: mark needed communication for given boundary, grid and neighbour 145 LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) :: lrecv_bdyint !: when searching towards the interior of the computational domain 146 LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) :: lsend_bdyext !: mark needed communication for given boundary, grid and neighbour 147 LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) :: lrecv_bdyext !: when searching towards the exterior of the computational domain 147 148 !!---------------------------------------------------------------------- 148 149 !! NEMO/OCE 4.0 , NEMO Consortium (2018) -
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 -
NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn3d.F90
r11071 r11191 42 42 INTEGER, INTENT(in) :: kt ! Main time step counter 43 43 ! 44 INTEGER :: ib_bdy ! loop index 44 INTEGER :: ib_bdy, ir ! BDY set index, rim index 45 LOGICAL :: llrim0 ! indicate if rim 0 is treated 45 46 LOGICAL, DIMENSION(4) :: llsend2, llrecv2, llsend3, llrecv3 ! indicate how communications are to be carried out 46 47 47 48 !!---------------------------------------------------------------------- 48 DO ib_bdy=1, nb_bdy 49 DO ir = 1, 0, -1 ! treat rim 1 before rim 0 50 IF( ir == 0 ) THEN ; llrim0 = .TRUE. 51 ELSE ; llrim0 = .FALSE. 52 END IF 53 DO ib_bdy=1, nb_bdy 54 ! 55 SELECT CASE( cn_dyn3d(ib_bdy) ) 56 CASE('none') ; CYCLE 57 CASE('frs' ) ! treat the whole boundary at once 58 IF( ir == 0) CALL bdy_dyn3d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 59 CASE('specified') ! treat the whole rim at once 60 IF( ir == 0) CALL bdy_dyn3d_spe( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 61 CASE('zero') ! treat the whole rim at once 62 IF( ir == 0) CALL bdy_dyn3d_zro( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 63 CASE('orlanski' ) ; CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, llrim0, ll_npo=.false. ) 64 CASE('orlanski_npo'); CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, llrim0, ll_npo=.true. ) 65 CASE('zerograd') ; CALL bdy_dyn3d_zgrad( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy, llrim0 ) 66 CASE('neumann') ; CALL bdy_dyn3d_nmn( idx_bdy(ib_bdy), ib_bdy, llrim0 ) 67 CASE DEFAULT ; CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' ) 68 END SELECT 69 END DO 49 70 ! 50 SELECT CASE( cn_dyn3d(ib_bdy) ) 51 CASE('none') ; CYCLE 52 CASE('frs' ) ; CALL bdy_dyn3d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 53 CASE('specified') ; CALL bdy_dyn3d_spe( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 54 CASE('zero') ; CALL bdy_dyn3d_zro( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 55 CASE('orlanski' ) ; CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.false. ) 56 CASE('orlanski_npo'); CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.true. ) 57 CASE('zerograd') ; CALL bdy_dyn3d_zgrad( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 58 CASE('neumann') ; CALL bdy_dyn3d_nmn( idx_bdy(ib_bdy), ib_bdy ) 59 CASE DEFAULT ; CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' ) 60 END SELECT 61 END DO 62 ! 63 llsend2(:) = .false. 64 llrecv2(:) = .false. 65 llsend3(:) = .false. 66 llrecv3(:) = .false. 67 DO ib_bdy=1, nb_bdy 68 SELECT CASE( cn_dyn3d(ib_bdy) ) 69 CASE('orlanski', 'orlanski_npo') 70 llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:) ! possibly every direction, U points 71 llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:) ! possibly every direction, U points 72 llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:) ! possibly every direction, V points 73 llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:) ! possibly every direction, V points 74 CASE('zerograd') 75 llsend2(3:4) = llsend2(3:4) .OR. lsend_bdyint(ib_bdy,2,3:4) ! north/south, U points 76 llrecv2(3:4) = llrecv2(3:4) .OR. lrecv_bdyint(ib_bdy,2,3:4) ! north/south, U points 77 llsend3(1:2) = llsend3(1:2) .OR. lsend_bdyint(ib_bdy,3,1:2) ! west/east, V points 78 llrecv3(1:2) = llrecv3(1:2) .OR. lrecv_bdyint(ib_bdy,3,1:2) ! west/east, V points 79 CASE('neumann') 80 llsend2(:) = llsend2(:) .OR. lsend_bdyint(ib_bdy,2,:) ! possibly every direction, U points 81 llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(ib_bdy,2,:) ! possibly every direction, U points 82 llsend3(:) = llsend3(:) .OR. lsend_bdyint(ib_bdy,3,:) ! possibly every direction, V points 83 llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(ib_bdy,3,:) ! possibly every direction, V points 84 END SELECT 85 END DO 86 ! 87 IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN ! if need to send/recv in at least one direction 88 CALL lbc_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, ua, 'U', -1. ) 89 END IF 90 IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN ! if need to send/recv in at least one direction 91 CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, va, 'V', -1. ) 92 END IF 71 llsend2(:) = .false. 72 llrecv2(:) = .false. 73 llsend3(:) = .false. 74 llrecv3(:) = .false. 75 DO ib_bdy=1, nb_bdy 76 SELECT CASE( cn_dyn3d(ib_bdy) ) 77 CASE('orlanski', 'orlanski_npo') 78 llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:,ir) ! possibly every direction, U points 79 llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:,ir) ! possibly every direction, U points 80 llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:,ir) ! possibly every direction, V points 81 llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:,ir) ! possibly every direction, V points 82 CASE('zerograd') 83 llsend2(3:4) = llsend2(3:4) .OR. lsend_bdyint(ib_bdy,2,3:4,ir) ! north/south, U points 84 llrecv2(3:4) = llrecv2(3:4) .OR. lrecv_bdyint(ib_bdy,2,3:4,ir) ! north/south, U points 85 llsend3(1:2) = llsend3(1:2) .OR. lsend_bdyint(ib_bdy,3,1:2,ir) ! west/east, V points 86 llrecv3(1:2) = llrecv3(1:2) .OR. lrecv_bdyint(ib_bdy,3,1:2,ir) ! west/east, V points 87 CASE('neumann') 88 llsend2(:) = llsend2(:) .OR. lsend_bdyint(ib_bdy,2,:,ir) ! possibly every direction, U points 89 llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(ib_bdy,2,:,ir) ! possibly every direction, U points 90 llsend3(:) = llsend3(:) .OR. lsend_bdyint(ib_bdy,3,:,ir) ! possibly every direction, V points 91 llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(ib_bdy,3,:,ir) ! possibly every direction, V points 92 END SELECT 93 END DO 94 ! 95 IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN ! if need to send/recv in at least one direction 96 CALL lbc_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, ua, 'U', -1. ) 97 END IF 98 IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN ! if need to send/recv in at least one direction 99 CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, va, 'V', -1. ) 100 END IF 101 END DO 93 102 ! 94 103 END SUBROUTINE bdy_dyn3d … … 133 142 134 143 135 SUBROUTINE bdy_dyn3d_zgrad( idx, dta, kt , ib_bdy)144 SUBROUTINE bdy_dyn3d_zgrad( idx, dta, kt, ib_bdy, llrim0 ) 136 145 !!---------------------------------------------------------------------- 137 146 !! *** SUBROUTINE bdy_dyn3d_zgrad *** … … 141 150 !!---------------------------------------------------------------------- 142 151 INTEGER :: kt 143 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 144 TYPE(OBC_DATA), INTENT(in) :: dta ! OBC external data 145 INTEGER, INTENT(in) :: ib_bdy ! BDY set index 152 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 153 TYPE(OBC_DATA), INTENT(in) :: dta ! OBC external data 154 INTEGER, INTENT(in) :: ib_bdy ! BDY set index 155 LOGICAL, INTENT(in) :: llrim0 ! indicate if rim 0 is treated 146 156 !! 147 157 INTEGER :: jb, jk ! dummy loop indices 148 158 INTEGER :: ii, ij, igrd ! local integers 149 159 INTEGER :: flagu, flagv ! short cuts 160 INTEGER :: ibeg, iend ! length of rim to be treated (rim 0 or rim 1 or both) 150 161 !!---------------------------------------------------------------------- 151 162 ! 152 163 igrd = 2 ! Copying tangential velocity into bdy points 153 DO jb = 1, idx%nblenrim(igrd) 164 IF( llrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(igrd) 165 ELSE ; ibeg = idx%nblenrim0(igrd)+1 ; iend = idx%nblenrim(igrd) 166 ENDIF 167 DO jb = ibeg, iend 154 168 ii = idx%nbi(jb,igrd) 155 169 ij = idx%nbj(jb,igrd) … … 168 182 ! 169 183 igrd = 3 ! Copying tangential velocity into bdy points 170 DO jb = 1, idx%nblenrim(igrd) 184 IF( llrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(igrd) 185 ELSE ; ibeg = idx%nblenrim0(igrd)+1 ; iend = idx%nblenrim(igrd) 186 ENDIF 187 DO jb = ibeg, iend 171 188 ii = idx%nbi(jb,igrd) 172 189 ij = idx%nbj(jb,igrd) … … 268 285 269 286 270 SUBROUTINE bdy_dyn3d_orlanski( idx, dta, ib_bdy, ll _npo )287 SUBROUTINE bdy_dyn3d_orlanski( idx, dta, ib_bdy, llrim0, ll_npo ) 271 288 !!---------------------------------------------------------------------- 272 289 !! *** SUBROUTINE bdy_dyn3d_orlanski *** … … 280 297 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 281 298 TYPE(OBC_DATA), INTENT(in) :: dta ! OBC external data 282 INTEGER, INTENT(in) :: ib_bdy ! BDY set index 283 LOGICAL, INTENT(in) :: ll_npo ! switch for NPO version 299 INTEGER, INTENT(in) :: ib_bdy ! BDY set index 300 LOGICAL, INTENT(in) :: llrim0 ! indicate if rim 0 is treated 301 LOGICAL, INTENT(in) :: ll_npo ! switch for NPO version 284 302 285 303 INTEGER :: jb, igrd ! dummy loop indices … … 290 308 igrd = 2 ! Orlanski bc on u-velocity; 291 309 ! 292 CALL bdy_orlanski_3d( idx, igrd, ub, ua, dta%u3d, ll_npo )310 CALL bdy_orlanski_3d( idx, igrd, ub, ua, dta%u3d, ll_npo, llrim0 ) 293 311 294 312 igrd = 3 ! Orlanski bc on v-velocity 295 313 ! 296 CALL bdy_orlanski_3d( idx, igrd, vb, va, dta%v3d, ll_npo )314 CALL bdy_orlanski_3d( idx, igrd, vb, va, dta%v3d, ll_npo, llrim0 ) 297 315 ! 298 316 END SUBROUTINE bdy_dyn3d_orlanski … … 347 365 348 366 349 SUBROUTINE bdy_dyn3d_nmn( idx, ib_bdy )367 SUBROUTINE bdy_dyn3d_nmn( idx, ib_bdy, llrim0 ) 350 368 !!---------------------------------------------------------------------- 351 369 !! *** SUBROUTINE bdy_dyn3d_nmn *** … … 356 374 !! 357 375 !!---------------------------------------------------------------------- 358 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 359 INTEGER, INTENT(in) :: ib_bdy ! BDY set index 360 INTEGER :: igrd ! dummy indice 376 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 377 INTEGER, INTENT(in) :: ib_bdy ! BDY set index 378 LOGICAL, INTENT(in) :: llrim0 ! indicate if rim 0 is treated 379 INTEGER :: igrd ! dummy indice 361 380 !!---------------------------------------------------------------------- 362 381 ! … … 365 384 igrd = 2 ! Neumann bc on u-velocity; 366 385 ! 367 CALL bdy_nmn( idx, igrd, ua ) ! ua is masked386 CALL bdy_nmn( idx, igrd, ua, llrim0 ) ! ua is masked 368 387 369 388 igrd = 3 ! Neumann bc on v-velocity 370 389 ! 371 CALL bdy_nmn( idx, igrd, va ) ! va is masked390 CALL bdy_nmn( idx, igrd, va, llrim0 ) ! va is masked 372 391 ! 373 392 END SUBROUTINE bdy_dyn3d_nmn -
NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdyice.F90
r11071 r11191 55 55 INTEGER, INTENT(in) :: kt ! Main time step counter 56 56 ! 57 INTEGER :: jbdy ! BDY set index 58 LOGICAL, DIMENSION(4) :: llsend1, llrecv1 ! indicate how communications are to be carried out 57 INTEGER :: jbdy, ir ! BDY set index, rim index 58 INTEGER :: ibeg, iend ! length of rim to be treated (rim 0 or rim 1) 59 LOGICAL :: llrim0 ! indicate if rim 0 is treated 60 LOGICAL, DIMENSION(4) :: llsend1, llrecv1 ! indicate how communications are to be carried out 59 61 !!---------------------------------------------------------------------- 60 62 ! controls … … 64 66 CALL ice_var_glo2eqv 65 67 ! 66 DO jbdy = 1, nb_bdy 68 DO ir = 1, 0, -1 ! treat rim 1 before rim 0 69 IF( ir == 0 ) THEN ; llrim0 = .TRUE. 70 ELSE ; llrim0 = .FALSE. 71 END IF 72 DO jbdy = 1, nb_bdy 73 ! 74 SELECT CASE( cn_ice(jbdy) ) 75 CASE('none') ; CYCLE 76 CASE('frs' ) ; CALL bdy_ice_frs( idx_bdy(jbdy), dta_bdy(jbdy), kt, jbdy, llrim0 ) 77 CASE DEFAULT 78 CALL ctl_stop( 'bdy_ice : unrecognised option for open boundaries for ice fields' ) 79 END SELECT 80 ! 81 END DO 67 82 ! 68 SELECT CASE( cn_ice(jbdy) ) 69 CASE('none') ; CYCLE 70 CASE('frs' ) ; CALL bdy_ice_frs( idx_bdy(jbdy), dta_bdy(jbdy), kt, jbdy ) 71 CASE DEFAULT 72 CALL ctl_stop( 'bdy_ice : unrecognised option for open boundaries for ice fields' ) 73 END SELECT 74 ! 75 END DO 76 ! 77 ! Update bdy points 78 llsend1(:) = .false. 79 llrecv1(:) = .false. 80 DO jbdy = 1, nb_bdy 81 IF( cn_ice(jbdy) == 'frs' ) THEN 82 llsend1(:) = llsend1(:) .OR. lsend_bdyint(jbdy,1,:) ! possibly every direction, T points 83 llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(jbdy,1,:) ! possibly every direction, T points 83 ! Update bdy points 84 llsend1(:) = .false. 85 llrecv1(:) = .false. 86 DO jbdy = 1, nb_bdy 87 IF( cn_ice(jbdy) == 'frs' ) THEN 88 llsend1(:) = llsend1(:) .OR. lsend_bdyint(jbdy,1,:,ir) ! possibly every direction, T points 89 llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(jbdy,1,:,ir) ! possibly every direction, T points 90 END IF 91 END DO ! jbdy 92 IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN ! if need to send/recv in at least one direction 93 ! exchange 3d arrays 94 CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, a_i , 'T', 1., h_i , 'T', 1., h_s , 'T', 1. & 95 & , oa_i, 'T', 1., a_ip, 'T', 1., v_ip, 'T', 1. & 96 & , s_i , 'T', 1., t_su, 'T', 1., v_i , 'T', 1. & 97 & , v_s , 'T', 1., sv_i, 'T', 1. ) 98 ! exchange 4d arrays 99 CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_s , 'T', 1., e_s , 'T', 1. ) ! third dimension = 1 100 CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_i , 'T', 1., e_i , 'T', 1. ) ! third dimension = jpk 84 101 END IF 85 END DO 86 IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN ! if need to send/recv in at least one direction 87 ! exchange 3d arrays 88 CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, a_i , 'T', 1., h_i , 'T', 1., h_s , 'T', 1. & 89 & , oa_i, 'T', 1., a_ip, 'T', 1., v_ip, 'T', 1. & 90 & , s_i , 'T', 1., t_su, 'T', 1., v_i , 'T', 1. & 91 & , v_s , 'T', 1., sv_i, 'T', 1. ) 92 ! exchange 4d arrays 93 CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_s , 'T', 1., e_s , 'T', 1. ) ! third dimension = 1 94 CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_i , 'T', 1., e_i , 'T', 1. ) ! third dimension = jpk 95 END IF 102 END DO ! ir 96 103 ! 97 104 CALL ice_cor( kt , 0 ) ! -- In case categories are out of bounds, do a remapping … … 108 115 109 116 110 SUBROUTINE bdy_ice_frs( idx, dta, kt, jbdy )117 SUBROUTINE bdy_ice_frs( idx, dta, kt, jbdy, llrim0 ) 111 118 !!------------------------------------------------------------------------------ 112 119 !! *** SUBROUTINE bdy_ice_frs *** … … 117 124 !! dimensional baroclinic ocean model with realistic topography. Tellus, 365-382. 118 125 !!------------------------------------------------------------------------------ 119 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 120 TYPE(OBC_DATA), INTENT(in) :: dta ! OBC external data 121 INTEGER, INTENT(in) :: kt ! main time-step counter 122 INTEGER, INTENT(in) :: jbdy ! BDY set index 126 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 127 TYPE(OBC_DATA), INTENT(in) :: dta ! OBC external data 128 INTEGER, INTENT(in) :: kt ! main time-step counter 129 INTEGER, INTENT(in) :: jbdy ! BDY set index 130 LOGICAL, INTENT(in) :: llrim0 ! indicate if rim 0 is treated 123 131 ! 124 132 INTEGER :: jpbound ! 0 = incoming ice 125 133 ! ! 1 = outgoing ice 134 INTEGER :: ibeg, iend ! length of rim to be treated (rim 0 or rim 1) 126 135 INTEGER :: i_bdy, jgrd ! dummy loop indices 127 136 INTEGER :: ji, jj, jk, jl, ib, jb … … 132 141 ! 133 142 jgrd = 1 ! Everything is at T-points here 143 IF( llrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(jgrd) 144 ELSE ; ibeg = idx%nblenrim0(jgrd)+1 ; iend = idx%nblenrim(jgrd) 145 END IF 134 146 ! 135 147 DO jl = 1, jpl 136 DO i_bdy = 1, idx%nblenrim(jgrd)148 DO i_bdy = ibeg, iend 137 149 ji = idx%nbi(i_bdy,jgrd) 138 150 jj = idx%nbj(i_bdy,jgrd) … … 162 174 163 175 DO jl = 1, jpl 164 DO i_bdy = 1, idx%nblenrim(jgrd)176 DO i_bdy = ibeg, iend 165 177 ji = idx%nbi(i_bdy,jgrd) 166 178 jj = idx%nbj(i_bdy,jgrd) … … 295 307 CHARACTER(len=1), INTENT(in) :: cd_type ! nature of velocity grid-points 296 308 ! 297 INTEGER :: i_bdy, jgrd ! dummy loop indices 298 INTEGER :: ji, jj ! local scalar 299 INTEGER :: jbdy ! BDY set index 309 INTEGER :: i_bdy, jgrd ! dummy loop indices 310 INTEGER :: ji, jj ! local scalar 311 INTEGER :: jbdy, ir ! BDY set index, rim index 312 INTEGER :: ibeg, iend ! length of rim to be treated (rim 0 or rim 1) 300 313 REAL(wp) :: zmsk1, zmsk2, zflag 301 314 LOGICAL, DIMENSION(4) :: llsend2, llrecv2, llsend3, llrecv3 ! indicate how communications are to be carried out … … 303 316 IF( ln_timing ) CALL timing_start('bdy_ice_dyn') 304 317 ! 305 DO jbdy=1, nb_bdy 318 DO ir = 1, 0, -1 319 DO jbdy = 1, nb_bdy 320 ! 321 SELECT CASE( cn_ice(jbdy) ) 322 ! 323 CASE('none') 324 CYCLE 325 ! 326 CASE('frs') 327 ! 328 IF( nn_ice_dta(jbdy) == 0 ) CYCLE ! case ice boundaries = initial conditions 329 ! ! do not change ice velocity (it is only computed by rheology) 330 SELECT CASE ( cd_type ) 331 ! 332 CASE ( 'U' ) 333 jgrd = 2 ! u velocity 334 IF( ir == 0 ) THEN ; ibeg = 1 ; iend = idx_bdy(jbdy)%nblenrim0(jgrd) 335 ELSE ; ibeg = idx_bdy(jbdy)%nblenrim0(jgrd)+1 ; iend = idx_bdy(jbdy)%nblenrim(jgrd) 336 END IF 337 DO i_bdy = ibeg, iend 338 ji = idx_bdy(jbdy)%nbi(i_bdy,jgrd) 339 jj = idx_bdy(jbdy)%nbj(i_bdy,jgrd) 340 zflag = idx_bdy(jbdy)%flagu(i_bdy,jgrd) 341 ! i-1 i i | ! i i i+1 | ! i i i+1 | 342 ! > ice > | ! o > ice | ! o > o | 343 ! => set at u_ice(i-1) ! => set to O ! => unchanged 344 IF( zflag == -1. .AND. ji > 1 .AND. ji < jpi ) THEN 345 IF ( vt_i(ji ,jj) > 0. ) THEN ; u_ice(ji,jj) = u_ice(ji-1,jj) 346 ELSEIF( vt_i(ji+1,jj) > 0. ) THEN ; u_ice(ji,jj) = 0._wp 347 END IF 348 END IF 349 ! | i i+1 i+1 ! | i i i+1 ! | i i i+1 350 ! | > ice > ! | ice > o ! | o > o 351 ! => set at u_ice(i+1) ! => set to O ! => unchanged 352 IF( zflag == 1. .AND. ji+1 < jpi+1 ) THEN 353 IF ( vt_i(ji+1,jj) > 0. ) THEN ; u_ice(ji,jj) = u_ice(ji+1,jj) 354 ELSEIF( vt_i(ji ,jj) > 0. ) THEN ; u_ice(ji,jj) = 0._wp 355 END IF 356 END IF 357 ! 358 IF( zflag == 0. ) u_ice(ji,jj) = 0._wp ! u_ice = 0 if north/south bdy 359 ! 360 END DO 361 ! 362 CASE ( 'V' ) 363 jgrd = 3 ! v velocity 364 IF( ir == 0 ) THEN ; ibeg = 1 ; iend = idx_bdy(jbdy)%nblenrim0(jgrd) 365 ELSE ; ibeg = idx_bdy(jbdy)%nblenrim0(jgrd)+1 ; iend = idx_bdy(jbdy)%nblenrim(jgrd) 366 END IF 367 DO i_bdy = ibeg, iend 368 ji = idx_bdy(jbdy)%nbi(i_bdy,jgrd) 369 jj = idx_bdy(jbdy)%nbj(i_bdy,jgrd) 370 zflag = idx_bdy(jbdy)%flagv(i_bdy,jgrd) 371 ! ¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨ ! ¨¨¨¨ïce¨¨¨(jj+1)¨¨ ! ¨¨¨¨¨¨ö¨¨¨¨(jj+1) 372 ! ^ (jj ) ! ^ (jj ) ! ^ (jj ) 373 ! ice (jj ) ! o (jj ) ! o (jj ) 374 ! ^ (jj-1) ! ! 375 ! => set to u_ice(jj-1) ! => set to 0 ! => unchanged 376 IF( zflag == -1. .AND. jj > 1 .AND. jj < jpj ) THEN 377 IF ( vt_i(ji,jj ) > 0. ) THEN ; v_ice(ji,jj) = v_ice(ji,jj-1) 378 ELSEIF( vt_i(ji,jj+1) > 0. ) THEN ; v_ice(ji,jj) = 0._wp 379 END IF 380 END IF 381 ! ^ (jj+1) ! ! 382 ! ice (jj+1) ! o (jj+1) ! o (jj+1) 383 ! ^ (jj ) ! ^ (jj ) ! ^ (jj ) 384 ! ________________ ! ____ice___(jj )_ ! _____o____(jj ) 385 ! => set to u_ice(jj+1) ! => set to 0 ! => unchanged 386 IF( zflag == 1. .AND. jj < jpj ) THEN 387 IF ( vt_i(ji,jj+1) > 0. ) THEN ; v_ice(ji,jj) = v_ice(ji,jj+1) 388 ELSEIF( vt_i(ji,jj ) > 0. ) THEN ; v_ice(ji,jj) = 0._wp 389 END IF 390 END IF 391 ! 392 IF( zflag == 0. ) v_ice(ji,jj) = 0._wp ! v_ice = 0 if west/east bdy 393 ! 394 END DO 395 ! 396 END SELECT 397 ! 398 CASE DEFAULT 399 CALL ctl_stop( 'bdy_ice_dyn : unrecognised option for open boundaries for ice fields' ) 400 END SELECT 401 ! 402 END DO ! jbdy 306 403 ! 307 SELECT CASE( cn_ice(jbdy) ) 308 ! 309 CASE('none') 310 CYCLE 311 ! 312 CASE('frs') 313 ! 314 IF( nn_ice_dta(jbdy) == 0 ) CYCLE ! case ice boundaries = initial conditions 315 ! ! do not change ice velocity (it is only computed by rheology) 316 SELECT CASE ( cd_type ) 317 ! 318 CASE ( 'U' ) 319 jgrd = 2 ! u velocity 320 DO i_bdy = 1, idx_bdy(jbdy)%nblenrim(jgrd) 321 ji = idx_bdy(jbdy)%nbi(i_bdy,jgrd) 322 jj = idx_bdy(jbdy)%nbj(i_bdy,jgrd) 323 zflag = idx_bdy(jbdy)%flagu(i_bdy,jgrd) 324 ! i-1 i i | ! i i i+1 | ! i i i+1 | 325 ! > ice > | ! o > ice | ! o > o | 326 ! => set at u_ice(i-1) ! => set to O ! => unchanged 327 IF( zflag == -1. .AND. ji > 1 .AND. ji < jpi ) THEN 328 IF ( vt_i(ji ,jj) > 0. ) THEN ; u_ice(ji,jj) = u_ice(ji-1,jj) 329 ELSEIF( vt_i(ji+1,jj) > 0. ) THEN ; u_ice(ji,jj) = 0._wp 330 END IF 331 END IF 332 ! | i i+1 i+1 ! | i i i+1 ! | i i i+1 333 ! | > ice > ! | ice > o ! | o > o 334 ! => set at u_ice(i+1) ! => set to O ! => unchanged 335 IF( zflag == 1. .AND. ji+1 < jpi+1 ) THEN 336 IF ( vt_i(ji+1,jj) > 0. ) THEN ; u_ice(ji,jj) = u_ice(ji+1,jj) 337 ELSEIF( vt_i(ji ,jj) > 0. ) THEN ; u_ice(ji,jj) = 0._wp 338 END IF 339 END IF 340 ! 341 IF( zflag == 0. ) u_ice(ji,jj) = 0._wp ! u_ice = 0 if north/south bdy 342 ! 343 END DO 344 ! 345 CASE ( 'V' ) 346 jgrd = 3 ! v velocity 347 DO i_bdy = 1, idx_bdy(jbdy)%nblenrim(jgrd) 348 ji = idx_bdy(jbdy)%nbi(i_bdy,jgrd) 349 jj = idx_bdy(jbdy)%nbj(i_bdy,jgrd) 350 zflag = idx_bdy(jbdy)%flagv(i_bdy,jgrd) 351 ! ¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨ ! ¨¨¨¨ïce¨¨¨(jj+1)¨¨ ! ¨¨¨¨¨¨ö¨¨¨¨(jj+1) 352 ! ^ (jj ) ! ^ (jj ) ! ^ (jj ) 353 ! ice (jj ) ! o (jj ) ! o (jj ) 354 ! ^ (jj-1) ! ! 355 ! => set to u_ice(jj-1) ! => set to 0 ! => unchanged 356 IF( zflag == -1. .AND. jj > 1 .AND. jj < jpj ) THEN 357 IF ( vt_i(ji,jj ) > 0. ) THEN ; v_ice(ji,jj) = v_ice(ji,jj-1) 358 ELSEIF( vt_i(ji,jj+1) > 0. ) THEN ; v_ice(ji,jj) = 0._wp 359 END IF 360 END IF 361 ! ^ (jj+1) ! ! 362 ! ice (jj+1) ! o (jj+1) ! o (jj+1) 363 ! ^ (jj ) ! ^ (jj ) ! ^ (jj ) 364 ! ________________ ! ____ice___(jj )_ ! _____o____(jj ) 365 ! => set to u_ice(jj+1) ! => set to 0 ! => unchanged 366 IF( zflag == 1. .AND. jj < jpj ) THEN 367 IF ( vt_i(ji,jj+1) > 0. ) THEN ; v_ice(ji,jj) = v_ice(ji,jj+1) 368 ELSEIF( vt_i(ji,jj ) > 0. ) THEN ; v_ice(ji,jj) = 0._wp 369 END IF 370 END IF 371 ! 372 IF( zflag == 0. ) v_ice(ji,jj) = 0._wp ! v_ice = 0 if west/east bdy 373 ! 374 END DO 375 ! 376 END SELECT 377 ! 378 CASE DEFAULT 379 CALL ctl_stop( 'bdy_ice_dyn : unrecognised option for open boundaries for ice fields' ) 404 SELECT CASE ( cd_type ) 405 CASE ( 'U' ) 406 llsend2(:) = .false. ; llrecv2(:) = .false. 407 DO jbdy = 1, nb_bdy 408 IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN 409 llsend2(:) = llsend2(:) .OR. lsend_bdyint(jbdy,2,:,ir) ! possibly every direction, U points 410 llsend2(1) = llsend2(1) .OR. lsend_bdyext(jbdy,2,1,ir) ! neighbour might search point towards its west bdy 411 llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(jbdy,2,:,ir) ! possibly every direction, U points 412 llrecv2(2) = llrecv2(2) .OR. lrecv_bdyext(jbdy,2,2,ir) ! might search point towards east bdy 413 END IF 414 END DO 415 IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN ! if need to send/recv in at least one direction 416 CALL lbc_bdy_lnk( 'bdyice', llsend2, llrecv2, u_ice, 'U', -1. ) 417 END IF 418 CASE ( 'V' ) 419 llsend3(:) = .false. ; llrecv3(:) = .false. 420 DO jbdy = 1, nb_bdy 421 IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN 422 llsend3(:) = llsend3(:) .OR. lsend_bdyint(jbdy,3,:,ir) ! possibly every direction, V points 423 llsend3(3) = llsend3(3) .OR. lsend_bdyext(jbdy,3,3,ir) ! neighbour might search point towards its south bdy 424 llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(jbdy,3,:,ir) ! possibly every direction, V points 425 llrecv3(4) = llrecv3(4) .OR. lrecv_bdyext(jbdy,3,4,ir) ! might search point towards north bdy 426 END IF 427 END DO 428 IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN ! if need to send/recv in at least one direction 429 CALL lbc_bdy_lnk( 'bdyice', llsend3, llrecv3, v_ice, 'V', -1. ) 430 END IF 380 431 END SELECT 381 ! 382 END DO ! jbdy 383 ! 384 SELECT CASE ( cd_type ) 385 CASE ( 'U' ) 386 llsend2(:) = .false. ; llrecv2(:) = .false. 387 DO jbdy = 1, nb_bdy 388 IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN 389 llsend2(:) = llsend2(:) .OR. lsend_bdyint(jbdy,2,:) ! possibly every direction, U points 390 llsend2(1) = llsend2(1) .OR. lsend_bdyext(jbdy,2,1) ! neighbour might search point towards bdy on its east 391 llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(jbdy,2,:) ! possibly every direction, U points 392 llrecv2(2) = llrecv2(2) .OR. lrecv_bdyext(jbdy,2,2) ! might search point towards bdy on the east 393 END IF 394 END DO 395 IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN ! if need to send/recv in at least one direction 396 CALL lbc_bdy_lnk( 'bdyice', llsend2, llrecv2, u_ice, 'U', -1. ) 397 END IF 398 CASE ( 'V' ) 399 llsend3(:) = .false. ; llrecv3(:) = .false. 400 DO jbdy = 1, nb_bdy 401 IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN 402 llsend3(:) = llsend3(:) .OR. lsend_bdyint(jbdy,3,:) ! possibly every direction, V points 403 llsend3(3) = llsend3(3) .OR. lsend_bdyext(jbdy,3,3) ! neighbour might search point towards bdy on its north 404 llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(jbdy,3,:) ! possibly every direction, V points 405 llrecv3(4) = llrecv3(4) .OR. lrecv_bdyext(jbdy,3,4) ! might search point towards bdy on the north 406 END IF 407 END DO 408 IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN ! if need to send/recv in at least one direction 409 CALL lbc_bdy_lnk( 'bdyice', llsend3, llrecv3, v_ice, 'V', -1. ) 410 END IF 411 END SELECT 432 END DO ! ir 412 433 ! 413 434 IF( ln_timing ) CALL timing_stop('bdy_ice_dyn') -
NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdyini.F90
r11183 r11191 33 33 PRIVATE 34 34 35 PUBLIC bdy_init ! routine called in nemo_init35 PUBLIC bdy_init ! routine called in nemo_init 36 36 PUBLIC find_neib ! routine called in bdy_nmn 37 37 … … 126 126 !! ** Input : bdy_init.nc, input file for unstructured open boundaries 127 127 !!---------------------------------------------------------------------- 128 INTEGER :: ib_bdy, ii, ij, ik, igrd, ib, ir, iseg ! dummy loop indices 129 INTEGER :: icount, icountr, ibr_max, ilen1, ibm1 ! local integers 128 INTEGER :: ib_bdy, ii, ij, igrd, ib, ir, iseg ! dummy loop indices 129 INTEGER :: icount, icountr, icountr0, ibr_max ! local integers 130 INTEGER :: ilen1, ibm1 ! - - 130 131 INTEGER :: iwe, ies, iso, ino, inum, id_dummy ! - - 131 INTEGER :: igrd_start, igrd_end, jpbdta ! - - 132 INTEGER :: jpbdtau, jpbdtas ! - - 132 INTEGER :: jpbdta, jpbdtau, jpbdtas ! - - 133 133 INTEGER :: ib_bdy1, ib_bdy2, ib1, ib2 ! - - 134 INTEGER :: i_offset, j_offset, inn ! - -135 134 INTEGER :: ii1, ii2, ii3, ij1, ij2, ij3 ! - - 136 135 INTEGER :: iibe, ijbe, iibi, ijbi ! - - 137 136 INTEGER :: flagu, flagv ! short cuts 138 INTEGER , POINTER :: nbi, nbj, nbr ! short cuts139 137 REAL(wp), POINTER, DIMENSION(:,:) :: pmask ! pointer to 2D mask fields 140 138 REAL(wp) :: zefl, zwfl, znfl, zsfl ! local scalars … … 145 143 CHARACTER(LEN=1),DIMENSION(jpbgrd) :: cgrid 146 144 INTEGER :: jpk_max 147 LOGICAL :: llsend_ea, llsend_we, llsend_so, llsend_no ! Flags for boundaries sending 148 LOGICAL :: llrecv_ea, llrecv_we, llrecv_so, llrecv_no ! Flags for boundaries receiving 149 REAL(wp), TARGET, DIMENSION(jpi,jpj) :: zfmask ! temporary fmask array excluding coastal boundary condition (shlat) 150 REAL(wp) , DIMENSION(jpi,jpj) :: ztmp 151 LOGICAL :: llnon, llson, llean, llwen ! local logicals 145 REAL(wp), DIMENSION(jpi,jpj) :: zfmask ! temporary fmask array excluding coastal boundary condition (shlat) 146 REAL(wp), DIMENSION(jpi,jpj) :: ztmask, zumask, zvmask ! temporary u/v mask array 152 147 !! 153 148 CHARACTER(LEN=1) :: ctypebdy ! - - … … 798 793 ! 799 794 iwe = mig(1) 800 ies = mig( nlci)795 ies = mig(jpi) 801 796 iso = mjg(1) 802 ino = mjg( nlcj)797 ino = mjg(jpj) 803 798 ! 804 799 DO ib_bdy = 1, nb_bdy 805 800 DO igrd = 1, jpbgrd 806 icount = 0 ! initialization of local bdy length 807 icountr = 0 ! initialization of local rim bdy length 808 idx_bdy(ib_bdy)%nblen(igrd) = 0 809 idx_bdy(ib_bdy)%nblenrim(igrd) = 0 801 icount = 0 ! initialization of local bdy length 802 icountr = 0 ! initialization of local rim 0 and rim 1 bdy length 803 icountr0 = 0 ! initialization of local rim 0 bdy length 804 idx_bdy(ib_bdy)%nblen(igrd) = 0 805 idx_bdy(ib_bdy)%nblenrim(igrd) = 0 806 idx_bdy(ib_bdy)%nblenrim0(igrd) = 0 810 807 DO ib = 1, nblendta(igrd,ib_bdy) 811 808 ! check that data is in correct order in file … … 823 820 ! 824 821 icount = icount + 1 825 IF( nbrdta(ib,igrd,ib_bdy) == 1 ) icountr = icountr + 1 822 IF( nbrdta(ib,igrd,ib_bdy) == 1 .OR. nbrdta(ib,igrd,ib_bdy) == 0 ) icountr = icountr + 1 823 IF( nbrdta(ib,igrd,ib_bdy) == 0 ) icountr0 = icountr0 + 1 826 824 ENDIF 827 825 END DO 828 idx_bdy(ib_bdy)%nblenrim(igrd) = icountr !: length of rim boundary data on each proc 829 idx_bdy(ib_bdy)%nblen (igrd) = icount !: length of boundary data on each proc 826 idx_bdy(ib_bdy)%nblen (igrd) = icount !: length of boundary data on each proc 827 idx_bdy(ib_bdy)%nblenrim (igrd) = icountr !: length of rim 0 and rim 1 boundary data on each proc 828 idx_bdy(ib_bdy)%nblenrim0(igrd) = icountr0 !: length of rim 0 boundary data on each proc 830 829 END DO ! igrd 831 830 … … 849 848 icount = 0 850 849 ! Outer loop on rimwidth to ensure outermost points come first in the local arrays. 851 DO ir =1, nn_rimwidth(ib_bdy)850 DO ir = 0, nn_rimwidth(ib_bdy) 852 851 DO ib = 1, nblendta(igrd,ib_bdy) 853 852 ! check if point is in local domain and equals ir … … 870 869 ! Initialize array indicating communications in bdy 871 870 ! ------------------------------------------------- 872 ALLOCATE( lsend_bdy(nb_bdy,jpbgrd,4 ), lrecv_bdy(nb_bdy,jpbgrd,4) )873 lsend_bdy(:,:,: ) = .false.874 lrecv_bdy(:,:,: ) = .false.871 ALLOCATE( lsend_bdy(nb_bdy,jpbgrd,4,0:1), lrecv_bdy(nb_bdy,jpbgrd,4,0:1) ) 872 lsend_bdy(:,:,:,:) = .false. 873 lrecv_bdy(:,:,:,:) = .false. 875 874 876 875 DO ib_bdy = 1, nb_bdy … … 879 878 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 880 879 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 880 IF( ib .LE. idx_bdy(ib_bdy)%nblenrim0(igrd) ) THEN ; ir = 0 881 ELSE ; ir = 1 882 END IF 881 883 ! 882 884 ! check if point has to be sent to a neighbour 883 885 ! W neighbour and on the inner left side 884 IF( ii == 2 .and. (nbondi == 0 .or. nbondi == 1) ) lsend_bdy(ib_bdy,igrd,1 ) = .true.886 IF( ii == 2 .and. (nbondi == 0 .or. nbondi == 1) ) lsend_bdy(ib_bdy,igrd,1,ir) = .true. 885 887 ! E neighbour and on the inner right side 886 IF( ii == jpi-1 .and. (nbondi == 0 .or. nbondi == -1) ) lsend_bdy(ib_bdy,igrd,2 ) = .true.888 IF( ii == jpi-1 .and. (nbondi == 0 .or. nbondi == -1) ) lsend_bdy(ib_bdy,igrd,2,ir) = .true. 887 889 ! S neighbour and on the inner down side 888 IF( ij == 2 .and. (nbondj == 0 .or. nbondj == 1) ) lsend_bdy(ib_bdy,igrd,3 ) = .true.890 IF( ij == 2 .and. (nbondj == 0 .or. nbondj == 1) ) lsend_bdy(ib_bdy,igrd,3,ir) = .true. 889 891 ! N neighbour and on the inner up side 890 IF( ij == jpj-1 .and. (nbondj == 0 .or. nbondj == -1) ) lsend_bdy(ib_bdy,igrd,4 ) = .true.892 IF( ij == jpj-1 .and. (nbondj == 0 .or. nbondj == -1) ) lsend_bdy(ib_bdy,igrd,4,ir) = .true. 891 893 ! 892 894 ! check if point has to be received from a neighbour 893 895 ! W neighbour and on the outter left side 894 IF( ii == 1 .and. (nbondi == 0 .or. nbondi == 1) ) lrecv_bdy(ib_bdy,igrd,1) = .true.896 IF( ii == 1 .and. (nbondi == 0 .or. nbondi == 1) ) lrecv_bdy(ib_bdy,igrd,1,ir) = .true. 895 897 ! E neighbour and on the outter right side 896 IF( ii == jpi .and. (nbondi == 0 .or. nbondi == -1) ) lrecv_bdy(ib_bdy,igrd,2) = .true.898 IF( ii == jpi .and. (nbondi == 0 .or. nbondi == -1) ) lrecv_bdy(ib_bdy,igrd,2,ir) = .true. 897 899 ! S neighbour and on the outter down side 898 IF( ij == 1 .and. (nbondj == 0 .or. nbondj == 1) ) lrecv_bdy(ib_bdy,igrd,3) = .true.900 IF( ij == 1 .and. (nbondj == 0 .or. nbondj == 1) ) lrecv_bdy(ib_bdy,igrd,3,ir) = .true. 899 901 ! N neighbour and on the outter up side 900 IF( ij == jpj .and. (nbondj == 0 .or. nbondj == -1) ) lrecv_bdy(ib_bdy,igrd,4) = .true.902 IF( ij == jpj .and. (nbondj == 0 .or. nbondj == -1) ) lrecv_bdy(ib_bdy,igrd,4,ir) = .true. 901 903 ! 902 904 END DO … … 907 909 DO igrd = 1, jpbgrd 908 910 DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd) 909 nbr => idx_bdy(ib_bdy)%nbr(ib,igrd)910 idx_bdy(ib_bdy)%nbw(ib,igrd) = 1.- TANH( REAL( nbr - 1 ) *0.5 ) ! tanh formulation911 ! idx_bdy(ib_bdy)%nbw(ib,igrd) = (REAL(nn_rimwidth(ib_bdy)+1- nbr)/REAL(nn_rimwidth(ib_bdy)))**2. ! quadratic912 ! idx_bdy(ib_bdy)%nbw(ib,igrd) = REAL(nn_rimwidth(ib_bdy)+1- nbr)/REAL(nn_rimwidth(ib_bdy)) ! linear911 ir = MAX( 1, idx_bdy(ib_bdy)%nbr(ib,igrd) ) ! both rim 0 and rim 1 have the same weights 912 idx_bdy(ib_bdy)%nbw(ib,igrd) = 1.- TANH( REAL( ir - 1 ) *0.5 ) ! tanh formulation 913 ! idx_bdy(ib_bdy)%nbw(ib,igrd) = (REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy)))**2. ! quadratic 914 ! idx_bdy(ib_bdy)%nbw(ib,igrd) = REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy)) ! linear 913 915 END DO 914 916 END DO … … 918 920 DO igrd = 1, jpbgrd 919 921 DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd) 920 nbr => idx_bdy(ib_bdy)%nbr(ib,igrd)922 ir = MAX( 1, idx_bdy(ib_bdy)%nbr(ib,igrd) ) ! both rim 0 and rim 1 have the same damping coefficients 921 923 idx_bdy(ib_bdy)%nbd(ib,igrd) = 1. / ( rn_time_dmp(ib_bdy) * rday ) & 922 & *(REAL(nn_rimwidth(ib_bdy)+1- nbr)/REAL(nn_rimwidth(ib_bdy)))**2. ! quadratic924 & *(REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy)))**2. ! quadratic 923 925 idx_bdy(ib_bdy)%nbdout(ib,igrd) = 1. / ( rn_time_dmp_out(ib_bdy) * rday ) & 924 & *(REAL(nn_rimwidth(ib_bdy)+1- nbr)/REAL(nn_rimwidth(ib_bdy)))**2. ! quadratic926 & *(REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy)))**2. ! quadratic 925 927 END DO 926 928 END DO … … 931 933 ! Initialise masks and find normal/tangential directions 932 934 ! ------------------------------------------------------ 935 936 ! ------------------------------------------ 937 ! handel rim0, do as if rim 1 was free ocean 938 ! ------------------------------------------ 939 940 ztmask(:,:) = tmask(:,:,1) ; zumask(:,:) = umask(:,:,1) ; zvmask(:,:) = vmask(:,:,1) 941 ! For the flagu/flagv calculation below we require a version of fmask without 942 ! the land boundary condition (shlat) included: 943 DO ij = 2, jpjm1 944 DO ii = 2, jpim1 945 zfmask(ii,ij) = ztmask(ii,ij ) * ztmask(ii+1,ij ) & 946 & * ztmask(ii,ij+1) * ztmask(ii+1,ij+1) 947 END DO 948 END DO 949 CALL lbc_lnk( 'bdyini', zfmask, 'F', 1. ) 933 950 934 951 ! Read global 2D mask at T-points: bdytmask … … 940 957 941 958 ! Derive mask on U and V grid from mask on T grid 942 943 bdyumask(:,:) = 0._wp944 bdyvmask(:,:) = 0._wp945 959 DO ij = 1, jpjm1 946 960 DO ii = 1, jpim1 947 bdyumask(ii,ij) = bdytmask(ii,ij) * bdytmask(ii+1, ij)961 bdyumask(ii,ij) = bdytmask(ii,ij) * bdytmask(ii+1,ij ) 948 962 bdyvmask(ii,ij) = bdytmask(ii,ij) * bdytmask(ii ,ij+1) 949 963 END DO 950 964 END DO 951 CALL lbc_lnk_multi( 'bdyini', bdyumask, 'U', 1. , bdyvmask, 'V', 1. ) ! Lateral boundary cond. 952 953 ! bdy masks are now set to zero on boundary points: 954 ! 955 igrd = 1 965 CALL lbc_lnk_multi( 'bdyini', bdyumask, 'U', 1., bdyvmask, 'V', 1. ) ! Lateral boundary cond. 966 967 ! bdy masks are now set to zero on rim 0 points: 956 968 DO ib_bdy = 1, nb_bdy 957 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 958 bdytmask(idx_bdy(ib_bdy)%nbi(ib,igrd), idx_bdy(ib_bdy)%nbj(ib,igrd)) = 0._wp 959 END DO 969 DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(1) ! extent of rim 0 970 bdytmask(idx_bdy(ib_bdy)%nbi(ib,1), idx_bdy(ib_bdy)%nbj(ib,1)) = 0._wp 971 END DO 972 DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(2) ! extent of rim 0 973 bdyumask(idx_bdy(ib_bdy)%nbi(ib,2), idx_bdy(ib_bdy)%nbj(ib,2)) = 0._wp 974 END DO 975 DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(3) ! extent of rim 0 976 bdyvmask(idx_bdy(ib_bdy)%nbi(ib,3), idx_bdy(ib_bdy)%nbj(ib,3)) = 0._wp 977 END DO 960 978 END DO 961 ! 962 igrd = 2 979 980 CALL bdy_rim_treat( zumask, zvmask, zfmask, .true. ) ! compute flagu, flagv, ntreat on rim 0 981 982 ! ------------------------------------ 983 ! handel rim1, do as if rim 0 was land 984 ! ------------------------------------ 985 963 986 DO ib_bdy = 1, nb_bdy 964 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 965 bdyumask(idx_bdy(ib_bdy)%nbi(ib,igrd), idx_bdy(ib_bdy)%nbj(ib,igrd)) = 0._wp 966 END DO 987 DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(1) ! extent of rim 0 988 ztmask(idx_bdy(ib_bdy)%nbi(ib,1), idx_bdy(ib_bdy)%nbj(ib,1)) = 0._wp 989 END DO 990 DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(2) ! extent of rim 0 991 zumask(idx_bdy(ib_bdy)%nbi(ib,2), idx_bdy(ib_bdy)%nbj(ib,2)) = 0._wp 992 END DO 993 DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(3) ! extent of rim 0 994 zvmask(idx_bdy(ib_bdy)%nbi(ib,3), idx_bdy(ib_bdy)%nbj(ib,3)) = 0._wp 995 END DO 967 996 END DO 968 ! 969 igrd = 3 970 DO ib_bdy = 1, nb_bdy 971 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 972 bdyvmask(idx_bdy(ib_bdy)%nbi(ib,igrd), idx_bdy(ib_bdy)%nbj(ib,igrd)) = 0._wp 973 END DO 974 END DO 975 976 ! For the flagu/flagv calculation below we require a version of fmask without 977 ! the land boundary condition (shlat) included: 978 zfmask(:,:) = 0 997 998 ! Recompute zfmask 979 999 DO ij = 2, jpjm1 980 1000 DO ii = 2, jpim1 981 zfmask(ii,ij) = tmask(ii,ij ,1) * tmask(ii+1,ij ,1) &982 & * tmask(ii,ij+1,1) * tmask(ii+1,ij+1,1)1001 zfmask(ii,ij) = ztmask(ii,ij ) * ztmask(ii+1,ij ) & 1002 & * ztmask(ii,ij+1) * ztmask(ii+1,ij+1) 983 1003 END DO 984 1004 END DO 985 986 ! Lateral boundary conditions 987 CALL lbc_lnk( 'bdyini', zfmask, 'F', 1. ) 988 CALL lbc_lnk_multi( 'bdyini', bdyumask, 'U', 1. , bdyvmask, 'V', 1., bdytmask, 'T', 1. ) 1005 CALL lbc_lnk( 'bdyini', zfmask, 'F', 1. ) 1006 1007 ! bdy masks are now set to zero on rim1 points: 1008 DO ib_bdy = 1, nb_bdy 1009 DO ib = idx_bdy(ib_bdy)%nblenrim0(1) + 1, idx_bdy(ib_bdy)%nblenrim(1) ! extent of rim 1 1010 bdytmask(idx_bdy(ib_bdy)%nbi(ib,1), idx_bdy(ib_bdy)%nbj(ib,1)) = 0._wp 1011 END DO 1012 DO ib = idx_bdy(ib_bdy)%nblenrim0(2) + 1, idx_bdy(ib_bdy)%nblenrim(2) ! extent of rim 1 1013 bdyumask(idx_bdy(ib_bdy)%nbi(ib,2), idx_bdy(ib_bdy)%nbj(ib,2)) = 0._wp 1014 END DO 1015 DO ib = idx_bdy(ib_bdy)%nblenrim0(3) + 1, idx_bdy(ib_bdy)%nblenrim(3) ! extent of rim 1 1016 bdyvmask(idx_bdy(ib_bdy)%nbi(ib,3), idx_bdy(ib_bdy)%nbj(ib,3)) = 0._wp 1017 END DO 1018 END DO 1019 1020 CALL bdy_rim_treat( zumask, zvmask, zfmask, .false. ) ! compute flagu, flagv, ntreat on rim 1 1021 1022 1023 ! 1024 ! Check which boundaries might need communication 1025 ALLOCATE( lsend_bdyint(nb_bdy,jpbgrd,4,0:1), lrecv_bdyint(nb_bdy,jpbgrd,4,0:1) ) 1026 lsend_bdyint(:,:,:,:) = .false. 1027 lrecv_bdyint(:,:,:,:) = .false. 1028 ALLOCATE( lsend_bdyext(nb_bdy,jpbgrd,4,0:1), lrecv_bdyext(nb_bdy,jpbgrd,4,0:1) ) 1029 lsend_bdyext(:,:,:,:) = .false. 1030 lrecv_bdyext(:,:,:,:) = .false. 1031 ! 1032 DO igrd = 1, jpbgrd 1033 DO ib_bdy = 1, nb_bdy 1034 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1035 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 1036 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 1037 IF( ib .LE. idx_bdy(ib_bdy)%nblenrim0(igrd) ) THEN ; ir = 0 1038 ELSE ; ir = 1 1039 END IF 1040 flagu = NINT(idx_bdy(ib_bdy)%flagu(ib,igrd)) 1041 flagv = NINT(idx_bdy(ib_bdy)%flagv(ib,igrd)) 1042 iibe = ii - flagu ! neighbouring point towards the exterior of the computational domain 1043 ijbe = ij - flagv 1044 iibi = ii + flagu ! neighbouring point towards the interior of the computational domain 1045 ijbi = ij + flagv 1046 CALL find_neib( ii, ij, idx_bdy(ib_bdy)%ntreat(ib,igrd), ii1, ij1, ii2, ij2, ii3, ij3 ) ! free ocean neighbours 1047 ! 1048 ! search neighbour in the west/east direction 1049 ! Rim is on the halo and computed ocean is towards exterior of mpi domain 1050 ! <-- (o exterior) --> 1051 ! (1) o|x OR (2) x|o 1052 ! |___ ___| 1053 IF( iibi == 0 .OR. ii1 == 0 .OR. ii2 == 0 .OR. ii3 == 0 ) lrecv_bdyint(ib_bdy,igrd,1,ir) = .true. 1054 IF( iibi == jpi+1 .OR. ii1 == jpi+1 .OR. ii2 == jpi+1 .OR. ii3 == jpi+1 ) lrecv_bdyint(ib_bdy,igrd,2,ir) = .true. 1055 IF( iibe == 0 ) lrecv_bdyext(ib_bdy,igrd,1,ir) = .true. 1056 IF( iibe == jpi+1 ) lrecv_bdyext(ib_bdy,igrd,2,ir) = .true. 1057 ! Check if neighbour has its rim parallel to its mpi subdomain border and located next to its halo 1058 ! :¨¨¨¨¨|¨¨--> | | <--¨¨|¨¨¨¨¨: 1059 ! : | x:o | neighbour limited by ... would need o | o:x | : 1060 ! :.....|_._:_____| (1) W neighbour E neighbour (2) |_____:_._|.....: 1061 IF( ii == 2 .AND. ( nbondi == 1 .OR. nbondi == 0 ) .AND. & 1062 & ( iibi == 3 .OR. ii1 == 3 .OR. ii2 == 3 .OR. ii3 == 3 ) ) lsend_bdyint(ib_bdy,igrd,1,ir) = .true. 1063 IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. & 1064 & ( iibi == jpi-2 .OR. ii1 == jpi-2 .OR. ii2 == jpi-2 .OR. ii3 == jpi-2) ) lsend_bdyint(ib_bdy,igrd,2,ir) = .true. 1065 IF( ii == 2 .AND. ( nbondi == 1 .OR. nbondi == 0 ) .AND. iibe == 3 ) lsend_bdyext(ib_bdy,igrd,1,ir) = .true. 1066 IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. iibe == jpi-2) lsend_bdyext(ib_bdy,igrd,2,ir) = .true. 1067 ! 1068 ! search neighbour in the north/south direction 1069 ! Rim is on the halo and computed ocean is towards exterior of mpi domain 1070 !(3) | | ^ ___o___ 1071 ! | |___x___| OR | | x | 1072 ! v o (4) | | 1073 IF( ijbi == 0 .OR. ij1 == 0 .OR. ij2 == 0 .OR. ij3 == 0 ) lrecv_bdyint(ib_bdy,igrd,3,ir) = .true. 1074 IF( ijbi == jpj+1 .OR. ij1 == jpj+1 .OR. ij2 == jpj+1 .OR. ij3 == jpj+1 ) lrecv_bdyint(ib_bdy,igrd,4,ir) = .true. 1075 IF( ijbe == 0 ) lrecv_bdyext(ib_bdy,igrd,3,ir) = .true. 1076 IF( ijbe == jpj+1 ) lrecv_bdyext(ib_bdy,igrd,4,ir) = .true. 1077 ! Check if neighbour has its rim parallel to its mpi subdomain _________ border and next to its halo 1078 ! ^ | o | : : 1079 ! | |¨¨¨¨x¨¨¨¨| neighbour limited by ... would need o | |....x....| 1080 ! :_________: (3) S neighbour N neighbour (4) v | o | 1081 IF( ij == 2 .AND. ( nbondj == 1 .OR. nbondj == 0 ) .AND. & 1082 & ( ijbi == 3 .OR. ij1 == 3 .OR. ij2 == 3 .OR. ij3 == 3 ) ) lsend_bdyint(ib_bdy,igrd,3,ir) = .true. 1083 IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. & 1084 & ( ijbi == jpj-2 .OR. ij1 == jpj-2 .OR. ij2 == jpj-2 .OR. ij3 == jpj-2) ) lsend_bdyint(ib_bdy,igrd,4,ir) = .true. 1085 IF( ij == 2 .AND. ( nbondj == 1 .OR. nbondj == 0 ) .AND. ijbe == 3 ) lsend_bdyext(ib_bdy,igrd,3,ir) = .true. 1086 IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. ijbe == jpj-2) lsend_bdyext(ib_bdy,igrd,4,ir) = .true. 1087 END DO 1088 END DO 1089 END DO 1090 1091 DO ib_bdy = 1,nb_bdy 1092 IF( cn_dyn2d(ib_bdy) == 'orlanski' .OR. cn_dyn2d(ib_bdy) == 'orlanski_npo' .OR. & 1093 & cn_dyn3d(ib_bdy) == 'orlanski' .OR. cn_dyn3d(ib_bdy) == 'orlanski_npo' .OR. & 1094 & cn_tra(ib_bdy) == 'orlanski' .OR. cn_tra(ib_bdy) == 'orlanski_npo' ) THEN 1095 DO igrd = 1, jpbgrd 1096 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1097 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 1098 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 1099 IF( mig(ii) > 2 .AND. mig(ii) < jpiglo-2 .AND. mjg(ij) > 2 .AND. mjg(ij) < jpjglo-2 ) THEN 1100 WRITE(ctmp1,*) ' Orlanski is not safe when the open boundaries are on the interior of the computational domain' 1101 WRITE(ctmp2,*) ' ========== ' 1102 CALL ctl_stop( ' ', ctmp1, ctmp2, ' ' ) 1103 END IF 1104 END DO 1105 END DO 1106 END IF 1107 END DO 1108 1109 ! 1110 ! Tidy up 1111 !-------- 1112 IF( nb_bdy>0 ) DEALLOCATE( nbidta, nbjdta, nbrdta ) 1113 ! 1114 END SUBROUTINE bdy_segs 1115 1116 1117 SUBROUTINE bdy_rim_treat( zumask, zvmask, zfmask, lrim0 ) 1118 !!---------------------------------------------------------------------- 1119 !! *** ROUTINE bdy_rim_treat *** 1120 !! 1121 !! ** Purpose : Initialize structures ( flagu, flagv, ntreat ) indicating how rim points 1122 !! are to be handeled in the boundary condition treatment 1123 !! 1124 !! ** Method : - to handel rim 0 zmasks must indicate ocean points (set at one on rim 0 and rim 1 and interior) 1125 !! and bdymasks must be set at 0 on rim 0 (set at one on rim 1 and interior) 1126 !! (as if rim 1 was free ocean) 1127 !! - to handel rim 1 zmasks must be set at 0 on rim 0 (set at one on rim 1 and interior) 1128 !! and bdymasks must indicate free ocean points (set at one on interior) 1129 !! (as if rim 0 was land) 1130 !! - we can then check in which direction the interior of the computational domain is with the difference 1131 !! mask array values on both sides to compute flagu and flagv 1132 !! - and look at the ocean neighbours to compute ntreat 1133 !!---------------------------------------------------------------------- 1134 REAL(wp), TARGET, DIMENSION(jpi,jpj), INTENT (in ) :: zfmask ! temporary fmask excluding coastal boundary condition (shlat) 1135 REAL(wp), TARGET, DIMENSION(jpi,jpj), INTENT (in ) :: zumask, zvmask ! temporary t/u/v mask array 1136 LOGICAL , INTENT (in ) :: lrim0 ! .true. -> rim 0 .false. -> rim 1 1137 INTEGER :: ib_bdy, ii, ij, igrd, ib, icount ! dummy loop indices 1138 INTEGER :: i_offset, j_offset, inn ! local integer 1139 INTEGER :: ibeg, iend ! local integer 1140 LOGICAL :: llnon, llson, llean, llwen ! local logicals indicating the presence of a ocean neighbour 1141 REAL(wp), POINTER, DIMENSION(:,:) :: pmask ! pointer to 2D mask fields 1142 REAL(wp) :: zefl, zwfl, znfl, zsfl ! local scalars 1143 CHARACTER(LEN=1), DIMENSION(jpbgrd) :: cgrid 1144 REAL(wp) , DIMENSION(jpi,jpj) :: ztmp 1145 !!---------------------------------------------------------------------- 1146 1147 cgrid = (/'t','u','v'/) 1148 989 1149 DO ib_bdy = 1, nb_bdy ! Indices and directions of rim velocity components 990 icount = 0991 1150 992 1151 ! Calculate relationship of U direction to the local orientation of the boundary … … 994 1153 ! flagu = 0 : u is tangential 995 1154 ! flagu = 1 : u is normal to the boundary and is direction is inward 996 997 1155 DO igrd = 1, jpbgrd 998 1156 SELECT CASE( igrd ) 999 CASE( 1 ) ; pmask => umask (:,:,1); i_offset = 01000 CASE( 2 ) ; pmask => bdytmask (:,:); i_offset = 11001 CASE( 3 ) ; pmask => zfmask (:,:); i_offset = 01157 CASE( 1 ) ; pmask => zumask ; i_offset = 0 1158 CASE( 2 ) ; pmask => bdytmask ; i_offset = 1 1159 CASE( 3 ) ; pmask => zfmask ; i_offset = 0 1002 1160 END SELECT 1003 1161 icount = 0 1004 1162 ztmp(:,:) = 0._wp 1005 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1163 IF( lrim0 ) THEN ! extent of rim 0 1164 ibeg = 1 ; iend = idx_bdy(ib_bdy)%nblenrim0(igrd) 1165 ELSE ! extent of rim 1 1166 ibeg = idx_bdy(ib_bdy)%nblenrim0(igrd) + 1 ; iend = idx_bdy(ib_bdy)%nblenrim(igrd) 1167 END IF 1168 DO ib = ibeg, iend 1006 1169 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 1007 1170 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 1008 1171 IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj ) CYCLE 1009 z efl = pmask(ii+i_offset-1,ij)1010 z wfl = pmask(ii+i_offset ,ij)1172 zwfl = pmask(ii+i_offset-1,ij) 1173 zefl = pmask(ii+i_offset ,ij) 1011 1174 ! This error check only works if you are using the bdyXmask arrays 1012 1175 IF( i_offset == 1 .and. zefl + zwfl == 2. ) THEN … … 1014 1177 IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(ii),mjg(ij) 1015 1178 ELSE 1016 ztmp(ii,ij) = -z efl + zwfl1179 ztmp(ii,ij) = -zwfl + zefl 1017 1180 ENDIF 1018 1181 END DO … … 1024 1187 ENDIF 1025 1188 CALL lbc_lnk( 'bdyini', ztmp, 'T', 1. ) 1026 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)1189 DO ib = ibeg, iend 1027 1190 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 1028 1191 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) … … 1035 1198 ! flagv = 0 : v is tangential 1036 1199 ! flagv = 1 : v is normal to the boundary and is direction is inward 1037 1038 1200 DO igrd = 1, jpbgrd 1039 1201 SELECT CASE( igrd ) 1040 CASE( 1 ) ; pmask => vmask (:,:,1); j_offset = 01041 CASE( 2 ) ; pmask => zfmask (:,:); j_offset = 01042 CASE( 3 ) ; pmask => bdytmask 1202 CASE( 1 ) ; pmask => zvmask ; j_offset = 0 1203 CASE( 2 ) ; pmask => zfmask ; j_offset = 0 1204 CASE( 3 ) ; pmask => bdytmask ; j_offset = 1 1043 1205 END SELECT 1044 1206 icount = 0 1045 1207 ztmp(:,:) = 0._wp 1046 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1208 IF( lrim0 ) THEN ! extent of rim 0 1209 ibeg = 1 ; iend = idx_bdy(ib_bdy)%nblenrim0(igrd) 1210 ELSE ! extent of rim 1 1211 ibeg = idx_bdy(ib_bdy)%nblenrim0(igrd) + 1 ; iend = idx_bdy(ib_bdy)%nblenrim(igrd) 1212 END IF 1213 DO ib = ibeg, iend 1047 1214 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 1048 1215 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 1049 1216 IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj ) CYCLE 1050 z nfl = pmask(ii,ij+j_offset-1)1051 z sfl = pmask(ii,ij+j_offset )1217 zsfl = pmask(ii,ij+j_offset-1) 1218 znfl = pmask(ii,ij+j_offset ) 1052 1219 ! This error check only works if you are using the bdyXmask arrays 1053 1220 IF( j_offset == 1 .and. znfl + zsfl == 2. ) THEN … … 1055 1222 icount = icount + 1 1056 1223 ELSE 1057 ztmp(ii,ij) = -z nfl + zsfl1224 ztmp(ii,ij) = -zsfl + znfl 1058 1225 END IF 1059 1226 END DO … … 1065 1232 ENDIF 1066 1233 CALL lbc_lnk( 'bdyini', ztmp, 'T', 1. ) 1067 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)1234 DO ib = ibeg, iend 1068 1235 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 1069 1236 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) … … 1072 1239 END DO 1073 1240 ! 1074 END DO 1241 END DO ! ib_bdy 1075 1242 1076 1243 DO ib_bdy = 1, nb_bdy … … 1082 1249 END SELECT 1083 1250 ztmp(:,:) = 0._wp 1084 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1251 IF( lrim0 ) THEN ! extent of rim 0 1252 ibeg = 1 ; iend = idx_bdy(ib_bdy)%nblenrim0(igrd) 1253 ELSE ! extent of rim 1 1254 ibeg = idx_bdy(ib_bdy)%nblenrim0(igrd) + 1 ; iend = idx_bdy(ib_bdy)%nblenrim(igrd) 1255 END IF 1256 DO ib = ibeg, iend 1085 1257 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 1086 1258 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 1087 IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj ) CYCLE 1088 ll ean = pmask(ii+1,ij ) == 1.1089 ll wen = pmask(ii-1,ij ) == 1.1090 ll non = pmask(ii ,ij+1) == 1.1091 ll son = pmask(ii ,ij-1) == 1.1259 IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj ) CYCLE 1260 llnon = pmask(ii ,ij+1) == 1. 1261 llson = pmask(ii ,ij-1) == 1. 1262 llean = pmask(ii+1,ij ) == 1. 1263 llwen = pmask(ii-1,ij ) == 1. 1092 1264 inn = COUNT( (/ llnon, llson, llean, llwen /) ) 1093 1265 IF( inn == 0 ) THEN ! no neighbours -> interior of a corner or cluster of rim points … … 1101 1273 ELSE 1102 1274 WRITE(ctmp1,*) ' E R R O R : Problem with ',cgrid(igrd) ,' grid point', ii, ij, & 1103 ' on boundary set ', ib_bdy, ' has no free ocean neighbour' 1104 WRITE(ctmp2,*) ' There seems to be a cluster of rim points.' 1275 ' on boundary set ', ib_bdy, ' has no free ocean neighbour' 1276 IF( lrim0 ) THEN 1277 WRITE(ctmp2,*) ' There seems to be a cluster of rim 0 points.' 1278 ELSE 1279 WRITE(ctmp2,*) ' There seems to be a cluster of rim 1 points.' 1280 END IF 1105 1281 WRITE(ctmp3,*) ' ========== ' 1106 1282 CALL ctl_stop( ' ', ctmp1, ctmp2, ctmp3, ' ' ) … … 1142 1318 END DO 1143 1319 CALL lbc_lnk( 'bdyini', ztmp, 'T', 1. ) 1144 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)1320 DO ib = ibeg, iend 1145 1321 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 1146 1322 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) … … 1150 1326 END DO 1151 1327 1152 1153 ALLOCATE( lsend_bdyint(nb_bdy,jpbgrd,4), lrecv_bdyint(nb_bdy,jpbgrd,4) ) 1154 lsend_bdyint(:,:,:) = .false. 1155 lrecv_bdyint(:,:,:) = .false. 1156 ALLOCATE( lsend_bdyext(nb_bdy,jpbgrd,4), lrecv_bdyext(nb_bdy,jpbgrd,4) ) 1157 lsend_bdyext(:,:,:) = .false. 1158 lrecv_bdyext(:,:,:) = .false. 1159 ! 1160 ! Check which boundaries might need communication 1161 DO igrd = 1, jpbgrd 1162 DO ib_bdy = 1, nb_bdy 1163 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1164 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 1165 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 1166 flagu = NINT(idx_bdy(ib_bdy)%flagu(ib,igrd)) 1167 flagv = NINT(idx_bdy(ib_bdy)%flagv(ib,igrd)) 1168 iibe = ii - flagu ! neighbouring point towards the exterior of the computational domain 1169 ijbe = ij - flagv 1170 iibi = ii + flagu ! neighbouring point towards the interior of the computational domain 1171 ijbi = ij + flagv 1172 CALL find_neib( ii, ij, idx_bdy(ib_bdy)%ntreat(ib,igrd), ii1, ij1, ii2, ij2, ii3, ij3 ) 1173 ! 1174 ! search neighbour in the west/east direction 1175 ! Rim is on the halo and computed ocean is towards exterior of mpi domain 1176 ! <-- (o exterior) --> 1177 ! (1) o|x OR (2) x|o 1178 ! |___ ___| 1179 IF( iibi == 0 .OR. ii1 == 0 .OR. ii2 == 0 .OR. ii3 == 0 ) lrecv_bdyint(ib_bdy,igrd,1) = .true. 1180 IF( iibi == jpi+1 .OR. ii1 == jpi+1 .OR. ii2 == jpi+1 .OR. ii3 == jpi+1 ) lrecv_bdyint(ib_bdy,igrd,2) = .true. 1181 IF( iibe == 0 ) lrecv_bdyext(ib_bdy,igrd,1) = .true. 1182 IF( iibe == jpi+1 ) lrecv_bdyext(ib_bdy,igrd,2) = .true. 1183 ! Check if neighbour has its rim parallel to its mpi subdomain border and located next to its halo 1184 ! :¨¨¨¨¨|¨¨--> | | <--¨¨|¨¨¨¨¨: 1185 ! : | x:o | neighbour limited by ... would need o | o:x | : 1186 ! :.....|_._:_____| (1) W neighbour E neighbour (2) |_____:_._|.....: 1187 IF( ii == 2 .AND. ( nbondi == 1 .OR. nbondi == 0 ) .AND. & 1188 & ( iibi == 3 .OR. ii1 == 3 .OR. ii2 == 3 .OR. ii3 == 3 ) ) lsend_bdyint(ib_bdy,igrd,1) = .true. 1189 IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. & 1190 & ( iibi == jpi-2 .OR. ii1 == jpi-2 .OR. ii2 == jpi-2 .OR. ii3 == jpi-2) ) lsend_bdyint(ib_bdy,igrd,2) = .true. 1191 IF( ii == 2 .AND. ( nbondi == 1 .OR. nbondi == 0 ) .AND. iibe == 3 ) lsend_bdyext(ib_bdy,igrd,1) = .true. 1192 IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. iibe == jpi-2 ) lsend_bdyext(ib_bdy,igrd,2) = .true. 1193 ! 1194 ! search neighbour in the north/south direction 1195 ! Rim is on the halo and computed ocean is towards exterior of mpi domain 1196 !(3) | | ^ ___o___ 1197 ! | |___x___| OR | | x | 1198 ! v o (4) | | 1199 IF( ijbi == 0 .OR. ij1 == 0 .OR. ij2 == 0 .OR. ij3 == 0 ) lrecv_bdyint(ib_bdy,igrd,3) = .true. 1200 IF( ijbi == jpj+1 .OR. ij1 == jpj+1 .OR. ij2 == jpj+1 .OR. ij3 == jpj+1 ) lrecv_bdyint(ib_bdy,igrd,4) = .true. 1201 IF( ijbe == 0 ) lrecv_bdyext(ib_bdy,igrd,3) = .true. 1202 IF( ijbe == jpj+1 ) lrecv_bdyext(ib_bdy,igrd,4) = .true. 1203 ! Check if neighbour has its rim parallel to its mpi subdomain _________ border and next to its halo 1204 ! ^ | o | : : 1205 ! | |¨¨¨¨x¨¨¨¨| neighbour limited by ... would need o | |....x....| 1206 ! :_________: (3) S neighbour N neighbour (4) v | o | 1207 IF( ij == 2 .AND. ( nbondj == 1 .OR. nbondj == 0 ) .AND. & 1208 & ( ijbi == 3 .OR. ij1 == 3 .OR. ij2 == 3 .OR. ij3 == 3 ) ) lsend_bdyint(ib_bdy,igrd,3) = .true. 1209 IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. & 1210 & ( ijbi == jpj-2 .OR. ij1 == jpj-2 .OR. ij2 == jpj-2 .OR. ij3 == jpj-2) ) lsend_bdyint(ib_bdy,igrd,4) = .true. 1211 IF( ij == 2 .AND. ( nbondj == 1 .OR. nbondj == 0 ) .AND. ijbe == 3 ) lsend_bdyext(ib_bdy,igrd,3) = .true. 1212 IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. ijbe == jpj-2 ) lsend_bdyext(ib_bdy,igrd,4) = .true. 1213 END DO 1214 END DO 1215 END DO 1216 1217 DO ib_bdy = 1,nb_bdy 1218 IF( cn_dyn2d(ib_bdy) == 'orlanski' .OR. cn_dyn2d(ib_bdy) == 'orlanski_npo' .OR. & 1219 & cn_dyn3d(ib_bdy) == 'orlanski' .OR. cn_dyn3d(ib_bdy) == 'orlanski_npo' .OR. & 1220 & cn_tra(ib_bdy) == 'orlanski' .OR. cn_tra(ib_bdy) == 'orlanski_npo' ) THEN 1221 DO igrd = 1, jpbgrd 1222 DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 1223 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 1224 ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 1225 IF( mig(ii) > 2 .AND. mig(ii) < jpiglo-2 .AND. mjg(ij) > 2 .AND. mjg(ij) < jpjglo-2 ) THEN 1226 WRITE(ctmp1,*) ' Orlanski can not be used when the open boundaries are on the interior of the computational domain' 1227 WRITE(ctmp2,*) ' ========== ' 1228 CALL ctl_stop( ' ', ctmp1, ctmp2, ' ' ) 1229 END IF 1230 END DO 1231 END DO 1232 END IF 1233 END DO 1234 1235 ! 1236 ! Tidy up 1237 !-------- 1238 IF( nb_bdy>0 ) DEALLOCATE( nbidta, nbjdta, nbrdta ) 1239 ! 1240 END SUBROUTINE bdy_segs 1241 1242 1243 SUBROUTINE find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 ) 1244 !!---------------------------------------------------------------------- 1245 !! *** ROUTINE find_neib *** 1246 !! 1247 !! ** Purpose : get ii1, ij1, ii2, ij2, ii3, ij3, the indices of 1248 !! the free ocean neighbours of (ii,ij) for bdy treatment 1249 !! 1250 !!---------------------------------------------------------------------- 1251 INTEGER, INTENT(in ) :: ii, ij, itreat 1252 INTEGER, INTENT( out) :: ii1, ij1, ii2, ij2, ii3, ij3 1253 !!---------------------------------------------------------------------- 1254 SELECT CASE( itreat ) ! points that will be used by bdy routines, -1 will be discarded 1255 ! ! ! _____ ! _____ 1256 ! 1 | o ! 2 o | ! 3 | x ! 4 x | 1257 ! |_x_ _ ! _ _x_| ! | o ! o | 1258 CASE( 1 ) ; ii1 = ii+1 ; ij1 = ij+1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1259 CASE( 2 ) ; ii1 = ii-1 ; ij1 = ij+1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1260 CASE( 3 ) ; ii1 = ii+1 ; ij1 = ij-1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1261 CASE( 4 ) ; ii1 = ii-1 ; ij1 = ij-1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1262 ! | ! | ! o ! ______ ! or incomplete corner 1263 ! 5 | x o ! 6 o x | ! 7 __x__ ! 8 x ! 7 ____ o 1264 ! | ! | ! ! o ! |x___ 1265 CASE( 5 ) ; ii1 = ii+1 ; ij1 = ij ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1266 CASE( 6 ) ; ii1 = ii-1 ; ij1 = ij ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1267 CASE( 7 ) ; ii1 = ii ; ij1 = ij+1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1268 CASE( 8 ) ; ii1 = ii ; ij1 = ij-1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1269 ! o ! o ! _____| ! |_____ 1270 ! 9 ____x o ! 10 o x___ ! 11 x o ! 12 o x 1271 ! | ! | ! o ! o 1272 CASE( 9 ) ; ii1 = ii ; ij1 = ij+1 ; ii2 = ii+1 ; ij2 = ij ; ii3 = -1 ; ij3 = -1 1273 CASE( 10 ) ; ii1 = ii ; ij1 = ij+1 ; ii2 = ii-1 ; ij2 = ij ; ii3 = -1 ; ij3 = -1 1274 CASE( 11 ) ; ii1 = ii ; ij1 = ij-1 ; ii2 = ii+1 ; ij2 = ij ; ii3 = -1 ; ij3 = -1 1275 CASE( 12 ) ; ii1 = ii ; ij1 = ij-1 ; ii2 = ii-1 ; ij2 = ij ; ii3 = -1 ; ij3 = -1 1276 ! |_ o ! o _| ! ¨¨|_|¨¨ ! o 1277 ! 13 _| x o ! 14 o x |_ ! 15 o x o ! 16 o x o 1278 ! | o ! o | ! o ! __|¨|__ 1279 CASE( 13 ) ; ii1 = ii ; ij1 = ij+1 ; ii2 = ii+1 ; ij2 = ij ; ii3 = ii ; ij3 = ij-1 1280 CASE( 14 ) ; ii1 = ii ; ij1 = ij+1 ; ii2 = ii-1 ; ij2 = ij ; ii3 = ii ; ij3 = ij-1 1281 CASE( 15 ) ; ii1 = ii-1 ; ij1 = ij ; ii2 = ii ; ij2 = ij-1 ; ii3 = ii+1 ; ij3 = ij 1282 CASE( 16 ) ; ii1 = ii-1 ; ij1 = ij ; ii2 = ii ; ij2 = ij+1 ; ii3 = ii+1 ; ij3 = ij 1283 END SELECT 1284 END SUBROUTINE find_neib 1285 1328 END SUBROUTINE bdy_rim_treat 1329 1330 1331 SUBROUTINE find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 ) 1332 !!---------------------------------------------------------------------- 1333 !! *** ROUTINE find_neib *** 1334 !! 1335 !! ** Purpose : get ii1, ij1, ii2, ij2, ii3, ij3, the indices of 1336 !! the free ocean neighbours of (ii,ij) for bdy treatment 1337 !! 1338 !! ** Method : use itreat input to select a case 1339 !! N.B. ntreat is defined for all bdy points in routine bdy_segs 1340 !! 1341 !!---------------------------------------------------------------------- 1342 INTEGER, INTENT(in ) :: ii, ij, itreat 1343 INTEGER, INTENT( out) :: ii1, ij1, ii2, ij2, ii3, ij3 1344 !!---------------------------------------------------------------------- 1345 SELECT CASE( itreat ) ! points that will be used by bdy routines, -1 will be discarded 1346 ! ! ! _____ ! _____ 1347 ! 1 | o ! 2 o | ! 3 | x ! 4 x | 1348 ! |_x_ _ ! _ _x_| ! | o ! o | 1349 CASE( 1 ) ; ii1 = ii+1 ; ij1 = ij+1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1350 CASE( 2 ) ; ii1 = ii-1 ; ij1 = ij+1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1351 CASE( 3 ) ; ii1 = ii+1 ; ij1 = ij-1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1352 CASE( 4 ) ; ii1 = ii-1 ; ij1 = ij-1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1353 ! | ! | ! o ! ______ ! or incomplete corner 1354 ! 5 | x o ! 6 o x | ! 7 __x__ ! 8 x ! 7 ____ o 1355 ! | ! | ! ! o ! |x___ 1356 CASE( 5 ) ; ii1 = ii+1 ; ij1 = ij ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1357 CASE( 6 ) ; ii1 = ii-1 ; ij1 = ij ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1358 CASE( 7 ) ; ii1 = ii ; ij1 = ij+1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1359 CASE( 8 ) ; ii1 = ii ; ij1 = ij-1 ; ii2 = -1 ; ij2 = -1 ; ii3 = -1 ; ij3 = -1 1360 ! o ! o ! _____| ! |_____ 1361 ! 9 ____x o ! 10 o x___ ! 11 x o ! 12 o x 1362 ! | ! | ! o ! o 1363 CASE( 9 ) ; ii1 = ii ; ij1 = ij+1 ; ii2 = ii+1 ; ij2 = ij ; ii3 = -1 ; ij3 = -1 1364 CASE( 10 ) ; ii1 = ii ; ij1 = ij+1 ; ii2 = ii-1 ; ij2 = ij ; ii3 = -1 ; ij3 = -1 1365 CASE( 11 ) ; ii1 = ii ; ij1 = ij-1 ; ii2 = ii+1 ; ij2 = ij ; ii3 = -1 ; ij3 = -1 1366 CASE( 12 ) ; ii1 = ii ; ij1 = ij-1 ; ii2 = ii-1 ; ij2 = ij ; ii3 = -1 ; ij3 = -1 1367 ! |_ o ! o _| ! ¨¨|_|¨¨ ! o 1368 ! 13 _| x o ! 14 o x |_ ! 15 o x o ! 16 o x o 1369 ! | o ! o | ! o ! __|¨|__ 1370 CASE( 13 ) ; ii1 = ii ; ij1 = ij+1 ; ii2 = ii+1 ; ij2 = ij ; ii3 = ii ; ij3 = ij-1 1371 CASE( 14 ) ; ii1 = ii ; ij1 = ij+1 ; ii2 = ii-1 ; ij2 = ij ; ii3 = ii ; ij3 = ij-1 1372 CASE( 15 ) ; ii1 = ii-1 ; ij1 = ij ; ii2 = ii ; ij2 = ij-1 ; ii3 = ii+1 ; ij3 = ij 1373 CASE( 16 ) ; ii1 = ii-1 ; ij1 = ij ; ii2 = ii ; ij2 = ij+1 ; ii3 = ii+1 ; ij3 = ij 1374 END SELECT 1375 END SUBROUTINE find_neib 1376 1286 1377 1287 1378 SUBROUTINE bdy_ctl_seg -
NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdylib.F90
r11183 r11191 92 92 93 93 94 SUBROUTINE bdy_orl( idx, ptb, pta, dta, l l_npo )94 SUBROUTINE bdy_orl( idx, ptb, pta, dta, lrim0, ll_npo ) 95 95 !!---------------------------------------------------------------------- 96 96 !! *** SUBROUTINE bdy_orl *** … … 104 104 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: ptb ! before tracer field 105 105 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pta ! tracer trend 106 LOGICAL , OPTIONAL, INTENT(in) :: lrim0 ! indicate if rim 0 is treated 106 107 LOGICAL, INTENT(in) :: ll_npo ! switch for NPO version 107 108 !! … … 111 112 igrd = 1 ! Everything is at T-points here 112 113 ! 113 CALL bdy_orlanski_3d( idx, igrd, ptb(:,:,:), pta(:,:,:), dta, l l_npo )114 CALL bdy_orlanski_3d( idx, igrd, ptb(:,:,:), pta(:,:,:), dta, lrim0, ll_npo ) 114 115 ! 115 116 END SUBROUTINE bdy_orl 116 117 117 118 118 SUBROUTINE bdy_orlanski_2d( idx, igrd, phib, phia, phi_ext, l l_npo )119 SUBROUTINE bdy_orlanski_2d( idx, igrd, phib, phia, phi_ext, lrim0, ll_npo ) 119 120 !!---------------------------------------------------------------------- 120 121 !! *** SUBROUTINE bdy_orlanski_2d *** … … 132 133 REAL(wp), DIMENSION(:,:), INTENT(inout) :: phia ! model after 2D field (to be updated) 133 134 REAL(wp), DIMENSION(:) , INTENT(in ) :: phi_ext ! external forcing data 135 LOGICAL, OPTIONAL, INTENT(in ) :: lrim0 ! indicate if rim 0 is treated 134 136 LOGICAL , INTENT(in ) :: ll_npo ! switch for NPO version 135 137 ! … … 140 142 INTEGER :: ii_offset, ij_offset ! offsets for mask indices 141 143 INTEGER :: flagu, flagv ! short cuts 144 INTEGER :: ibeg, iend ! length of rim to be treated (rim 0 or rim 1 or both) 142 145 REAL(wp) :: zmask_x, zmask_y1, zmask_y2 143 146 REAL(wp) :: zex1, zex2, zey, zey1, zey2 … … 185 188 END SELECT 186 189 ! 187 DO jb = 1, idx%nblenrim(igrd) 190 IF( PRESENT(lrim0) ) THEN 191 IF( lrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(igrd) ! rim 0 192 ELSE ; ibeg = idx%nblenrim0(igrd)+1 ; iend = idx%nblenrim(igrd) ! rim 1 193 END IF 194 ELSE ; ibeg = 1 ; iend = idx%nblenrim(igrd) ! both 195 END IF 196 ! 197 DO jb = ibeg, iend 188 198 ii = idx%nbi(jb,igrd) 189 199 ij = idx%nbj(jb,igrd) … … 272 282 273 283 274 SUBROUTINE bdy_orlanski_3d( idx, igrd, phib, phia, phi_ext, l l_npo )284 SUBROUTINE bdy_orlanski_3d( idx, igrd, phib, phia, phi_ext, lrim0, ll_npo ) 275 285 !!---------------------------------------------------------------------- 276 286 !! *** SUBROUTINE bdy_orlanski_3d *** … … 288 298 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: phia ! model after 3D field (to be updated) 289 299 REAL(wp), DIMENSION(:,:) , INTENT(in ) :: phi_ext ! external forcing data 300 LOGICAL, OPTIONAL, INTENT(in ) :: lrim0 ! indicate if rim 0 is treated 290 301 LOGICAL , INTENT(in ) :: ll_npo ! switch for NPO version 291 302 ! … … 296 307 INTEGER :: ii_offset, ij_offset ! offsets for mask indices 297 308 INTEGER :: flagu, flagv ! short cuts 309 INTEGER :: ibeg, iend ! length of rim to be treated (rim 0 or rim 1 or both) 298 310 REAL(wp) :: zmask_x, zmask_y1, zmask_y2 299 311 REAL(wp) :: zex1, zex2, zey, zey1, zey2 … … 340 352 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for igrd in bdy_orlanksi_2d' ) 341 353 END SELECT 342 354 ! 355 IF( PRESENT(lrim0) ) THEN 356 IF( lrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(igrd) ! rim 0 357 ELSE ; ibeg = idx%nblenrim0(igrd)+1 ; iend = idx%nblenrim(igrd) ! rim 1 358 END IF 359 ELSE ; ibeg = 1 ; iend = idx%nblenrim(igrd) ! both 360 END IF 361 ! 343 362 DO jk = 1, jpk 344 363 ! 345 DO jb = 1, idx%nblenrim(igrd)364 DO jb = ibeg, iend 346 365 ii = idx%nbi(jb,igrd) 347 366 ij = idx%nbj(jb,igrd) … … 430 449 END SUBROUTINE bdy_orlanski_3d 431 450 432 SUBROUTINE bdy_nmn( idx, igrd, phia )451 SUBROUTINE bdy_nmn( idx, igrd, phia, lrim0 ) 433 452 !!---------------------------------------------------------------------- 434 453 !! *** SUBROUTINE bdy_nmn *** … … 444 463 !! ! o 445 464 !!---------------------------------------------------------------------- 446 INTEGER, INTENT(in ):: igrd ! grid index465 INTEGER, INTENT(in ) :: igrd ! grid index 447 466 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: phia ! model after 3D field (to be updated), must be masked 448 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 467 TYPE(OBC_INDEX), INTENT(in ) :: idx ! OBC indices 468 LOGICAL, OPTIONAL, INTENT(in ) :: lrim0 ! indicate if rim 0 is treated 449 469 !! 450 470 REAL(wp) :: zweight … … 453 473 INTEGER :: ii, ij ! 2D addresses 454 474 INTEGER :: ipkm1 ! size of phia third dimension minus 1 475 INTEGER :: ibeg, iend ! length of rim to be treated (rim 0 or rim 1 or both) 455 476 INTEGER :: ii1, ii2, ii3, ij1, ij2, ij3, itreat 456 477 !!---------------------------------------------------------------------- … … 465 486 END SELECT 466 487 ! 467 DO ib = 1, idx%nblenrim(igrd) 488 IF( PRESENT(lrim0) ) THEN 489 IF( lrim0 ) THEN ; ibeg = 1 ; iend = idx%nblenrim0(igrd) ! rim 0 490 ELSE ; ibeg = idx%nblenrim0(igrd)+1 ; iend = idx%nblenrim(igrd) ! rim 1 491 END IF 492 ELSE ; ibeg = 1 ; iend = idx%nblenrim(igrd) ! both 493 END IF 494 ! 495 DO ib = ibeg, iend 468 496 ii = idx%nbi(ib,igrd) 469 497 ij = idx%nbj(ib,igrd) 470 498 itreat = idx%ntreat(ib,igrd) 471 CALL find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 ) 499 CALL find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 ) ! find free ocean neighbours 472 500 SELECT CASE( itreat ) 473 501 CASE( 1:8 ) -
NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdytra.F90
r11071 r11191 49 49 INTEGER, INTENT(in) :: kt ! Main time step counter 50 50 ! 51 INTEGER :: ib_bdy, jn, igrd ! Loop indeces 52 TYPE(ztrabdy), DIMENSION(jpts) :: zdta ! Temporary data structure 53 LOGICAL, DIMENSION(4) :: llsend1, llrecv1 ! indicate how communications are to be carried out 51 INTEGER :: ib_bdy, jn, igrd, ir ! Loop indeces 52 TYPE(ztrabdy), DIMENSION(jpts) :: zdta ! Temporary data structure 53 LOGICAL :: llrim0 ! indicate if rim 0 is treated 54 LOGICAL, DIMENSION(4) :: llsend1, llrecv1 ! indicate how communications are to be carried out 54 55 !!---------------------------------------------------------------------- 55 56 igrd = 1 56 57 57 DO ib_bdy=1, nb_bdy 58 DO ir = 1, 0, -1 ! treat rim 1 before rim 0 59 IF( ir == 0 ) THEN ; llrim0 = .TRUE. 60 ELSE ; llrim0 = .FALSE. 61 END IF 62 DO ib_bdy=1, nb_bdy 63 ! 64 zdta(1)%tra => dta_bdy(ib_bdy)%tem 65 zdta(2)%tra => dta_bdy(ib_bdy)%sal 66 ! 67 DO jn = 1, jpts 68 ! 69 SELECT CASE( TRIM(cn_tra(ib_bdy)) ) 70 CASE('none' ) ; CYCLE 71 CASE('frs' ) ! treat the whole boundary at once 72 IF( ir == 0 ) CALL bdy_frs ( idx_bdy(ib_bdy), tsa(:,:,:,jn), zdta(jn)%tra ) 73 CASE('specified' ) ! treat the whole rim at once 74 IF( ir == 0 ) CALL bdy_spe ( idx_bdy(ib_bdy), tsa(:,:,:,jn), zdta(jn)%tra ) 75 CASE('neumann' ) ; CALL bdy_nmn ( idx_bdy(ib_bdy), igrd , tsa(:,:,:,jn), llrim0 ) ! tsa masked 76 CASE('orlanski' ) ; CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), & 77 & zdta(jn)%tra, llrim0, ll_npo=.false. ) 78 CASE('orlanski_npo') ; CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), & 79 & zdta(jn)%tra, llrim0, ll_npo=.true. ) 80 CASE('runoff' ) ; CALL bdy_rnf ( idx_bdy(ib_bdy), tsa(:,:,:,jn), jn, llrim0 ) 81 CASE DEFAULT ; CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' ) 82 END SELECT 83 ! 84 END DO 85 END DO 58 86 ! 59 zdta(1)%tra => dta_bdy(ib_bdy)%tem 60 zdta(2)%tra => dta_bdy(ib_bdy)%sal 61 ! 62 DO jn = 1, jpts 63 ! 87 llsend1(:) = .false. 88 llrecv1(:) = .false. 89 DO ib_bdy=1, nb_bdy 64 90 SELECT CASE( TRIM(cn_tra(ib_bdy)) ) 65 CASE('none' ) ; CYCLE 66 CASE('frs' ) ; CALL bdy_frs ( idx_bdy(ib_bdy), tsa(:,:,:,jn), zdta(jn)%tra ) 67 CASE('specified' ) ; CALL bdy_spe ( idx_bdy(ib_bdy), tsa(:,:,:,jn), zdta(jn)%tra ) 68 CASE('neumann' ) ; CALL bdy_nmn ( idx_bdy(ib_bdy), igrd , tsa(:,:,:,jn) ) ! tsa masked 69 CASE('orlanski' ) ; CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), zdta(jn)%tra, ll_npo=.false. ) 70 CASE('orlanski_npo') ; CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), zdta(jn)%tra, ll_npo=.true. ) 71 CASE('runoff' ) ; CALL bdy_rnf ( idx_bdy(ib_bdy), tsa(:,:,:,jn), jn ) 72 CASE DEFAULT ; CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' ) 91 CASE('neumann','runoff') 92 llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:,ir) ! possibly every direction, T points 93 llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:,ir) ! possibly every direction, T points 94 CASE('orlanski', 'orlanski_npo') 95 llsend1(:) = llsend1(:) .OR. lsend_bdy(ib_bdy,1,:,ir) ! possibly every direction, T points 96 llrecv1(:) = llrecv1(:) .OR. lrecv_bdy(ib_bdy,1,:,ir) ! possibly every direction, T points 73 97 END SELECT 74 !75 98 END DO 99 IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN ! if need to send/recv in at least one direction 100 CALL lbc_bdy_lnk( 'bdytra', llsend1, llrecv1, tsa, 'T', 1. ) 101 END IF 76 102 END DO 77 !78 llsend1(:) = .false.79 llrecv1(:) = .false.80 DO ib_bdy=1, nb_bdy81 SELECT CASE( TRIM(cn_tra(ib_bdy)) )82 CASE('neumann','runoff')83 llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:) ! possibly every direction, T points84 llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:) ! possibly every direction, T points85 CASE('orlanski', 'orlanski_npo')86 llsend1(:) = llsend1(:) .OR. lsend_bdy(ib_bdy,1,:) ! possibly every direction, T points87 llrecv1(:) = llrecv1(:) .OR. lrecv_bdy(ib_bdy,1,:) ! possibly every direction, T points88 END SELECT89 END DO90 IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN ! if need to send/recv in at least one direction91 CALL lbc_bdy_lnk( 'bdytra', llsend1, llrecv1, tsa, 'T', 1. )92 END IF93 103 ! 94 104 END SUBROUTINE bdy_tra 95 105 96 106 97 SUBROUTINE bdy_rnf( idx, pta, jpa )107 SUBROUTINE bdy_rnf( idx, pta, jpa, llrim0 ) 98 108 !!---------------------------------------------------------------------- 99 109 !! *** SUBROUTINE bdy_rnf *** … … 104 114 !! 105 115 !!---------------------------------------------------------------------- 106 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 107 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pta ! tracer trend 108 INTEGER, INTENT(in) :: jpa ! TRA index 116 TYPE(OBC_INDEX), INTENT(in) :: idx ! OBC indices 117 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pta ! tracer trend 118 INTEGER, INTENT(in) :: jpa ! TRA index 119 LOGICAL, INTENT(in) :: llrim0 ! indicate if rim 0 is treated 109 120 ! 110 121 INTEGER :: ib, ii, ij, igrd ! dummy loop indices … … 113 124 ! 114 125 igrd = 1 ! Everything is at T-points here 115 IF( jpa == jp_tem ) THEN 116 CALL bdy_nmn( idx, igrd, pta ) 117 ELSE IF( jpa == jp_sal ) THEN 118 DO ib = 1, idx%nblenrim(igrd) 126 IF( jpa == jp_tem ) THEN 127 CALL bdy_nmn( idx, igrd, pta, llrim0 ) 128 ELSE IF( jpa == jp_sal ) THEN 129 IF( .NOT. llrim0 ) RETURN 130 DO ib = 1, idx%nblenrim(igrd) ! if llrim0 then treat the whole rim 119 131 ii = idx%nbi(ib,igrd) 120 132 ij = idx%nbj(ib,igrd) -
NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/DYN/dynkeg.F90
r11071 r11191 140 140 llrecv1(:) = .false. 141 141 DO ib_bdy = 1, nb_bdy 142 llsend1(2) = llsend1(2) .OR. lsend_bdy(ib_bdy,igrd,2) ! send east143 llsend1(4) = llsend1(4) .OR. lsend_bdy(ib_bdy,igrd,4) ! send north144 llrecv1(1) = llrecv1(1) .OR. lrecv_bdy(ib_bdy,igrd,1) ! receive west145 llrecv1(3) = llrecv1(3) .OR. lrecv_bdy(ib_bdy,igrd,3) ! receive south142 llsend1(2) = llsend1(2) .OR. ANY(lsend_bdy(ib_bdy,igrd,2,:)) ! send east 143 llsend1(4) = llsend1(4) .OR. ANY(lsend_bdy(ib_bdy,igrd,4,:)) ! send north 144 llrecv1(1) = llrecv1(1) .OR. ANY(lrecv_bdy(ib_bdy,igrd,1,:)) ! receive west 145 llrecv1(3) = llrecv1(3) .OR. ANY(lrecv_bdy(ib_bdy,igrd,3,:)) ! receive south 146 146 END DO 147 147
Note: See TracChangeset
for help on using the changeset viewer.
