- Timestamp:
- 2017-12-01T18:44:09+01:00 (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/NST_SRC/agrif_opa_sponge.F90
r7646 r8882 3 3 MODULE agrif_opa_sponge 4 4 !!====================================================================== 5 !! *** MODULE agrif_opa_update***6 !! AGRIF :5 !! *** MODULE agrif_opa_interp *** 6 !! AGRIF: sponge package for the ocean dynamics (OPA) 7 7 !!====================================================================== 8 !! History : 8 !! History : 2.0 ! 2002-06 (XXX) Original cade 9 !! - ! 2005-11 (XXX) 10 !! 3.2 ! 2009-04 (R. Benshila) 11 !! 3.6 ! 2014-09 (R. Benshila) 9 12 !!---------------------------------------------------------------------- 10 13 #if defined key_agrif 14 !!---------------------------------------------------------------------- 15 !! 'key_agrif' AGRIF zoom 16 !!---------------------------------------------------------------------- 11 17 USE par_oce 12 18 USE oce 13 19 USE dom_oce 20 ! 14 21 USE in_out_manager 15 22 USE agrif_oce 16 USE wrk_nemo17 23 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 18 24 … … 24 30 25 31 !!---------------------------------------------------------------------- 26 !! NEMO/NST 3.7 , NEMO Consortium (2015)32 !! NEMO/NST 4.0 , NEMO Consortium (2017) 27 33 !! $Id$ 28 34 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 31 37 32 38 SUBROUTINE Agrif_Sponge_Tra 33 !!--------------------------------------------- 34 !! *** ROUTINE Agrif_Sponge_Tra ***35 !!--------------------------------------------- 36 REAL(wp) :: timecoeff37 !!--------------------------------------------- 39 !!---------------------------------------------------------------------- 40 !! *** ROUTINE Agrif_Sponge_Tra *** 41 !!---------------------------------------------------------------------- 42 REAL(wp) :: zcoef ! local scalar 43 !!---------------------------------------------------------------------- 38 44 ! 39 45 #if defined SPONGE 40 timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot()41 46 zcoef = REAL( Agrif_NbStepint(), wp ) / Agrif_rhot() 47 ! 42 48 CALL Agrif_Sponge 43 Agrif_SpecialValue =0.49 Agrif_SpecialValue = 0._wp 44 50 Agrif_UseSpecialValue = .TRUE. 45 tabspongedone_tsn = .FALSE.46 47 CALL Agrif_Bc_Variable( tsn_sponge_id,calledweight=timecoeff,procname=interptsn_sponge)48 51 tabspongedone_tsn = .FALSE. 52 ! 53 CALL Agrif_Bc_Variable( tsn_sponge_id, calledweight=zcoef, procname=interptsn_sponge ) 54 ! 49 55 Agrif_UseSpecialValue = .FALSE. 50 56 #endif … … 54 60 55 61 SUBROUTINE Agrif_Sponge_dyn 56 !!--------------------------------------------- 57 !! *** ROUTINE Agrif_Sponge_dyn ***58 !!--------------------------------------------- 59 REAL(wp) :: timecoeff60 !!--------------------------------------------- 61 62 !!---------------------------------------------------------------------- 63 !! *** ROUTINE Agrif_Sponge_dyn *** 64 !!---------------------------------------------------------------------- 65 REAL(wp) :: zcoef ! local scalar 66 !!---------------------------------------------------------------------- 67 ! 62 68 #if defined SPONGE 63 timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot()64 65 Agrif_SpecialValue =0.69 zcoef = REAL( Agrif_NbStepint(), wp ) / Agrif_rhot() 70 ! 71 Agrif_SpecialValue = 0._wp 66 72 Agrif_UseSpecialValue = ln_spc_dyn 67 73 ! 68 74 tabspongedone_u = .FALSE. 69 75 tabspongedone_v = .FALSE. 70 CALL Agrif_Bc_Variable( un_sponge_id,calledweight=timecoeff,procname=interpun_sponge)71 76 CALL Agrif_Bc_Variable( un_sponge_id, calledweight=zcoef, procname=interpun_sponge ) 77 ! 72 78 tabspongedone_u = .FALSE. 73 79 tabspongedone_v = .FALSE. 74 CALL Agrif_Bc_Variable( vn_sponge_id,calledweight=timecoeff,procname=interpvn_sponge)75 80 CALL Agrif_Bc_Variable( vn_sponge_id, calledweight=zcoef, procname=interpvn_sponge ) 81 ! 76 82 Agrif_UseSpecialValue = .FALSE. 77 83 #endif … … 81 87 82 88 SUBROUTINE Agrif_Sponge 83 !!--------------------------------------------- 84 !! *** ROUTINE Agrif_Sponge ***85 !!--------------------------------------------- 86 INTEGER :: ji,jj,jk87 INTEGER :: ispongearea, ilci, ilcj88 LOGICAL :: ll_spdone89 REAL(wp) :: z1spongearea, zramp90 REAL(wp), POINTER, DIMENSION(:,:) :: ztabramp91 89 !!---------------------------------------------------------------------- 90 !! *** ROUTINE Agrif_Sponge *** 91 !!---------------------------------------------------------------------- 92 INTEGER :: ji, jj, ind1, ind2 93 INTEGER :: ispongearea 94 REAL(wp) :: z1_spongearea 95 REAL(wp), DIMENSION(jpi,jpj) :: ztabramp 96 !!---------------------------------------------------------------------- 97 ! 92 98 #if defined SPONGE || defined SPONGE_TOP 93 ll_spdone=.TRUE.94 99 IF (( .NOT. spongedoneT ).OR.( .NOT. spongedoneU )) THEN 95 ! Define ramp from boundaries towards domain interior 96 ! at T-points 100 ! Define ramp from boundaries towards domain interior at T-points 97 101 ! Store it in ztabramp 98 ll_spdone=.FALSE.99 100 CALL wrk_alloc( jpi, jpj, ztabramp )101 102 102 103 ispongearea = 2 + nn_sponge_len * Agrif_irhox() 103 ilci = nlci - ispongearea 104 ilcj = nlcj - ispongearea 105 z1spongearea = 1._wp / REAL( ispongearea - 2 ) 106 104 z1_spongearea = 1._wp / REAL( ispongearea - 1 ) 105 107 106 ztabramp(:,:) = 0._wp 108 107 108 ! --- West --- ! 109 109 IF( (nbondi == -1) .OR. (nbondi == 2) ) THEN 110 ind1 = 1+nbghostcells 111 ind2 = 1+nbghostcells + (ispongearea-1) 110 112 DO jj = 1, jpj 111 IF ( umask(2,jj,1) == 1._wp ) THEN 112 DO ji = 2, ispongearea 113 ztabramp(ji,jj) = ( ispongearea-ji ) * z1spongearea 114 END DO 115 ENDIF 113 DO ji = ind1, ind2 114 ztabramp(ji,jj) = REAL( ind2 - ji ) * z1_spongearea * umask(ind1,jj,1) 115 END DO 116 116 ENDDO 117 117 ENDIF 118 118 119 ! --- East --- ! 119 120 IF( (nbondi == 1) .OR. (nbondi == 2) ) THEN 121 ind1 = nlci - (1+nbghostcells) - (ispongearea-1) 122 ind2 = nlci - (1+nbghostcells) 120 123 DO jj = 1, jpj 121 IF ( umask(nlci-2,jj,1) == 1._wp ) THEN 122 DO ji = ilci+1,nlci-1 123 zramp = (ji - (ilci+1) ) * z1spongearea 124 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp ) 125 ENDDO 126 ENDIF 124 DO ji = ind1, ind2 125 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( ji - ind2 ) * z1_spongearea * umask(ind2-1,jj,1) ) 126 ENDDO 127 127 ENDDO 128 128 ENDIF 129 129 130 ! --- South --- ! 130 131 IF( (nbondj == -1) .OR. (nbondj == 2) ) THEN 131 DO ji = 1, jpi 132 IF ( vmask(ji,2,1) == 1._wp ) THEN 133 DO jj = 2, ispongearea 134 zramp = ( ispongearea-jj ) * z1spongearea 135 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp ) 136 END DO 137 ENDIF 132 ind1 = 1+nbghostcells 133 ind2 = 1+nbghostcells + (ispongearea-1) 134 DO jj = ind1, ind2 135 DO ji = 1, jpi 136 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( ind2 - jj ) * z1_spongearea * vmask(ji,ind1,1) ) 137 END DO 138 138 ENDDO 139 139 ENDIF 140 140 141 ! --- North --- ! 141 142 IF( (nbondj == 1) .OR. (nbondj == 2) ) THEN 142 DO ji = 1, jpi 143 IF ( vmask(ji,nlcj-2,1) == 1._wp ) THEN 144 DO jj = ilcj+1,nlcj-1 145 zramp = (jj - (ilcj+1) ) * z1spongearea 146 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp ) 147 END DO 148 ENDIF 143 ind1 = nlcj - (1+nbghostcells) - (ispongearea-1) 144 ind2 = nlcj - (1+nbghostcells) 145 DO jj = ind1, ind2 146 DO ji = 1, jpi 147 ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( jj - ind2 ) * z1_spongearea * vmask(ji,ind2-1,1) ) 148 END DO 149 149 ENDDO 150 150 ENDIF … … 158 158 DO jj = 2, jpjm1 159 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 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 165 164 CALL lbc_lnk( fsaht_spu, 'U', 1. ) ! Lateral boundary conditions 166 165 CALL lbc_lnk( fsaht_spv, 'V', 1. ) 166 167 167 spongedoneT = .TRUE. 168 168 ENDIF … … 179 179 END DO 180 180 END DO 181 182 181 CALL lbc_lnk( fsahm_spt, 'T', 1. ) ! Lateral boundary conditions 183 182 CALL lbc_lnk( fsahm_spf, 'F', 1. ) 183 184 184 spongedoneU = .TRUE. 185 185 ENDIF 186 186 ! 187 IF (.NOT.ll_spdone) CALL wrk_dealloc( jpi, jpj, ztabramp )188 !189 187 #endif 190 188 ! … … 192 190 193 191 194 SUBROUTINE interptsn_sponge( tabres,i1,i2,j1,j2,k1,k2,n1,n2,before)195 !!--------------------------------------------- 196 !! *** ROUTINE interptsn_sponge ***197 !!--------------------------------------------- 198 INTEGER , INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2199 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres200 LOGICAL , INTENT(in) ::before192 SUBROUTINE interptsn_sponge( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before ) 193 !!---------------------------------------------------------------------- 194 !! *** ROUTINE interptsn_sponge *** 195 !!---------------------------------------------------------------------- 196 INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2 197 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres 198 LOGICAL , INTENT(in ) :: before 201 199 ! 202 200 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 205 203 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ztu, ztv 206 204 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2) ::tsbdiff 205 !!---------------------------------------------------------------------- 207 206 ! 208 207 IF( before ) THEN … … 241 240 zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 242 241 ! horizontal diffusive trends 243 ztsa = zbtr * ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) + ztv(ji,jj,jk) - ztv(ji 242 ztsa = zbtr * ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) + ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) 244 243 ! add it to the general tracer trends 245 244 tsa(ji,jj,jk,jn) = tsa(ji,jj,jk,jn) + ztsa … … 258 257 259 258 260 SUBROUTINE interpun_sponge(tabres,i1,i2,j1,j2,k1,k2, before) 261 !!--------------------------------------------- 262 !! *** ROUTINE interpun_sponge *** 263 !!--------------------------------------------- 264 INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2 265 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 266 LOGICAL, INTENT(in) :: before 267 268 INTEGER :: ji,jj,jk 269 259 SUBROUTINE interpun_sponge( tabres, i1, i2, j1, j2, k1, k2, before ) 260 !!---------------------------------------------------------------------- 261 !! *** ROUTINE interpun_sponge *** 262 !!---------------------------------------------------------------------- 263 INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2 264 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 265 LOGICAL , INTENT(in ) :: before 266 !! 267 INTEGER :: ji, jj, jk 270 268 ! sponge parameters 269 INTEGER :: jmax 271 270 REAL(wp) :: ze2u, ze1v, zua, zva, zbtr 272 271 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ubdiff 273 272 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: rotdiff, hdivdiff 274 INTEGER :: jmax 275 !!--------------------------------------------- 273 !!---------------------------------------------------------------------- 276 274 ! 277 275 IF( before ) THEN 278 276 tabres = un(i1:i2,j1:j2,:) 279 277 ELSE 280 ubdiff(i1:i2,j1:j2,:) = ( ub(i1:i2,j1:j2,:) - tabres(:,:,:))*umask(i1:i2,j1:j2,:)278 ubdiff(i1:i2,j1:j2,:) = ( ub(i1:i2,j1:j2,:) - tabres(:,:,:) )*umask(i1:i2,j1:j2,:) 281 279 ! 282 280 DO jk = 1, jpkm1 ! Horizontal slab … … 297 295 DO ji = i1,i2 ! vector opt. 298 296 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 297 rotdiff(ji,jj,jk) = ( -e1u(ji,jj+1) * ubdiff(ji,jj+1,jk) & 298 & +e1u(ji,jj ) * ubdiff(ji,jj ,jk) ) * fmask(ji,jj,jk) * zbtr 302 299 END DO 303 300 END DO … … 312 309 ze1v = hdivdiff(ji,jj,jk) 313 310 ! 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)311 zua = - ( ze2u - rotdiff (ji,jj-1,jk) ) / ( e2u(ji,jj) * e3u_n(ji,jj,jk) ) & 312 + ( hdivdiff(ji+1,jj,jk) - ze1v ) * r1_e1u(ji,jj) 316 313 317 314 ! add it to the general momentum trends … … 327 324 328 325 jmax = j2-1 329 IF ((nbondj == 1).OR.(nbondj == 2)) jmax = MIN(jmax,nlcj- 3)326 IF ((nbondj == 1).OR.(nbondj == 2)) jmax = MIN(jmax,nlcj-nbghostcells-2) ! North 330 327 331 328 DO jj = j1+1, jmax … … 338 335 339 336 ! 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)337 zva = + ( ze2u - rotdiff (ji-1,jj,jk) ) / ( e1v(ji,jj) * e3v_n(ji,jj,jk) ) & 338 + ( hdivdiff(ji,jj+1,jk) - ze1v ) * r1_e2v(ji,jj) 342 339 343 340 ! add it to the general momentum trends … … 356 353 357 354 358 SUBROUTINE interpvn_sponge(tabres,i1,i2,j1,j2,k1,k2, before,nb,ndir) 359 !!--------------------------------------------- 360 !! *** ROUTINE interpvn_sponge *** 361 !!--------------------------------------------- 362 INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2 363 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 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 !!--------------------------------------------- 355 SUBROUTINE interpvn_sponge( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir ) 356 !!---------------------------------------------------------------------- 357 !! *** ROUTINE interpvn_sponge *** 358 !!---------------------------------------------------------------------- 359 INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2 360 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 361 LOGICAL , INTENT(in ) :: before 362 INTEGER , INTENT(in ) :: nb , ndir 363 ! 364 INTEGER :: ji, jj, jk 365 INTEGER :: imax 366 REAL(wp):: ze2u, ze1v, zua, zva, zbtr 367 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: vbdiff, rotdiff, hdivdiff 368 !!---------------------------------------------------------------------- 373 369 374 370 IF( before ) THEN … … 376 372 ELSE 377 373 ! 378 vbdiff(i1:i2,j1:j2,:) = ( vb(i1:i2,j1:j2,:) - tabres(:,:,:))*vmask(i1:i2,j1:j2,:)374 vbdiff(i1:i2,j1:j2,:) = ( vb(i1:i2,j1:j2,:) - tabres(:,:,:) ) * vmask(i1:i2,j1:j2,:) 379 375 ! 380 376 DO jk = 1, jpkm1 ! Horizontal slab … … 403 399 ! 404 400 405 imax = i2 -1406 IF ((nbondi == 1).OR.(nbondi == 2)) imax = MIN(imax,nlci- 3)401 imax = i2 - 1 402 IF ((nbondi == 1).OR.(nbondi == 2)) imax = MIN(imax,nlci-nbghostcells-2) ! East 407 403 408 404 DO jj = j1+1, j2 … … 437 433 438 434 #else 435 !!---------------------------------------------------------------------- 436 !! Empty module no AGRIF zoom 437 !!---------------------------------------------------------------------- 439 438 CONTAINS 440 439 SUBROUTINE agrif_opa_sponge_empty 441 !!---------------------------------------------442 !! *** ROUTINE agrif_OPA_sponge_empty ***443 !!---------------------------------------------444 440 WRITE(*,*) 'agrif_opa_sponge : You should not have seen this print! error?' 445 441 END SUBROUTINE agrif_opa_sponge_empty
Note: See TracChangeset
for help on using the changeset viewer.