- Timestamp:
- 2016-07-19T10:38:35+02:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/NERC/dev_r5549_BDY_ZEROGRAD/NEMOGCM/NEMO/NST_SRC/agrif_opa_sponge.F90
r4153 r6808 1 1 #define SPONGE && define SPONGE_TOP 2 2 3 Module agrif_opa_sponge 3 MODULE agrif_opa_sponge 4 !!====================================================================== 5 !! *** MODULE agrif_opa_update *** 6 !! AGRIF : 7 !!====================================================================== 8 !! History : 9 !!---------------------------------------------------------------------- 4 10 #if defined key_agrif && ! defined key_offline 5 11 USE par_oce … … 9 15 USE agrif_oce 10 16 USE wrk_nemo 17 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 11 18 12 19 IMPLICIT NONE 13 20 PRIVATE 14 21 15 PUBLIC Agrif_Sponge, Agrif_Sponge_Tra, Agrif_Sponge_Dyn, interptsn, interpun, interpvn 16 17 !! * Substitutions 18 # include "domzgr_substitute.h90" 22 PUBLIC Agrif_Sponge, Agrif_Sponge_Tra, Agrif_Sponge_Dyn 23 PUBLIC interptsn_sponge, interpun_sponge, interpvn_sponge 24 19 25 !!---------------------------------------------------------------------- 20 !! NEMO/NST 3. 3 , NEMO Consortium (2010)26 !! NEMO/NST 3.7 , NEMO Consortium (2015) 21 27 !! $Id$ 22 28 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 23 29 !!---------------------------------------------------------------------- 24 25 CONTAINS 30 CONTAINS 26 31 27 32 SUBROUTINE Agrif_Sponge_Tra … … 29 34 !! *** ROUTINE Agrif_Sponge_Tra *** 30 35 !!--------------------------------------------- 31 !!32 INTEGER :: ji,jj,jk,jn33 36 REAL(wp) :: timecoeff 34 REAL(wp) :: ztsa, zabe1, zabe2, zbtr 35 REAL(wp), POINTER, DIMENSION(:,: ) :: ztu, ztv 36 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: ztab 37 REAL(wp), POINTER, DIMENSION(:,:,:,:) :: tsbdiff 38 37 !!--------------------------------------------- 38 ! 39 39 #if defined SPONGE 40 CALL wrk_alloc( jpi, jpj, ztu, ztv )41 CALL wrk_alloc( jpi, jpj, jpk, jpts, ztab, tsbdiff )42 43 40 timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot() 44 41 42 CALL Agrif_Sponge 45 43 Agrif_SpecialValue=0. 46 44 Agrif_UseSpecialValue = .TRUE. 47 ztab = 0.e0 48 CALL Agrif_Bc_Variable(ztab, tsa_id,calledweight=timecoeff,procname=interptsn) 45 tabspongedone_tsn = .FALSE. 46 47 CALL Agrif_Bc_Variable(tsn_sponge_id,calledweight=timecoeff,procname=interptsn_sponge) 48 49 49 Agrif_UseSpecialValue = .FALSE. 50 51 tsbdiff(:,:,:,:) = tsb(:,:,:,:) - ztab(:,:,:,:)52 53 CALL Agrif_Sponge54 55 DO jn = 1, jpts56 DO jk = 1, jpkm157 !58 DO jj = 1, jpjm159 DO ji = 1, jpim160 zabe1 = umask(ji,jj,jk) * spe1ur(ji,jj) * fse3u(ji,jj,jk)61 zabe2 = vmask(ji,jj,jk) * spe2vr(ji,jj) * fse3v(ji,jj,jk)62 ztu(ji,jj) = zabe1 * ( tsbdiff(ji+1,jj ,jk,jn) - tsbdiff(ji,jj,jk,jn) )63 ztv(ji,jj) = zabe2 * ( tsbdiff(ji ,jj+1,jk,jn) - tsbdiff(ji,jj,jk,jn) )64 ENDDO65 ENDDO66 67 DO jj = 2, jpjm168 DO ji = 2, jpim169 zbtr = spbtr2(ji,jj) / fse3t(ji,jj,jk)70 ! horizontal diffusive trends71 ztsa = zbtr * ( ztu(ji,jj) - ztu(ji-1,jj ) &72 & + ztv(ji,jj) - ztv(ji ,jj-1) )73 ! add it to the general tracer trends74 tsa(ji,jj,jk,jn) = tsa(ji,jj,jk,jn) + ztsa75 END DO76 END DO77 !78 ENDDO79 ENDDO80 81 CALL wrk_dealloc( jpi, jpj, ztu, ztv )82 CALL wrk_dealloc( jpi, jpj, jpk, jpts, ztab, tsbdiff )83 50 #endif 84 51 ! 85 52 END SUBROUTINE Agrif_Sponge_Tra 86 53 54 87 55 SUBROUTINE Agrif_Sponge_dyn 88 56 !!--------------------------------------------- 89 57 !! *** ROUTINE Agrif_Sponge_dyn *** 90 58 !!--------------------------------------------- 91 !!92 INTEGER :: ji,jj,jk93 59 REAL(wp) :: timecoeff 94 REAL(wp) :: ze2u, ze1v, zua, zva, zbtr 95 REAL(wp), POINTER, DIMENSION(:,:,:) :: ubdiff, vbdiff 96 REAL(wp), POINTER, DIMENSION(:,:,:) :: rotdiff, hdivdiff 97 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztab 60 !!--------------------------------------------- 98 61 99 62 #if defined SPONGE 100 CALL wrk_alloc( jpi, jpj, jpk, ztab, ubdiff, vbdiff, rotdiff, hdivdiff )101 102 63 timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot() 103 64 104 65 Agrif_SpecialValue=0. 105 66 Agrif_UseSpecialValue = ln_spc_dyn 106 ztab = 0.e0 107 CALL Agrif_Bc_Variable(ztab, ua_id,calledweight=timecoeff,procname=interpun) 67 68 tabspongedone_u = .FALSE. 69 tabspongedone_v = .FALSE. 70 CALL Agrif_Bc_Variable(un_sponge_id,calledweight=timecoeff,procname=interpun_sponge) 71 72 tabspongedone_u = .FALSE. 73 tabspongedone_v = .FALSE. 74 CALL Agrif_Bc_Variable(vn_sponge_id,calledweight=timecoeff,procname=interpvn_sponge) 75 108 76 Agrif_UseSpecialValue = .FALSE. 109 110 ubdiff(:,:,:) = ( ub(:,:,:) - ztab(:,:,:) ) * umask(:,:,:)111 112 ztab = 0.e0113 Agrif_SpecialValue=0.114 Agrif_UseSpecialValue = ln_spc_dyn115 CALL Agrif_Bc_Variable(ztab, va_id,calledweight=timecoeff,procname=interpvn)116 Agrif_UseSpecialValue = .FALSE.117 118 vbdiff(:,:,:) = ( vb(:,:,:) - ztab(:,:,:) ) * vmask(:,:,:)119 120 CALL Agrif_Sponge121 122 DO jk = 1,jpkm1123 ubdiff(:,:,jk) = ubdiff(:,:,jk) * spe1ur2(:,:)124 vbdiff(:,:,jk) = vbdiff(:,:,jk) * spe2vr2(:,:)125 ENDDO126 127 hdivdiff = 0.128 rotdiff = 0.129 130 DO jk = 1, jpkm1 ! Horizontal slab131 ! ! ===============132 133 ! ! --------134 ! Horizontal divergence ! div135 ! ! --------136 DO jj = 2, jpjm1137 DO ji = 2, jpim1 ! vector opt.138 zbtr = spbtr2(ji,jj) / fse3t(ji,jj,jk)139 hdivdiff(ji,jj,jk) = ( e2u(ji ,jj ) * fse3u(ji ,jj ,jk) * ubdiff(ji ,jj ,jk) &140 & - e2u(ji-1,jj ) * fse3u(ji-1,jj ,jk) * ubdiff(ji-1,jj ,jk) &141 & + e1v(ji ,jj ) * fse3v(ji ,jj ,jk) * vbdiff(ji ,jj ,jk) &142 & - e1v(ji ,jj-1) * fse3v(ji ,jj-1,jk) * vbdiff(ji ,jj-1,jk) ) * zbtr143 END DO144 END DO145 146 DO jj = 1, jpjm1147 DO ji = 1, jpim1 ! vector opt.148 zbtr = spbtr3(ji,jj) * fse3f(ji,jj,jk)149 rotdiff(ji,jj,jk) = ( e2v(ji+1,jj ) * vbdiff(ji+1,jj ,jk) - e2v(ji,jj) * vbdiff(ji,jj,jk) &150 & - e1u(ji ,jj+1) * ubdiff(ji ,jj+1,jk) + e1u(ji,jj) * ubdiff(ji,jj,jk) ) &151 & * fmask(ji,jj,jk) * zbtr152 END DO153 END DO154 155 ENDDO156 157 ! ! ===============158 DO jk = 1, jpkm1 ! Horizontal slab159 ! ! ===============160 DO jj = 2, jpjm1161 DO ji = 2, jpim1 ! vector opt.162 ! horizontal diffusive trends163 zua = - ( rotdiff (ji ,jj,jk) - rotdiff (ji,jj-1,jk) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) ) &164 + ( hdivdiff(ji+1,jj,jk) - hdivdiff(ji,jj ,jk) ) / e1u(ji,jj)165 166 zva = + ( rotdiff (ji,jj ,jk) - rotdiff (ji-1,jj,jk) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) ) &167 + ( hdivdiff(ji,jj+1,jk) - hdivdiff(ji ,jj,jk) ) / e2v(ji,jj)168 ! add it to the general momentum trends169 ua(ji,jj,jk) = ua(ji,jj,jk) + zua170 va(ji,jj,jk) = va(ji,jj,jk) + zva171 END DO172 END DO173 ! ! ===============174 END DO ! End of slab175 ! ! ===============176 CALL wrk_dealloc( jpi, jpj, jpk, ztab, ubdiff, vbdiff, rotdiff, hdivdiff )177 77 #endif 178 78 ! 179 79 END SUBROUTINE Agrif_Sponge_dyn 80 180 81 181 82 SUBROUTINE Agrif_Sponge … … 199 100 CALL wrk_alloc( jpi, jpj, ztabramp ) 200 101 201 ispongearea = 2 + 2* Agrif_irhox()102 ispongearea = 2 + nn_sponge_len * Agrif_irhox() 202 103 ilci = nlci - ispongearea 203 104 ilcj = nlcj - ispongearea 204 105 z1spongearea = 1._wp / REAL( ispongearea - 2 ) 205 spbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:) ) 206 207 ztabramp(:,:) = 0. 106 107 ztabramp(:,:) = 0._wp 208 108 209 109 IF( (nbondi == -1) .OR. (nbondi == 2) ) THEN … … 254 154 ! Tracers 255 155 IF( .NOT. spongedoneT ) THEN 256 spe1ur(:,:) = 0. 257 spe2vr(:,:) = 0. 258 259 IF( (nbondi == -1) .OR. (nbondi == 2) ) THEN 260 spe1ur(2:ispongearea-1,: ) = visc_tra & 261 & * 0.5 * ( ztabramp(2:ispongearea-1,: ) & 262 & + ztabramp(3:ispongearea ,: ) ) & 263 & * e2u(2:ispongearea-1,:) / e1u(2:ispongearea-1,:) 264 265 spe2vr(2:ispongearea ,1:jpjm1 ) = visc_tra & 266 & * 0.5 * ( ztabramp(2:ispongearea ,1:jpjm1) & 267 & + ztabramp(2:ispongearea,2 :jpj ) ) & 268 & * e1v(2:ispongearea,1:jpjm1) / e2v(2:ispongearea,1:jpjm1) 269 ENDIF 270 271 IF( (nbondi == 1) .OR. (nbondi == 2) ) THEN 272 spe1ur(ilci+1:nlci-2,: ) = visc_tra & 273 & * 0.5 * ( ztabramp(ilci+1:nlci-2,: ) & 274 & + ztabramp(ilci+2:nlci-1,: ) ) & 275 & * e2u(ilci+1:nlci-2,:) / e1u(ilci+1:nlci-2,:) 276 277 spe2vr(ilci+1:nlci-1,1:jpjm1 ) = visc_tra & 278 & * 0.5 * ( ztabramp(ilci+1:nlci-1,1:jpjm1) & 279 & + ztabramp(ilci+1:nlci-1,2:jpj ) ) & 280 & * e1v(ilci+1:nlci-1,1:jpjm1) / e2v(ilci+1:nlci-1,1:jpjm1) 281 ENDIF 282 283 IF( (nbondj == -1) .OR. (nbondj == 2) ) THEN 284 spe1ur(1:jpim1,2:ispongearea ) = visc_tra & 285 & * 0.5 * ( ztabramp(1:jpim1,2:ispongearea ) & 286 & + ztabramp(2:jpi ,2:ispongearea ) ) & 287 & * e2u(1:jpim1,2:ispongearea) / e1u(1:jpim1,2:ispongearea) 288 289 spe2vr(: ,2:ispongearea-1) = visc_tra & 290 & * 0.5 * ( ztabramp(: ,2:ispongearea-1) & 291 & + ztabramp(: ,3:ispongearea ) ) & 292 & * e1v(:,2:ispongearea-1) / e2v(:,2:ispongearea-1) 293 ENDIF 294 295 IF( (nbondj == 1) .OR. (nbondj == 2) ) THEN 296 spe1ur(1:jpim1,ilcj+1:nlcj-1) = visc_tra & 297 & * 0.5 * ( ztabramp(1:jpim1,ilcj+1:nlcj-1) & 298 & + ztabramp(2:jpi ,ilcj+1:nlcj-1) ) & 299 & * e2u(1:jpim1,ilcj+1:nlcj-1) / e1u(1:jpim1,ilcj+1:nlcj-1) 300 301 spe2vr(: ,ilcj+1:nlcj-2) = visc_tra & 302 & * 0.5 * ( ztabramp(: ,ilcj+1:nlcj-2) & 303 & + ztabramp(: ,ilcj+2:nlcj-1) ) & 304 & * e1v(:,ilcj+1:nlcj-2) / e2v(:,ilcj+1:nlcj-2) 305 ENDIF 156 fsaht_spu(:,:) = 0._wp 157 fsaht_spv(:,:) = 0._wp 158 DO jj = 2, jpjm1 159 DO ji = 2, jpim1 ! vector opt. 160 fsaht_spu(ji,jj) = 0.5_wp * visc_tra * (ztabramp(ji,jj) + ztabramp(ji+1,jj )) 161 fsaht_spv(ji,jj) = 0.5_wp * visc_tra * (ztabramp(ji,jj) + ztabramp(ji ,jj+1)) 162 END DO 163 END DO 164 165 CALL lbc_lnk( fsaht_spu, 'U', 1. ) ! Lateral boundary conditions 166 CALL lbc_lnk( fsaht_spv, 'V', 1. ) 306 167 spongedoneT = .TRUE. 307 168 ENDIF … … 309 170 ! Dynamics 310 171 IF( .NOT. spongedoneU ) THEN 311 spe1ur2(:,:) = 0. 312 spe2vr2(:,:) = 0. 313 314 IF( (nbondi == -1) .OR. (nbondi == 2) ) THEN 315 spe1ur2(2:ispongearea-1,: ) = visc_dyn & 316 & * 0.5 * ( ztabramp(2:ispongearea-1,: ) & 317 & + ztabramp(3:ispongearea ,: ) ) 318 spe2vr2(2:ispongearea ,1:jpjm1) = visc_dyn & 319 & * 0.5 * ( ztabramp(2:ispongearea ,1:jpjm1) & 320 & + ztabramp(2:ispongearea ,2:jpj ) ) 321 ENDIF 322 323 IF( (nbondi == 1) .OR. (nbondi == 2) ) THEN 324 spe1ur2(ilci+1:nlci-2 ,: ) = visc_dyn & 325 & * 0.5 * ( ztabramp(ilci+1:nlci-2, : ) & 326 & + ztabramp(ilci+2:nlci-1, : ) ) 327 spe2vr2(ilci+1:nlci-1 ,1:jpjm1) = visc_dyn & 328 & * 0.5 * ( ztabramp(ilci+1:nlci-1,1:jpjm1 ) & 329 & + ztabramp(ilci+1:nlci-1,2:jpj ) ) 330 ENDIF 331 332 IF( (nbondj == -1) .OR. (nbondj == 2) ) THEN 333 spe1ur2(1:jpim1,2:ispongearea ) = visc_dyn & 334 & * 0.5 * ( ztabramp(1:jpim1,2:ispongearea ) & 335 & + ztabramp(2:jpi ,2:ispongearea ) ) 336 spe2vr2(: ,2:ispongearea-1) = visc_dyn & 337 & * 0.5 * ( ztabramp(: ,2:ispongearea-1) & 338 & + ztabramp(: ,3:ispongearea ) ) 339 ENDIF 340 341 IF( (nbondj == 1) .OR. (nbondj == 2) ) THEN 342 spe1ur2(1:jpim1,ilcj+1:nlcj-1 ) = visc_dyn & 343 & * 0.5 * ( ztabramp(1:jpim1,ilcj+1:nlcj-1 ) & 344 & + ztabramp(2:jpi ,ilcj+1:nlcj-1 ) ) 345 spe2vr2(: ,ilcj+1:nlcj-2 ) = visc_dyn & 346 & * 0.5 * ( ztabramp(: ,ilcj+1:nlcj-2 ) & 347 & + ztabramp(: ,ilcj+2:nlcj-1 ) ) 348 ENDIF 172 fsahm_spt(:,:) = 0._wp 173 fsahm_spf(:,:) = 0._wp 174 DO jj = 2, jpjm1 175 DO ji = 2, jpim1 ! vector opt. 176 fsahm_spt(ji,jj) = visc_dyn * ztabramp(ji,jj) 177 fsahm_spf(ji,jj) = 0.25_wp * visc_dyn * ( ztabramp(ji,jj) + ztabramp(ji ,jj+1) & 178 & +ztabramp(ji,jj) + ztabramp(ji+1,jj ) ) 179 END DO 180 END DO 181 182 CALL lbc_lnk( fsahm_spt, 'T', 1. ) ! Lateral boundary conditions 183 CALL lbc_lnk( fsahm_spf, 'F', 1. ) 349 184 spongedoneU = .TRUE. 350 spbtr3(:,:) = 1. / ( e1f(:,:) * e2f(:,:) )351 185 ENDIF 352 186 ! … … 354 188 ! 355 189 #endif 356 190 ! 357 191 END SUBROUTINE Agrif_Sponge 358 192 359 SUBROUTINE interptsn(tabres,i1,i2,j1,j2,k1,k2,n1,n2) 360 !!--------------------------------------------- 361 !! *** ROUTINE interptsn *** 193 194 SUBROUTINE interptsn_sponge(tabres,i1,i2,j1,j2,k1,k2,n1,n2,before) 195 !!--------------------------------------------- 196 !! *** ROUTINE interptsn_sponge *** 362 197 !!--------------------------------------------- 363 198 INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2 364 199 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres 365 366 tabres(i1:i2,j1:j2,k1:k2,n1:n2) = tsn(i1:i2,j1:j2,k1:k2,n1:n2) 367 368 END SUBROUTINE interptsn 369 370 SUBROUTINE interpun(tabres,i1,i2,j1,j2,k1,k2) 371 !!--------------------------------------------- 372 !! *** ROUTINE interpun *** 373 !!--------------------------------------------- 200 LOGICAL, INTENT(in) :: before 201 ! 202 INTEGER :: ji, jj, jk, jn ! dummy loop indices 203 INTEGER :: iku, ikv 204 REAL(wp) :: ztsa, zabe1, zabe2, zbtr 205 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ztu, ztv 206 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2) ::tsbdiff 207 ! 208 IF( before ) THEN 209 tabres(i1:i2,j1:j2,k1:k2,n1:n2) = tsn(i1:i2,j1:j2,k1:k2,n1:n2) 210 ELSE 211 ! 212 tsbdiff(:,:,:,:) = tsb(i1:i2,j1:j2,:,:) - tabres(:,:,:,:) 213 DO jn = 1, jpts 214 DO jk = 1, jpkm1 215 DO jj = j1,j2-1 216 DO ji = i1,i2-1 217 zabe1 = fsaht_spu(ji,jj) * umask(ji,jj,jk) * e2_e1u(ji,jj) * e3u_n(ji,jj,jk) 218 zabe2 = fsaht_spv(ji,jj) * vmask(ji,jj,jk) * e1_e2v(ji,jj) * e3v_n(ji,jj,jk) 219 ztu(ji,jj,jk) = zabe1 * ( tsbdiff(ji+1,jj ,jk,jn) - tsbdiff(ji,jj,jk,jn) ) 220 ztv(ji,jj,jk) = zabe2 * ( tsbdiff(ji ,jj+1,jk,jn) - tsbdiff(ji,jj,jk,jn) ) 221 END DO 222 END DO 223 ! 224 IF( ln_zps ) THEN ! set gradient at partial step level 225 DO jj = j1,j2-1 226 DO ji = i1,i2-1 227 ! last level 228 iku = mbku(ji,jj) 229 ikv = mbkv(ji,jj) 230 IF( iku == jk ) ztu(ji,jj,jk) = 0._wp 231 IF( ikv == jk ) ztv(ji,jj,jk) = 0._wp 232 END DO 233 END DO 234 ENDIF 235 END DO 236 ! 237 DO jk = 1, jpkm1 238 DO jj = j1+1,j2-1 239 DO ji = i1+1,i2-1 240 IF (.NOT. tabspongedone_tsn(ji,jj)) THEN 241 zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 242 ! horizontal diffusive trends 243 ztsa = zbtr * ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) + ztv(ji,jj,jk) - ztv(ji ,jj-1,jk) ) 244 ! add it to the general tracer trends 245 tsa(ji,jj,jk,jn) = tsa(ji,jj,jk,jn) + ztsa 246 ENDIF 247 END DO 248 END DO 249 END DO 250 ! 251 END DO 252 ! 253 tabspongedone_tsn(i1+1:i2-1,j1+1:j2-1) = .TRUE. 254 ! 255 ENDIF 256 ! 257 END SUBROUTINE interptsn_sponge 258 259 260 SUBROUTINE interpun_sponge(tabres,i1,i2,j1,j2,k1,k2, before) 261 !!--------------------------------------------- 262 !! *** ROUTINE interpun_sponge *** 263 !!--------------------------------------------- 374 264 INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2 375 265 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 376 377 tabres(i1:i2,j1:j2,k1:k2) = un(i1:i2,j1:j2,k1:k2) 378 379 END SUBROUTINE interpun 380 381 SUBROUTINE interpvn(tabres,i1,i2,j1,j2,k1,k2) 382 !!--------------------------------------------- 383 !! *** ROUTINE interpvn *** 384 !!--------------------------------------------- 266 LOGICAL, INTENT(in) :: before 267 268 INTEGER :: ji,jj,jk 269 270 ! sponge parameters 271 REAL(wp) :: ze2u, ze1v, zua, zva, zbtr 272 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ubdiff 273 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: rotdiff, hdivdiff 274 INTEGER :: jmax 275 !!--------------------------------------------- 276 ! 277 IF( before ) THEN 278 tabres = un(i1:i2,j1:j2,:) 279 ELSE 280 ubdiff(i1:i2,j1:j2,:) = (ub(i1:i2,j1:j2,:) - tabres(:,:,:))*umask(i1:i2,j1:j2,:) 281 ! 282 DO jk = 1, jpkm1 ! Horizontal slab 283 ! ! =============== 284 285 ! ! -------- 286 ! Horizontal divergence ! div 287 ! ! -------- 288 DO jj = j1,j2 289 DO ji = i1+1,i2 ! vector opt. 290 zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) * fsahm_spt(ji,jj) 291 hdivdiff(ji,jj,jk) = ( e2u(ji ,jj)*e3u_n(ji ,jj,jk) * ubdiff(ji ,jj,jk) & 292 & -e2u(ji-1,jj)*e3u_n(ji-1,jj,jk) * ubdiff(ji-1,jj,jk) ) * zbtr 293 END DO 294 END DO 295 296 DO jj = j1,j2-1 297 DO ji = i1,i2 ! vector opt. 298 zbtr = r1_e1e2f(ji,jj) * e3f_n(ji,jj,jk) * fsahm_spf(ji,jj) 299 rotdiff(ji,jj,jk) = (-e1u(ji,jj+1) * ubdiff(ji,jj+1,jk) & 300 +e1u(ji,jj ) * ubdiff(ji,jj ,jk) & 301 & ) * fmask(ji,jj,jk) * zbtr 302 END DO 303 END DO 304 END DO 305 ! 306 DO jj = j1+1, j2-1 307 DO ji = i1+1, i2-1 ! vector opt. 308 309 IF (.NOT. tabspongedone_u(ji,jj)) THEN 310 DO jk = 1, jpkm1 ! Horizontal slab 311 ze2u = rotdiff (ji,jj,jk) 312 ze1v = hdivdiff(ji,jj,jk) 313 ! horizontal diffusive trends 314 zua = - ( ze2u - rotdiff (ji,jj-1,jk)) / ( e2u(ji,jj) * e3u_n(ji,jj,jk) ) & 315 + ( hdivdiff(ji+1,jj,jk) - ze1v ) / e1u(ji,jj) 316 317 ! add it to the general momentum trends 318 ua(ji,jj,jk) = ua(ji,jj,jk) + zua 319 320 END DO 321 ENDIF 322 323 END DO 324 END DO 325 326 tabspongedone_u(i1+1:i2-1,j1+1:j2-1) = .TRUE. 327 328 jmax = j2-1 329 IF ((nbondj == 1).OR.(nbondj == 2)) jmax = MIN(jmax,nlcj-3) 330 331 DO jj = j1+1, jmax 332 DO ji = i1+1, i2 ! vector opt. 333 334 IF (.NOT. tabspongedone_v(ji,jj)) THEN 335 DO jk = 1, jpkm1 ! Horizontal slab 336 ze2u = rotdiff (ji,jj,jk) 337 ze1v = hdivdiff(ji,jj,jk) 338 339 ! horizontal diffusive trends 340 zva = + ( ze2u - rotdiff (ji-1,jj,jk)) / ( e1v(ji,jj) * e3v_n(ji,jj,jk) ) & 341 + ( hdivdiff(ji,jj+1,jk) - ze1v ) / e2v(ji,jj) 342 343 ! add it to the general momentum trends 344 va(ji,jj,jk) = va(ji,jj,jk) + zva 345 END DO 346 ENDIF 347 ! 348 END DO 349 END DO 350 ! 351 tabspongedone_v(i1+1:i2,j1+1:jmax) = .TRUE. 352 ! 353 ENDIF 354 ! 355 END SUBROUTINE interpun_sponge 356 357 358 SUBROUTINE interpvn_sponge(tabres,i1,i2,j1,j2,k1,k2, before,nb,ndir) 359 !!--------------------------------------------- 360 !! *** ROUTINE interpvn_sponge *** 361 !!--------------------------------------------- 385 362 INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2 386 363 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 387 388 tabres(i1:i2,j1:j2,k1:k2) = vn(i1:i2,j1:j2,k1:k2) 389 390 END SUBROUTINE interpvn 364 LOGICAL, INTENT(in) :: before 365 INTEGER, INTENT(in) :: nb , ndir 366 ! 367 INTEGER :: ji, jj, jk 368 REAL(wp) :: ze2u, ze1v, zua, zva, zbtr 369 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: vbdiff 370 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: rotdiff, hdivdiff 371 INTEGER :: imax 372 !!--------------------------------------------- 373 374 IF( before ) THEN 375 tabres = vn(i1:i2,j1:j2,:) 376 ELSE 377 ! 378 vbdiff(i1:i2,j1:j2,:) = (vb(i1:i2,j1:j2,:) - tabres(:,:,:))*vmask(i1:i2,j1:j2,:) 379 ! 380 DO jk = 1, jpkm1 ! Horizontal slab 381 ! ! =============== 382 383 ! ! -------- 384 ! Horizontal divergence ! div 385 ! ! -------- 386 DO jj = j1+1,j2 387 DO ji = i1,i2 ! vector opt. 388 zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) * fsahm_spt(ji,jj) 389 hdivdiff(ji,jj,jk) = ( e1v(ji,jj ) * e3v_n(ji,jj ,jk) * vbdiff(ji,jj ,jk) & 390 & -e1v(ji,jj-1) * e3v_n(ji,jj-1,jk) * vbdiff(ji,jj-1,jk) ) * zbtr 391 END DO 392 END DO 393 DO jj = j1,j2 394 DO ji = i1,i2-1 ! vector opt. 395 zbtr = r1_e1e2f(ji,jj) * e3f_n(ji,jj,jk) * fsahm_spf(ji,jj) 396 rotdiff(ji,jj,jk) = ( e2v(ji+1,jj) * vbdiff(ji+1,jj,jk) & 397 & -e2v(ji ,jj) * vbdiff(ji ,jj,jk) ) * fmask(ji,jj,jk) * zbtr 398 END DO 399 END DO 400 END DO 401 402 ! ! =============== 403 ! 404 405 imax = i2-1 406 IF ((nbondi == 1).OR.(nbondi == 2)) imax = MIN(imax,nlci-3) 407 408 DO jj = j1+1, j2 409 DO ji = i1+1, imax ! vector opt. 410 IF( .NOT. tabspongedone_u(ji,jj) ) THEN 411 DO jk = 1, jpkm1 412 ua(ji,jj,jk) = ua(ji,jj,jk) & 413 & - ( rotdiff (ji ,jj,jk) - rotdiff (ji,jj-1,jk)) / ( e2u(ji,jj) * e3u_n(ji,jj,jk) ) & 414 & + ( hdivdiff(ji+1,jj,jk) - hdivdiff(ji,jj ,jk)) * r1_e1u(ji,jj) 415 END DO 416 ENDIF 417 END DO 418 END DO 419 ! 420 tabspongedone_u(i1+1:imax,j1+1:j2) = .TRUE. 421 ! 422 DO jj = j1+1, j2-1 423 DO ji = i1+1, i2-1 ! vector opt. 424 IF( .NOT. tabspongedone_v(ji,jj) ) THEN 425 DO jk = 1, jpkm1 426 va(ji,jj,jk) = va(ji,jj,jk) & 427 & + ( rotdiff (ji,jj ,jk) - rotdiff (ji-1,jj,jk) ) / ( e1v(ji,jj) * e3v_n(ji,jj,jk) ) & 428 & + ( hdivdiff(ji,jj+1,jk) - hdivdiff(ji ,jj,jk) ) * r1_e2v(ji,jj) 429 END DO 430 ENDIF 431 END DO 432 END DO 433 tabspongedone_v(i1+1:i2-1,j1+1:j2-1) = .TRUE. 434 ENDIF 435 ! 436 END SUBROUTINE interpvn_sponge 391 437 392 438 #else 393 439 CONTAINS 394 395 440 SUBROUTINE agrif_opa_sponge_empty 396 441 !!--------------------------------------------- … … 401 446 #endif 402 447 448 !!====================================================================== 403 449 END MODULE agrif_opa_sponge
Note: See TracChangeset
for help on using the changeset viewer.
