Changeset 7514 for branches/2016/dev_merge_2016
- Timestamp:
- 2016-12-20T15:40:04+01:00 (8 years ago)
- Location:
- branches/2016/dev_merge_2016/NEMOGCM
- Files:
-
- 8 deleted
- 11 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_merge_2016/NEMOGCM/CONFIG/WAD_TEST_CASES/EXP00/README.py
r7412 r7514 1 ncks -d y,17,17,1 WAD_1ts_00010101_00010101_grid_T.nc sshtime.nc 2 ncks -4 -A -v gdepw_0,bathy -d y,17,17,1 -d z,9,9,1 -C mesh_mask_wad1.nc sshtime.nc 1 ncks -O -d y,17,17,1 WAD_1ts_00010101_00010101_grid_T.nc sshtime.nc 2 ncks -4 -A -v gdepw_0,ht_wd -d y,17,17,1 -d z,10,10,1 -C mesh_mask.nc sshtime.nc 3 4 python2.7 matpoly2.py sshtime.nc 5 animate wadfr*.png -
branches/2016/dev_merge_2016/NEMOGCM/CONFIG/WAD_TEST_CASES/EXP00/matpoly2.py
r7412 r7514 21 21 fw = Dataset(pathin) 22 22 ssh = fw.variables['sossheig'][:,:,:] 23 #bat = fw.variables['gdepw_0'][0,0,:,:]24 bat = fw.variables['bathy'][0,:,:]23 bat = fw.variables['gdepw_0'][0,0,:,:] 24 #bat = fw.variables['ht_wd'][0,:,:] 25 25 #vot = fw.variables['votemper'][:,0,:,:] 26 26 vot = fw.variables['sosaline'][:,:,:] … … 49 49 dy = np.int(t/24.0) 50 50 51 tfac = 1 2.0/3600.051 tfac = 18.0/3600.0 52 52 t24 = np.int(np.mod(t*tfac,24)) 53 53 dy = np.int(t*tfac/24.0) … … 112 112 p.set_array(np.array(colors)) 113 113 114 ax.set_ylim([-10., 4.0]) 115 ax.set_xlim([0., 50.0]) 114 116 ax.add_collection(p) 115 117 ax.plot(xs,ys, '--', color='black', ms=10) -
branches/2016/dev_merge_2016/NEMOGCM/CONFIG/WAD_TEST_CASES/MY_SRC/usrdef_istate.F90
r7467 r7514 14 14 !!---------------------------------------------------------------------- 15 15 USE par_oce ! ocean space and time domain 16 USE dom_oce , ONLY : mi0, mig, mjg, glamt, ht_016 USE dom_oce , ONLY : mi0, mig, mjg, glamt, gphit, ht_0 17 17 USE phycst ! physical constants 18 18 USE wet_dry ! Wetting and drying … … 74 74 !!---------------------------------------------------------------------- 75 75 ! 76 ! Uniform T & S in alltest cases76 ! Uniform T & S in most test cases 77 77 pts(:,:,:,jp_tem) = 10._wp 78 !tsb(:,:,:,jp_tem) = 10._wp79 78 pts(:,:,:,jp_sal) = 35._wp 80 !tsb(:,:,:,jp_sal) = 35._wp81 79 SELECT CASE ( nn_cfg ) 82 80 ! ! ==================== … … 90 88 do ji = 1,jpi 91 89 pssh(ji,:) = ( -5.5_wp + 5.5_wp*glamt(ji,1)/50._wp)*ptmask(ji,:,1) 92 !write(*,*) 'SSH ',jpk,mig(ji),pssh(ji,1),ht_0(ji,1),-5.5_wp+5.5_wp*glamt(ji,1)/50._wp,ptmask(ji,1,1)93 90 end do 94 91 ! ! ==================== … … 101 98 ! 102 99 do ji = 1,jpi 103 pssh(ji,:) = ( -5.5_wp + 3.9_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1)100 pssh(ji,:) = ( -5.5_wp + 6.5_wp*glamt(ji,1)/50._wp)*ptmask(ji,:,1) 104 101 end do 105 102 ! ! ==================== … … 112 109 ! 113 110 do ji = 1,jpi 114 pssh(ji,:) = ( - 7.5_wp + 6.9_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1)111 pssh(ji,:) = ( -5.5_wp + 5.5_wp*glamt(ji,1)/50._wp)*ptmask(ji,:,1) 115 112 end do 116 113 … … 125 122 ! 126 123 DO ji = 1, jpi 127 zi = MAX(1.0- FLOAT((mig(ji)-25)**2)/400.0, 0.0 )124 zi = MAX(1.0-((glamt(ji,1)-25._wp)**2)/400.0, 0.0 ) 128 125 DO jj = 1, jpj 129 zj = MAX(1.0- FLOAT((mjg(jj)-17)**2)/144.0, 0.0 )126 zj = MAX(1.0-((gphit(1,jj)-17._wp)**2)/144.0, 0.0 ) 130 127 pssh(ji,jj) = -8.5_wp + 8.5_wp*zi*zj 131 128 END DO … … 141 138 IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 142 139 ! 143 ! Needed rn_wdmin2 increased to 0.01 for this case?144 140 do ji = 1,jpi 145 pssh(ji,:) = ( -5.5_wp + 9.0_wp*FLOAT(mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1)141 pssh(ji,:) = ( -5.5_wp + 5.5_wp*glamt(ji,1)/50._wp)*ptmask(ji,:,1) 146 142 end do 147 143 … … 166 162 !pssh(ji,:) = ( -4.5_wp + 8.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1) 167 163 !6e 168 pssh(ji,:) = ( -3.5_wp + 7.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1)164 pssh(ji,:) = ( -3.5_wp + 5.5_wp*(50._wp-glamt(ji,1))/50._wp)*ptmask(ji,:,1) 169 165 !6f 170 166 !pssh(ji,:) = ( 0.5_wp + 3.75_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1) … … 173 169 do ji = mi0(jpiglo/2), mi0(jpiglo) 174 170 pts(ji,:,:,jp_sal) = 30._wp 175 !tsb(ji,:,:,jp_sal) = 30._wp176 171 end do 177 172 ! … … 191 186 IF( ht_wd(ji,jj) + pssh(ji,jj) < rn_wdmin1 ) THEN 192 187 pssh(ji,jj) = ptmask(ji,jj,1)*( rn_wdmin1 - ht_wd(ji,jj) ) 193 !IF (jj.eq.4) write(*,*) 'SSH2 ',mig(ji),pssh(ji,4),ht_0(ji,4),-5.5_wp+5.5_wp*glamt(ji,4)/50._wp,ptmask(ji,4,1)194 188 ENDIF 195 189 end do 196 190 end do 197 !sshb = sshn198 !ssha = sshn199 191 ! 200 ! CALL wad_istate ! WAD test configuration : start from pre-defined T-S fields and initial ssh slope201 !202 192 END SUBROUTINE usr_def_istate 203 193 -
branches/2016/dev_merge_2016/NEMOGCM/CONFIG/WAD_TEST_CASES/MY_SRC/usrdef_zgr.F90
r7467 r7514 106 106 zi = MIN((glamt(ji,1) - 10.0)/40.0, 1.0 ) 107 107 zht(ji,:) = MAX(zbathy*zi, 1.5) 108 IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zht(ji,1)109 108 END DO 110 109 zht(mi0(1):mi1(1),:) = -4._wp … … 122 121 DO ji = 1, jpi 123 122 zi = MAX(1.0-((glamt(ji,1)-25.0)**2)/484.0, 0.0 ) 124 zi = MAX(1.0-((glamt(ji,1)-25.0)**2)/484.0, -2.0 ) 125 zht(ji,:) = MAX(zbathy*zi, -20.0) 126 IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zht(ji,1) 127 END DO 128 zht(mi0(1):mi1(2),:) = -4._wp 129 zht(mi0(jpiglo-1):mi1(jpiglo),:) = -4._wp 123 zht(ji,:) = MAX(zbathy*zi, 1.5) 124 END DO 125 zht(mi0(1):mi1(1),:) = -4._wp 126 zht(mi0(jpiglo):mi1(jpiglo),:) = -4._wp 130 127 zht(:,mj0(1):mj1(1)) = -4._wp 131 128 zht(:,mj0(jpjglo):mj1(jpjglo)) = -4._wp … … 142 139 DO jj = 1, jpj 143 140 zj = MAX(1.0-((gphit(1,jj)-17.0)**2)/196.0, -2.0 ) 144 zht(ji,jj) = MAX(zbathy*zi*zj, -2.0) 145 END DO 146 IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zht(ji,1) 141 zht(ji,jj) = MAX(zbathy*zi*zj, 1.5) 142 END DO 147 143 END DO 148 144 zht(mi0(1):mi1(1),:) = -4._wp 149 zht(mi0(jpiglo -1):mi1(jpiglo),:) = -4._wp145 zht(mi0(jpiglo):mi1(jpiglo),:) = -4._wp 150 146 zht(:,mj0(1):mj1(1)) = -4._wp 151 147 zht(:,mj0(jpjglo):mj1(jpjglo)) = -4._wp … … 160 156 DO ji = 1, jpi 161 157 zi = MIN(glamt(ji,1)/45.0, 1.0 ) 162 zht(ji,:) = MAX(zbathy*zi, 0.5)158 zht(ji,:) = MAX(zbathy*zi, 1.5) 163 159 IF( glamt(ji,1) >= 46.0 ) THEN 164 160 zht(ji,:) = 10.0 … … 170 166 ELSE IF( glamt(ji,1) >= 4.0 .AND. glamt(ji,1) < 15.0 ) THEN 171 167 zi = 2.0/11.0 172 zht(ji,:) = MAX(2.5 - zi*(16.0-glamt(ji,1)), 0.5)168 zht(ji,:) = MAX(2.5 - zi*(16.0-glamt(ji,1)), 1.5) 173 169 ELSE IF( glamt(ji,1) >= 1.0 .AND. glamt(ji,1) < 4.0 ) THEN 174 zht(ji,:) = 0.5170 zht(ji,:) = 1.5 175 171 ENDIF 176 IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zht(ji,1)177 172 END DO 178 173 ! ! =========================== … … 192 187 zi = MAX(1.0-((glamt(ji,1)-25.0)**2)/484.0, -2.0 ) 193 188 zj = 0.95*MAX(EXP(-1.0*((glamt(ji,1)-25.0)**2)/32.0) , 0.0 ) 194 zht(ji,:) = MAX(zbathy*(zi-zj), -2.0) 195 IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zht(ji,1) 189 zht(ji,:) = MAX(zbathy*(zi-zj), 1.0) 196 190 END DO 197 191 zht(mi0(1):mi1(1),:) = -4._wp … … 214 208 zhu(ji,:) = 0.5_wp * ( zht(ji,:) + zht(ji+1,:) ) 215 209 END DO 216 IF(lwp) write(numout,*) 'E3a zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))217 210 CALL lbc_lnk( zhu, 'U', 1. ) ! boundary condition: this mask the surrounding grid-points 218 211 ! ! ==>>> set by hand non-zero value on first/last columns & rows 219 IF(lwp) write(numout,*) 'E3b zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))220 212 DO ji = mi0(1), mi1(1) ! first row of global domain only 221 213 zhu(ji,:) = zht(1,:) … … 224 216 zhu(ji,:) = zht(jpi,:) 225 217 END DO 226 IF(lwp) write(numout,*) 'E3c zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))227 218 ! at v-point: averaging zht 228 219 zhv = 0._wp … … 231 222 END DO 232 223 CALL lbc_lnk( zhv, 'V', 1. ) ! boundary condition: this mask the surrounding grid-points 233 IF(lwp) write(numout,*) 'E3d zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))234 224 !avoid the zero depth on T- (U-,V-,F-) points 235 225 DO jj = 1, jpj … … 243 233 END DO 244 234 END DO 245 IF(lwp) write(numout,*) 'E3e zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))246 235 ! 247 236 CALL zgr_z1d( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d ) ! Reference z-coordinate system … … 274 263 DO ji = 1, jpi 275 264 IF( zht(ji,jj) <= -(rn_wdld - rn_wdmin2)) THEN 276 !WRITE(numout,*) 'KBOT ',ji,jj,zht(ji,jj),rn_wdld,rn_wdmin2277 265 k_bot(ji,jj) = 0 278 266 zht(ji,jj) = rn_wdmin1 * REAL( jpkm1 , wp) * 0.5_wp … … 280 268 k_bot(ji,jj) = 2 281 269 zht(ji,jj) = rn_wdmin1 * REAL( jpkm1 , wp) * 0.5_wp 282 !WRITE(numout,*) 'KBOT1 ',ji,jj,zht(ji,jj),rn_wdmin1283 270 ENDIF 284 271 END DO … … 287 274 zhu(ji,:) = 0.5_wp * ( zht(ji,:) + zht(ji+1,:) ) 288 275 END DO 289 IF(lwp) write(numout,*) 'E3f zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))290 276 ! ! ==>>> set by hand non-zero value on 291 277 ! first/last columns & rows 292 IF(lwp) write(numout,*) 'E3f Zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))293 278 CALL lbc_lnk( zhu, 'U', 1. ) ! boundary condition: this mask the surrounding grid-points 294 IF(lwp) write(numout,*) 'E3F zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))295 279 DO ji = mi0(1), mi1(1) ! first row of global domain only 296 280 zhu(ji,:) = zht(ji,:) … … 300 284 END DO 301 285 ! at v-point: averaging zht 302 IF(lwp) write(numout,*) 'E3g zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))303 286 zhv = 0._wp 304 287 DO jj = 1, jpjm1 305 288 zhv(:,jj) = 0.5_wp * ( zht(:,jj) + zht(:,jj+1) ) 306 289 END DO 307 IF(lwp) write(numout,*) 'E3h zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))308 290 CALL lbc_lnk( zhv, 'V', 1. ) 309 IF(lwp) write(numout,*) 'E3H zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:))310 291 DO jj = mj0(1), mj1(1) ! first row of global domain only 311 292 zhv(:,jj) = zht(:,jj) … … 359 340 CALL lbc_lnk( pe3vw, 'V', 1. ) 360 341 CALL lbc_lnk( pe3f , 'F', 1. ) 361 WHERE( pe3t (:,:,:) == 0._wp ) pe3t (:,:,:) = 1._wp 362 WHERE( pe3u (:,:,:) == 0._wp ) pe3u (:,:,:) = 1._wp 363 WHERE( pe3v (:,:,:) == 0._wp ) pe3v (:,:,:) = 1._wp 364 WHERE( pe3f (:,:,:) == 0._wp ) pe3f (:,:,:) = 1._wp 365 WHERE( pe3w (:,:,:) == 0._wp ) pe3w (:,:,:) = 1._wp 366 WHERE( pe3uw(:,:,:) == 0._wp ) pe3uw(:,:,:) = 1._wp 367 WHERE( pe3vw(:,:,:) == 0._wp ) pe3vw(:,:,:) = 1._wp 368 IF(lwp) write(numout,*) 'E3 TU ',SUM(0._wp/pe3t(:,:,:)), MINVAL(pe3t(:,:,:)), rn_wdmin1 ; CALL FLUSH(numout) 369 IF(lwp) write(numout,*) 'E3 tu ',SUM(0._wp/pe3u(:,:,:)), MINVAL(pe3u(:,:,:)) 370 IF(lwp) write(numout,*) 'E3 tv ',SUM(0._wp/pe3v(:,:,:)), MINVAL(pe3v(:,:,:)) 371 IF(lwp) write(numout,*) 'E3 tf ',SUM(0._wp/pe3f(:,:,:)), MINVAL(pe3f(:,:,:)) 372 IF(lwp) write(numout,*) 'E3 zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 342 WHERE( pe3t (:,:,:) == 0._wp ) pe3t (:,:,:) = 1._wp 343 WHERE( pe3u (:,:,:) == 0._wp ) pe3u (:,:,:) = 1._wp 344 WHERE( pe3v (:,:,:) == 0._wp ) pe3v (:,:,:) = 1._wp 345 WHERE( pe3f (:,:,:) == 0._wp ) pe3f (:,:,:) = 1._wp 346 WHERE( pe3w (:,:,:) == 0._wp ) pe3w (:,:,:) = 1._wp 347 WHERE( pe3uw(:,:,:) == 0._wp ) pe3uw(:,:,:) = 1._wp 348 WHERE( pe3vw(:,:,:) == 0._wp ) pe3vw(:,:,:) = 1._wp 373 349 ENDIF 374 350 ! -
branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90
r7421 r7514 41 41 USE zdfddm ! vertical physics: double diffusion 42 42 USE diahth ! thermocline diagnostics 43 USE wet_dry ! wetting and drying 43 44 ! 44 45 USE lbclnk ! ocean lateral boundary conditions (or mpp link) … … 153 154 154 155 CALL iom_put( "ssh" , sshn ) ! sea surface height 156 IF( iom_use("wetdep") ) & ! wet depth 157 CALL iom_put( "wetdep" , ht_wd(:,:) + sshn(:,:) ) 155 158 156 159 CALL iom_put( "toce", tsn(:,:,:,jp_tem) ) ! 3D temperature -
branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90
r7490 r7514 39 39 USE domc1d ! 1D configuration: column location 40 40 USE dyncor_c1d ! 1D configuration: Coriolis term (cor_c1d routine) 41 !41 USE wet_dry ! wetting and drying 42 42 ! 43 43 USE in_out_manager ! I/O manager … … 664 664 ENDIF 665 665 ! 666 IF( ln_wd ) THEN ! wetting and drying domain 667 CALL iom_rstput( 0, 0, inum, 'ht_0' , ht_0 , ktype = jp_r8 ) 668 CALL iom_rstput( 0, 0, inum, 'ht_wd' , ht_wd , ktype = jp_r8 ) 669 ENDIF 670 ! 666 671 ! ! ============================ 667 672 ! ! close the files -
branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
r7430 r7514 24 24 USE sbc_oce ! ocean surface boundary condition 25 25 USE wet_dry ! wetting and drying 26 USE usrdef_istate ! user defined initial state (wad only) 26 27 USE restart ! ocean restart 27 28 ! … … 876 877 ! 877 878 IF( ln_wd .AND. ( cn_cfg == 'wad' ) ) THEN 878 ! 879 CALL wad_istate ! WAD test configuration : start from 879 ! Wetting and drying test case 880 CALL usr_def_istate( gdept_b, tmask, tsb, ub, vb, sshb ) 881 tsn (:,:,:,:) = tsb (:,:,:,:) ! set now values from to before ones 882 sshn (:,:) = sshb(:,:) 883 un (:,:,:) = ub (:,:,:) 884 vn (:,:,:) = vb (:,:,:) 880 885 ! uniform T-S fields and initial ssh slope 881 886 ! needs to be called here and in istate which is called later. … … 893 898 DO ji = 1, jpi 894 899 IF( e3t_0(ji,jj,1) <= 0.5_wp * rn_wdmin1 ) THEN 900 ! potential bug 901 ! Warning this assumes 2 layers only over wetting locations. needs investigating 895 902 e3t_b(ji,jj,:) = 0.5_wp * rn_wdmin1 896 903 e3t_n(ji,jj,:) = 0.5_wp * rn_wdmin1 897 904 e3t_a(ji,jj,:) = 0.5_wp * rn_wdmin1 898 sshb(ji,jj) = rn_wdmin1 - ht_ 0(ji,jj) !!gm I don't understand that !899 sshn(ji,jj) = rn_wdmin1 - ht_ 0(ji,jj)900 ssha(ji,jj) = rn_wdmin1 - ht_ 0(ji,jj)905 sshb(ji,jj) = rn_wdmin1 - ht_wd(ji,jj) !!gm I don't understand that ! 906 sshn(ji,jj) = rn_wdmin1 - ht_wd(ji,jj) 907 ssha(ji,jj) = rn_wdmin1 - ht_wd(ji,jj) 901 908 ENDIF 902 909 ENDDO -
branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domwri.F90
r7421 r7514 18 18 USE dom_oce ! ocean space and time domain 19 19 USE phycst , ONLY : rsmall 20 USE wet_dry, ONLY : ln_wd, ht_wd 20 21 ! 21 22 USE in_out_manager ! I/O manager … … 194 195 IF( ln_sco ) THEN ! s-coordinate stiffness 195 196 CALL dom_stiff( zprt ) 196 CALL iom_rstput( 0, 0, inum, 'stiffness', zprt ) ! ! Max. grid stiffness ratio 197 CALL iom_rstput( 0, 0, inum, 'stiffness', zprt ) ! Max. grid stiffness ratio 198 ENDIF 199 ! 200 IF( ln_wd ) THEN ! wetting and drying domain 201 CALL iom_rstput( 0, 0, inum, 'ht_0' , ht_0 , ktype = jp_r8 ) 202 CALL iom_rstput( 0, 0, inum, 'ht_wd' , ht_wd , ktype = jp_r8 ) 197 203 ENDIF 198 204 ! ! ============================ -
branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/dynhpg.F90
r7437 r7514 432 432 INTEGER :: ji, jj, jk, jii, jjj ! dummy loop indices 433 433 REAL(wp) :: zcoef0, zuap, zvap, znad, ztmp ! temporary scalars 434 LOGICAL :: ll_tmp1, ll_tmp2 , ll_tmp3! local logical variables434 LOGICAL :: ll_tmp1, ll_tmp2 ! local logical variables 435 435 REAL(wp), POINTER, DIMENSION(:,:,:) :: zhpi, zhpj 436 436 REAL(wp), POINTER, DIMENSION(:,:) :: zcpx, zcpy !W/D pressure filter … … 438 438 ! 439 439 CALL wrk_alloc( jpi,jpj,jpk, zhpi, zhpj ) 440 IF( ln_wd) CALL wrk_alloc( jpi,jpj, zcpx, zcpy )440 IF( ln_wd ) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 441 441 ! 442 442 IF( kt == nit000 ) THEN … … 451 451 ENDIF 452 452 ! 453 IF( ln_wd) THEN453 IF( ln_wd ) THEN 454 454 DO jj = 2, jpjm1 455 455 DO ji = 2, jpim1 456 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) 457 ll_tmp2 = MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji+1,jj) + ht_0(ji+1,jj)) > rn_wdmin1 + rn_wdmin2 458 ll_tmp3 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) + & 459 & rn_wdmin1 + rn_wdmin2 460 461 IF(ll_tmp1.AND.ll_tmp2) THEN 456 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji+1,jj) ) > & 457 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) .AND. & 458 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 459 & > rn_wdmin1 + rn_wdmin2 460 ll_tmp2 = MAX( sshn(ji,jj) , sshn(ji+1,jj) ) > & 461 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 462 463 IF(ll_tmp1) THEN 462 464 zcpx(ji,jj) = 1.0_wp 463 ELSE IF(ll_tmp 3) THEN464 ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen ! here465 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_0(ji+1,jj) - sshn(ji,jj) - ht_0(ji,jj)) /&466 & (sshn(ji+1,jj) - sshn(ji,jj)))465 ELSE IF(ll_tmp2) THEN 466 ! no worries about sshn(ji+1,jj) - sshn(ji ,jj) = 0, it won't happen ! here 467 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 468 & / (sshn(ji+1,jj) - sshn(ji ,jj)) ) 467 469 ELSE 468 470 zcpx(ji,jj) = 0._wp 469 471 END IF 470 472 471 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) 472 ll_tmp2 = MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji,jj+1) + ht_0(ji,jj+1)) > rn_wdmin1 + rn_wdmin2 473 ll_tmp3 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) + & 474 & rn_wdmin1 + rn_wdmin2 475 476 IF(ll_tmp1.AND.ll_tmp2) THEN 473 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji,jj+1) ) > & 474 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) .AND. & 475 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji,jj+1) + ht_wd(ji,jj+1) ) & 476 & > rn_wdmin1 + rn_wdmin2 477 ll_tmp2 = MAX( sshn(ji,jj) , sshn(ji,jj+1) ) > & 478 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 479 480 IF(ll_tmp1) THEN 477 481 zcpy(ji,jj) = 1.0_wp 478 ELSE IF(ll_tmp 3) THEN479 ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen ! here480 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_0(ji,jj+1) - sshn(ji,jj) - ht_0(ji,jj)) /&481 & (sshn(ji,jj+1) - sshn(ji,jj)))482 ELSE IF(ll_tmp2) THEN 483 ! no worries about sshn(ji,jj+1) - sshn(ji,jj ) = 0, it won't happen ! here 484 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 485 & / (sshn(ji,jj+1) - sshn(ji,jj )) ) 482 486 ELSE 483 487 zcpy(ji,jj) = 0._wp … … 486 490 END DO 487 491 CALL lbc_lnk( zcpx, 'U', 1._wp ) ; CALL lbc_lnk( zcpy, 'V', 1._wp ) 488 ENDIF 489 492 END IF 490 493 491 494 ! Surface value … … 504 507 505 508 506 IF( ln_wd) THEN509 IF( ln_wd ) THEN 507 510 508 511 zhpi(ji,jj,1) = zhpi(ji,jj,1) * zcpx(ji,jj) … … 535 538 & * ( gde3w_n(ji ,jj+1,jk) - gde3w_n(ji,jj,jk) ) * r1_e2v(ji,jj) 536 539 537 IF( ln_wd) THEN540 IF( ln_wd ) THEN 538 541 zhpi(ji,jj,jk) = zhpi(ji,jj,jk) * zcpx(ji,jj) 539 542 zhpj(ji,jj,jk) = zhpj(ji,jj,jk) * zcpy(ji,jj) … … 550 553 ! 551 554 CALL wrk_dealloc( jpi,jpj,jpk, zhpi, zhpj ) 552 IF( ln_wd) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy )555 IF( ln_wd ) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 553 556 ! 554 557 END SUBROUTINE hpg_sco … … 695 698 CALL wrk_alloc( jpi, jpj, jpk, drhox, drhoy, drhoz, drhou, drhov, drhow ) 696 699 CALL wrk_alloc( jpi, jpj, jpk, rho_i, rho_j, rho_k, zhpi, zhpj ) 697 IF( ln_wd) CALL wrk_alloc( jpi,jpj, zcpx, zcpy )698 ! 699 ! 700 IF( ln_wd) THEN700 IF( ln_wd ) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 701 ! 702 ! 703 IF( ln_wd ) THEN 701 704 DO jj = 2, jpjm1 702 705 DO ji = 2, jpim1 703 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) & 704 & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji+1,jj) + ht_0(ji+1,jj)) & 705 & > rn_wdmin1 + rn_wdmin2 706 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) +& 707 & rn_wdmin1 + rn_wdmin2 706 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji+1,jj) ) > & 707 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) .AND. & 708 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 709 & > rn_wdmin1 + rn_wdmin2 710 ll_tmp2 = MAX( sshn(ji,jj) , sshn(ji+1,jj) ) > & 711 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 708 712 709 713 IF(ll_tmp1) THEN 710 714 zcpx(ji,jj) = 1.0_wp 711 715 ELSE IF(ll_tmp2) THEN 712 ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen ! here713 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_0(ji+1,jj) - sshn(ji,jj) - ht_0(ji,jj)) /&714 & (sshn(ji+1,jj) - sshn(ji,jj)))716 ! no worries about sshn(ji+1,jj) - sshn(ji ,jj) = 0, it won't happen ! here 717 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 718 & / (sshn(ji+1,jj) - sshn(ji ,jj)) ) 715 719 ELSE 716 720 zcpx(ji,jj) = 0._wp 717 721 END IF 718 722 719 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) & 720 & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji,jj+1) + ht_0(ji,jj+1)) & 721 & > rn_wdmin1 + rn_wdmin2 722 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) +& 723 & rn_wdmin1 + rn_wdmin2 723 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji,jj+1) ) > & 724 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) .AND. & 725 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji,jj+1) + ht_wd(ji,jj+1) ) & 726 & > rn_wdmin1 + rn_wdmin2 727 ll_tmp2 = MAX( sshn(ji,jj) , sshn(ji,jj+1) ) > & 728 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 724 729 725 730 IF(ll_tmp1) THEN 726 731 zcpy(ji,jj) = 1.0_wp 727 732 ELSE IF(ll_tmp2) THEN 728 ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen ! here729 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_0(ji,jj+1) - sshn(ji,jj) - ht_0(ji,jj)) /&730 & (sshn(ji,jj+1) - sshn(ji,jj)))733 ! no worries about sshn(ji,jj+1) - sshn(ji,jj ) = 0, it won't happen ! here 734 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 735 & / (sshn(ji,jj+1) - sshn(ji,jj )) ) 731 736 ELSE 732 737 zcpy(ji,jj) = 0._wp … … 735 740 END DO 736 741 CALL lbc_lnk( zcpx, 'U', 1._wp ) ; CALL lbc_lnk( zcpy, 'V', 1._wp ) 737 ENDIF 738 742 END IF 739 743 740 744 IF( kt == nit000 ) THEN … … 907 911 zhpi(ji,jj,1) = ( rho_k(ji+1,jj ,1) - rho_k(ji,jj,1) - rho_i(ji,jj,1) ) * r1_e1u(ji,jj) 908 912 zhpj(ji,jj,1) = ( rho_k(ji ,jj+1,1) - rho_k(ji,jj,1) - rho_j(ji,jj,1) ) * r1_e2v(ji,jj) 909 IF( ln_wd) THEN913 IF( ln_wd ) THEN 910 914 zhpi(ji,jj,1) = zhpi(ji,jj,1) * zcpx(ji,jj) 911 915 zhpj(ji,jj,1) = zhpj(ji,jj,1) * zcpy(ji,jj) … … 930 934 & + ( ( rho_k(ji,jj+1,jk) - rho_k(ji,jj,jk ) ) & 931 935 & -( rho_j(ji,jj ,jk) - rho_j(ji,jj,jk-1) ) ) * r1_e2v(ji,jj) 932 IF( ln_wd) THEN936 IF( ln_wd ) THEN 933 937 zhpi(ji,jj,jk) = zhpi(ji,jj,jk) * zcpx(ji,jj) 934 938 zhpj(ji,jj,jk) = zhpj(ji,jj,jk) * zcpy(ji,jj) … … 944 948 CALL wrk_dealloc( jpi, jpj, jpk, drhox, drhoy, drhoz, drhou, drhov, drhow ) 945 949 CALL wrk_dealloc( jpi, jpj, jpk, rho_i, rho_j, rho_k, zhpi, zhpj ) 946 IF( ln_wd) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy )950 IF( ln_wd ) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 947 951 ! 948 952 END SUBROUTINE hpg_djc … … 981 985 CALL wrk_alloc( jpi,jpj,jpk, zdept, zrhh ) 982 986 CALL wrk_alloc( jpi,jpj, zsshu_n, zsshv_n ) 983 IF( ln_wd) CALL wrk_alloc( jpi,jpj, zcpx, zcpy )987 IF( ln_wd ) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 984 988 ! 985 989 IF( kt == nit000 ) THEN … … 994 998 IF( ln_linssh ) znad = 0._wp 995 999 996 IF( ln_wd) THEN1000 IF( ln_wd ) THEN 997 1001 DO jj = 2, jpjm1 998 1002 DO ji = 2, jpim1 999 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) & 1000 & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji+1,jj) + ht_0(ji+1,jj)) & 1001 & > rn_wdmin1 + rn_wdmin2 1002 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) +& 1003 & rn_wdmin1 + rn_wdmin2 1003 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji+1,jj) ) > & 1004 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) .AND. & 1005 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 1006 & > rn_wdmin1 + rn_wdmin2 1007 ll_tmp2 = MAX( sshn(ji,jj) , sshn(ji+1,jj) ) > & 1008 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 1004 1009 1005 1010 IF(ll_tmp1) THEN 1006 1011 zcpx(ji,jj) = 1.0_wp 1007 1012 ELSE IF(ll_tmp2) THEN 1008 ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen ! here1009 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_0(ji+1,jj) - sshn(ji,jj) - ht_0(ji,jj)) /&1010 & (sshn(ji+1,jj) - sshn(ji,jj)))1013 ! no worries about sshn(ji+1,jj) - sshn(ji ,jj) = 0, it won't happen ! here 1014 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 1015 & / (sshn(ji+1,jj) - sshn(ji ,jj)) ) 1011 1016 ELSE 1012 1017 zcpx(ji,jj) = 0._wp 1013 1018 END IF 1014 1019 1015 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) & 1016 & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji,jj+1) + ht_0(ji,jj+1)) & 1017 & > rn_wdmin1 + rn_wdmin2 1018 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) +& 1019 & rn_wdmin1 + rn_wdmin2 1020 1021 IF(ll_tmp1.OR.ll_tmp2) THEN 1020 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji,jj+1) ) > & 1021 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) .AND. & 1022 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji,jj+1) + ht_wd(ji,jj+1) ) & 1023 & > rn_wdmin1 + rn_wdmin2 1024 ll_tmp2 = MAX( sshn(ji,jj) , sshn(ji,jj+1) ) > & 1025 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 1026 1027 IF(ll_tmp1) THEN 1022 1028 zcpy(ji,jj) = 1.0_wp 1023 1029 ELSE IF(ll_tmp2) THEN 1024 ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen ! here1025 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_0(ji,jj+1) - sshn(ji,jj) - ht_0(ji,jj)) /&1026 & (sshn(ji,jj+1) - sshn(ji,jj)))1030 ! no worries about sshn(ji,jj+1) - sshn(ji,jj ) = 0, it won't happen ! here 1031 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 1032 & / (sshn(ji,jj+1) - sshn(ji,jj )) ) 1027 1033 ELSE 1028 1034 zcpy(ji,jj) = 0._wp … … 1031 1037 END DO 1032 1038 CALL lbc_lnk( zcpx, 'U', 1._wp ) ; CALL lbc_lnk( zcpy, 'V', 1._wp ) 1033 END IF1039 END IF 1034 1040 1035 1041 ! Clean 3-D work arrays … … 1215 1221 zdpdx2 = zcoef0 * r1_e1u(ji,jj) * REAL(jis-jid, wp) * (zpwes + zpwed) 1216 1222 ENDIF 1217 IF( ln_wd) THEN1223 IF( ln_wd ) THEN 1218 1224 zdpdx1 = zdpdx1 * zcpx(ji,jj) 1219 1225 zdpdx2 = zdpdx2 * zcpx(ji,jj) … … 1274 1280 zdpdy2 = zcoef0 * r1_e2v(ji,jj) * REAL(jjs-jjd, wp) * (zpnss + zpnsd ) 1275 1281 ENDIF 1276 IF( ln_wd) THEN1282 IF( ln_wd ) THEN 1277 1283 zdpdy1 = zdpdy1 * zcpy(ji,jj) 1278 1284 zdpdy2 = zdpdy2 * zcpy(ji,jj) … … 1289 1295 CALL wrk_dealloc( jpi,jpj,jpk, zdept, zrhh ) 1290 1296 CALL wrk_dealloc( jpi,jpj, zsshu_n, zsshv_n ) 1291 IF( ln_wd) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy )1297 IF( ln_wd ) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 1292 1298 ! 1293 1299 END SUBROUTINE hpg_prj -
branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90
r7482 r7514 261 261 !!gm 262 262 !! 263 !! IF ( .not. ln_sco ) THEN 264 !! 265 !! !!gm agree the JC comment : this should be done in a much clear way 266 !! 267 !! ! JC: It not clear yet what should be the depth at f-points over land in z-coordinate case 268 !! ! Set it to zero for the time being 269 !! ! IF( rn_hmin < 0._wp ) THEN ; jk = - INT( rn_hmin ) ! from a nb of level 270 !! ! ELSE ; jk = MINLOC( gdepw_0, mask = gdepw_0 > rn_hmin, dim = 1 ) ! from a depth 271 !! ! ENDIF 272 !! ! zhf(:,:) = gdepw_0(:,:,jk+1) 273 !! ELSE 274 !! zhf(:,:) = hbatf(:,:) 275 !! END IF 276 !! 277 !! DO jj = 1, jpjm1 278 !! zhf(:,jj) = zhf(:,jj) * (1._wp- umask(:,jj,1) * umask(:,jj+1,1)) 279 !! END DO 263 IF ( .not. ln_sco ) THEN 264 265 !!gm agree the JC comment : this should be done in a much clear way 266 267 ! JC: It not clear yet what should be the depth at f-points over land in z-coordinate case 268 ! Set it to zero for the time being 269 ! IF( rn_hmin < 0._wp ) THEN ; jk = - INT( rn_hmin ) ! from a nb of level 270 ! ELSE ; jk = MINLOC( gdepw_0, mask = gdepw_0 > rn_hmin, dim = 1 ) ! from a depth 271 ! ENDIF 272 ! zhf(:,:) = gdepw_0(:,:,jk+1) 273 ELSE 274 !zhf(:,:) = hbatf(:,:) 275 DO jj = 1, jpjm1 276 DO ji = 1, jpim1 277 zhf(ji,jj) = MAX( 0._wp, & 278 & ( ht_0(ji ,jj )*tmask(ji ,jj ,1) + & 279 & ht_0(ji+1,jj )*tmask(ji+1,jj ,1) + & 280 & ht_0(ji ,jj+1)*tmask(ji ,jj+1,1) + & 281 & ht_0(ji+1,jj+1)*tmask(ji+1,jj+1,1) ) / & 282 & ( tmask(ji ,jj ,1) + tmask(ji+1,jj ,1) +& 283 & tmask(ji ,jj+1,1) + tmask(ji+1,jj+1,1) +& 284 & rsmall ) ) 285 END DO 286 END DO 287 END IF 288 289 DO jj = 1, jpjm1 290 zhf(:,jj) = zhf(:,jj) * (1._wp- umask(:,jj,1) * umask(:,jj+1,1)) 291 END DO 280 292 !!gm end 281 293 … … 381 393 IF( .NOT.ln_linssh ) THEN ! Variable volume : remove surface pressure gradient 382 394 IF( ln_wd ) THEN ! Calculating and applying W/D gravity filters 383 DO jj = 2, jpjm1 384 DO ji = 2, jpim1 385 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) & 386 & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji+1,jj) + ht_0(ji+1,jj)) & 387 & > rn_wdmin1 + rn_wdmin2 388 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) & 389 & + rn_wdmin1 + rn_wdmin2 395 DO jj = 2, jpjm1 396 DO ji = 2, jpim1 397 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji+1,jj) ) > & 398 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) .AND. & 399 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 400 & > rn_wdmin1 + rn_wdmin2 401 ll_tmp2 = MAX( sshn(ji,jj) , sshn(ji+1,jj) ) > & 402 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 403 390 404 IF(ll_tmp1) THEN 391 zcpx(ji,jj) 392 ELSE IF(ll_tmp2) THEN393 ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happenhere394 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_0(ji+1,jj) - sshn(ji,jj) - ht_0(ji,jj)) &395 & /(sshn(ji+1,jj) - sshn(ji,jj)))405 zcpx(ji,jj) = 1.0_wp 406 ELSE IF(ll_tmp2) THEN 407 ! no worries about sshn(ji+1,jj) - sshn(ji ,jj) = 0, it won't happen ! here 408 zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 409 & / (sshn(ji+1,jj) - sshn(ji ,jj)) ) 396 410 ELSE 397 zcpx(ji,jj) 411 zcpx(ji,jj) = 0._wp 398 412 END IF 399 400 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) & 401 & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji,jj+1) + ht_0(ji,jj+1)) & 402 & > rn_wdmin1 + rn_wdmin2 403 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) & 404 & + rn_wdmin1 + rn_wdmin2 413 414 ll_tmp1 = MIN( sshn(ji,jj) , sshn(ji,jj+1) ) > & 415 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) .AND. & 416 & MAX( sshn(ji,jj) + ht_wd(ji,jj), sshn(ji,jj+1) + ht_wd(ji,jj+1) ) & 417 & > rn_wdmin1 + rn_wdmin2 418 ll_tmp2 = MAX( sshn(ji,jj) , sshn(ji,jj+1) ) > & 419 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 420 405 421 IF(ll_tmp1) THEN 406 zcpy(ji,jj)= 1.0_wp407 ELSE IF(ll_tmp2) THEN408 ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happenhere409 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_0(ji,jj+1) - sshn(ji,jj) - ht_0(ji,jj)) &410 & /(sshn(ji,jj+1) - sshn(ji,jj)))422 zcpy(ji,jj) = 1.0_wp 423 ELSE IF(ll_tmp2) THEN 424 ! no worries about sshn(ji,jj+1) - sshn(ji,jj ) = 0, it won't happen ! here 425 zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 426 & / (sshn(ji,jj+1) - sshn(ji,jj )) ) 411 427 ELSE 412 zcpy(ji,jj) = 0._wp 413 ENDIF 414 415 END DO 428 zcpy(ji,jj) = 0._wp 429 END IF 430 END DO 416 431 END DO 417 418 CALL lbc_lnk( zcpx, 'U', 1._wp ) ; CALL lbc_lnk( zcpy, 'V', 1._wp ) 419 432 420 433 DO jj = 2, jpjm1 421 434 DO ji = 2, jpim1 422 zu_trd(ji,jj) = (zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) &423 & * r1_e1u(ji,jj) ) * zcpx(ji,jj)424 zv_trd(ji,jj) = (zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) &425 & * r1_e2v(ji,jj) ) * zcpy(ji,jj)435 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj ) - sshn(ji ,jj ) ) & 436 & * r1_e1u(ji,jj) * zcpx(ji,jj) 437 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji ,jj+1) - sshn(ji ,jj ) ) & 438 & * r1_e2v(ji,jj) * zcpy(ji,jj) 426 439 END DO 427 440 END DO … … 570 583 ENDIF 571 584 572 IF( ln_wd ) THEN !preserve the positivity of water depth573 !ssh[b,n,a] should have already been processed for this574 sshbb_e(:,:) = MAX(sshbb_e(:,:), rn_wdmin1 - ht_0(:,:))575 sshb_e(:,:) = MAX(sshb_e(:,:) , rn_wdmin1 - ht_0(:,:))576 ENDIF577 585 ! 578 586 IF (ln_bt_fw) THEN ! FORWARD integration: start from NOW fields … … 649 657 zhup2_e (:,:) = hu_0(:,:) + zwx(:,:) ! Ocean depth at U- and V-points 650 658 zhvp2_e (:,:) = hv_0(:,:) + zwy(:,:) 651 IF( ln_wd ) THEN652 zhup2_e(:,:) = MAX(zhup2_e (:,:), rn_wdmin1)653 zhvp2_e(:,:) = MAX(zhvp2_e (:,:), rn_wdmin1)654 END IF655 659 ELSE 656 660 zhup2_e (:,:) = hu_n(:,:) … … 705 709 END DO 706 710 ssha_e(:,:) = ( sshn_e(:,:) - rdtbt * ( zssh_frc(:,:) + zhdiv(:,:) ) ) * ssmask(:,:) 707 IF( ln_wd ) ssha_e(:,:) = MAX(ssha_e(:,:), rn_wdmin1 - ht_0(:,:))711 708 712 CALL lbc_lnk( ssha_e, 'T', 1._wp ) 709 713 … … 752 756 IF( ln_wd ) THEN ! Calculating and applying W/D gravity filters 753 757 DO jj = 2, jpjm1 754 DO ji = 2, jpim1 755 ll_tmp1 = MIN( zsshp2_e(ji,jj), zsshp2_e(ji+1,jj) ) > MAX( -ht_0(ji,jj), -ht_0(ji+1,jj) ) & 756 & .AND. MAX( zsshp2_e(ji,jj) + ht_0(ji,jj), zsshp2_e(ji+1,jj) + ht_0(ji+1,jj) ) & 757 & > rn_wdmin1 + rn_wdmin2 758 ll_tmp2 = MAX( zsshp2_e(ji,jj), zsshp2_e(ji+1,jj) ) > MAX( -ht_0(ji,jj), -ht_0(ji+1,jj) ) & 759 & + rn_wdmin1 + rn_wdmin2 760 IF(ll_tmp1) THEN 761 zcpx(ji,jj) = 1._wp 762 ELSE IF(ll_tmp2) THEN 763 ! no worries about zsshp2_e(ji+1,jj)-zsshp2_e(ji,jj) = 0, it won't happen here 764 zcpx(ji,jj) = ABS( (zsshp2_e(ji+1,jj) + ht_0(ji+1,jj) - zsshp2_e(ji,jj) - ht_0(ji,jj)) & 765 & / (zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj)) ) 766 ELSE 767 zcpx(ji,jj) = 0._wp 768 END IF 769 770 ll_tmp1 = MIN( zsshp2_e(ji,jj), zsshp2_e(ji,jj+1) ) > MAX( -ht_0(ji,jj), -ht_0(ji,jj+1) ) & 771 & .AND. MAX( zsshp2_e(ji,jj) + ht_0(ji,jj), zsshp2_e(ji,jj+1) + ht_0(ji,jj+1) ) & 772 & > rn_wdmin1 + rn_wdmin2 773 ll_tmp2 = MAX( zsshp2_e(ji,jj), zsshp2_e(ji,jj+1) ) > MAX( -ht_0(ji,jj), -ht_0(ji,jj+1) ) & 774 & + rn_wdmin1 + rn_wdmin2 775 IF(ll_tmp1) THEN 776 zcpy(ji,jj) = 1._wp 777 ELSE IF(ll_tmp2) THEN 778 ! no worries about zsshp2_e(ji,jj+1)-zsshp2_e(ji,jj) = 0, it won't happen here 779 zcpy(ji,jj) = ABS( (zsshp2_e(ji,jj+1) + ht_0(ji,jj+1) - zsshp2_e(ji,jj) - ht_0(ji,jj)) & 780 & / (zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj)) ) 781 ELSE 782 zcpy(ji,jj) = 0._wp 783 END IF 758 DO ji = 2, jpim1 759 ll_tmp1 = MIN( zsshp2_e(ji,jj) , zsshp2_e(ji+1,jj) ) > & 760 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) .AND. & 761 & MAX( zsshp2_e(ji,jj) + ht_wd(ji,jj), zsshp2_e(ji+1,jj) + ht_wd(ji+1,jj) ) & 762 & > rn_wdmin1 + rn_wdmin2 763 ll_tmp2 = MAX( zsshp2_e(ji,jj) , zsshp2_e(ji+1,jj) ) > & 764 & MAX( -ht_wd(ji,jj) , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 765 766 IF(ll_tmp1) THEN 767 zcpx(ji,jj) = 1.0_wp 768 ELSE IF(ll_tmp2) THEN 769 ! no worries about zsshp2_e(ji+1,jj) - zsshp2_e(ji ,jj) = 0, it won't happen ! here 770 zcpx(ji,jj) = ABS( (zsshp2_e(ji+1,jj) + ht_wd(ji+1,jj) - zsshp2_e(ji,jj) - ht_wd(ji,jj)) & 771 & / (zsshp2_e(ji+1,jj) - zsshp2_e(ji ,jj)) ) 772 ELSE 773 zcpx(ji,jj) = 0._wp 774 END IF 775 776 ll_tmp1 = MIN( zsshp2_e(ji,jj) , zsshp2_e(ji,jj+1) ) > & 777 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) .AND. & 778 & MAX( zsshp2_e(ji,jj) + ht_wd(ji,jj), zsshp2_e(ji,jj+1) + ht_wd(ji,jj+1) ) & 779 & > rn_wdmin1 + rn_wdmin2 780 ll_tmp2 = MAX( zsshp2_e(ji,jj) , zsshp2_e(ji,jj+1) ) > & 781 & MAX( -ht_wd(ji,jj) , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 782 783 IF(ll_tmp1) THEN 784 zcpy(ji,jj) = 1.0_wp 785 ELSE IF(ll_tmp2) THEN 786 ! no worries about zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj ) = 0, it won't happen ! here 787 zcpy(ji,jj) = ABS( (zsshp2_e(ji,jj+1) + ht_wd(ji,jj+1) - zsshp2_e(ji,jj) - ht_wd(ji,jj)) & 788 & / (zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj )) ) 789 ELSE 790 zcpy(ji,jj) = 0._wp 791 END IF 784 792 END DO 785 END DO 786 CALL lbc_lnk( zcpx, 'U', 1._wp ) ; CALL lbc_lnk( zcpy, 'V', 1._wp ) 787 ENDIF 793 END DO 794 END IF 788 795 ! 789 796 ! Compute associated depths at U and V points: … … 802 809 END DO 803 810 END DO 804 805 IF( ln_wd ) THEN806 zhust_e(:,:) = MAX(zhust_e (:,:), rn_wdmin1 )807 zhvst_e(:,:) = MAX(zhvst_e (:,:), rn_wdmin1 )808 END IF809 811 810 812 ENDIF … … 1243 1245 zxr2 = r1_e1t(ji,jj) * r1_e1t(ji,jj) 1244 1246 zyr2 = r1_e2t(ji,jj) * r1_e2t(ji,jj) 1245 zcu(ji,jj) = SQRT( grav * ht_0(ji,jj) * (zxr2 + zyr2) )1247 zcu(ji,jj) = SQRT( grav * MAX(ht_0(ji,jj),0._wp) * (zxr2 + zyr2) ) 1246 1248 END DO 1247 1249 END DO -
branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/wet_dry.F90
r7436 r7514 32 32 33 33 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: wdmask !: u- and v- limiter 34 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: ht_wd !: wetting and drying t-pt depths 35 ! (can include negative depths) 34 36 35 37 LOGICAL, PUBLIC :: ln_wd !: Wetting/drying activation switch (T:on,F:off) … … 85 87 ! 86 88 IF(ln_wd) THEN 87 ALLOCATE( wdmask(jpi,jpj), STAT=ierr )89 ALLOCATE( wdmask(jpi,jpj), ht_wd(jpi,jpj), STAT=ierr ) 88 90 IF( ierr /= 0 ) CALL ctl_stop('STOP', 'wad_init : Array allocation error') 89 91 ENDIF … … 123 125 IF(ln_wd) THEN 124 126 125 CALL wrk_alloc( jpi, jpj,zflxp, zflxn, zflxu, zflxv, zflxu1, zflxv1 )126 CALL wrk_alloc( jpi, jpj,zwdlmtu, zwdlmtv)127 CALL wrk_alloc( jpi, jpj, zflxp, zflxn, zflxu, zflxv, zflxu1, zflxv1 ) 128 CALL wrk_alloc( jpi, jpj, zwdlmtu, zwdlmtv) 127 129 ! 128 130 … … 159 161 DO ji = 2, jpi 160 162 161 IF( tmask(ji, jj,1) == 0._wp )CYCLE ! we don't care about land cells162 IF( ht_ 0 (ji,jj) > zdepwd ) CYCLE ! and cells which will unlikely go dried out163 IF( tmask(ji, jj, 1) < 0.5_wp ) CYCLE ! we don't care about land cells 164 IF( ht_wd(ji,jj) > zdepwd ) CYCLE ! and cells which are unlikely to dry 163 165 164 166 zflxp(ji,jj) = max(zflxu(ji,jj), 0._wp) - min(zflxu(ji-1,jj), 0._wp) + & … … 167 169 & min(zflxv(ji,jj), 0._wp) - max(zflxv(ji, jj-1), 0._wp) 168 170 169 zdep2 = ht_ 0(ji,jj) + sshb1(ji,jj) - rn_wdmin1170 IF(zdep2 <0._wp) THEN !add more safty, but not necessary171 zdep2 = ht_wd(ji,jj) + sshb1(ji,jj) - rn_wdmin1 172 IF(zdep2 .le. 0._wp) THEN !add more safty, but not necessary 171 173 !zdep2 = 0._wp 172 sshb1(ji,jj) = rn_wdmin1 - ht_0(ji,jj) 174 sshb1(ji,jj) = rn_wdmin1 - ht_wd(ji,jj) 175 wdmask(ji,jj) = 0._wp 173 176 END IF 174 177 ENDDO … … 186 189 DO ji = 2, jpi 187 190 188 wdmask(ji,jj) = 0 189 IF( tmask(ji,jj,1) < 0.5_wp) CYCLE 190 IF( ht_0(ji,jj) > zdepwd) CYCLE 191 IF( tmask(ji, jj, 1) < 0.5_wp ) CYCLE 192 IF( ht_wd(ji,jj) > zdepwd ) CYCLE 191 193 192 194 ztmp = e1e2t(ji,jj) … … 198 200 199 201 zdep1 = (zzflxp + zzflxn) * z2dt / ztmp 200 zdep2 = ht_ 0(ji,jj) + sshb1(ji,jj) - rn_wdmin1 - z2dt * sshemp(ji,jj) ! this one can be moved out of the loop202 zdep2 = ht_wd(ji,jj) + sshb1(ji,jj) - rn_wdmin1 - z2dt * sshemp(ji,jj) 201 203 202 204 IF(zdep1 > zdep2) THEN 203 205 zflag = 1 204 206 wdmask(ji, jj) = 0 205 zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 207 !zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 208 zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt ) 206 209 zcoef = max(zcoef, 0._wp) 207 210 IF(jk1 > nn_wdit) zcoef = 0._wp … … 307 310 DO ji = 2, jpi 308 311 309 IF( tmask(ji,jj,1) < 0.5_wp) CYCLE ! we don't care about land cells310 IF( ht_0 (ji,jj) > zdepwd) CYCLE ! and cells which will unlikely go dried out312 IF( tmask(ji, jj, 1 ) < 0.5_wp) CYCLE ! we don't care about land cells 313 IF( ht_wd(ji,jj) > zdepwd ) CYCLE ! and cells which are unlikely to dry 311 314 312 315 zflxp(ji,jj) = max(zflxu(ji,jj), 0._wp) - min(zflxu(ji-1,jj), 0._wp) + & … … 315 318 & min(zflxv(ji,jj), 0._wp) - max(zflxv(ji, jj-1), 0._wp) 316 319 317 zdep2 = ht_ 0(ji,jj) + sshn_e(ji,jj) - rn_wdmin1320 zdep2 = ht_wd(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 318 321 ENDDO 319 322 END DO … … 330 333 DO ji = 2, jpi 331 334 332 IF( tmask(ji,jj,1) < 0.5_wp) CYCLE333 IF( ht_0 (ji,jj) > zdepwd)CYCLE335 IF( tmask(ji, jj, 1 ) < 0.5_wp) CYCLE 336 IF( ht_wd(ji,jj) > zdepwd ) CYCLE 334 337 335 338 ztmp = e1e2t(ji,jj) … … 341 344 342 345 zdep1 = (zzflxp + zzflxn) * z2dt / ztmp 343 zdep2 = ht_ 0(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 ! this one can be moved out of the loop346 zdep2 = ht_wd(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 344 347 zdep2 = zdep2 - z2dt * zssh_frc(ji,jj) 345 348 346 349 IF(zdep1 > zdep2) THEN 347 350 zflag = 1 348 ! wdmask(ji, jj) = 1349 zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( z flxp(ji,jj)* z2dt )351 !zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 352 zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt ) 350 353 zcoef = max(zcoef, 0._wp) 351 354 IF(jk1 > nn_wdit) zcoef = 0._wp … … 487 490 do ji = 1,jpi 488 491 !6a 489 sshn(ji,:) = ( -5.5_wp + 9.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*tmask(ji,:,1)492 !sshn(ji,:) = ( -5.5_wp + 9.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*tmask(ji,:,1) 490 493 !Some variations in initial slope that have been tested 491 494 !6b … … 495 498 !6d 496 499 !sshn(ji,:) = ( -4.5_wp + 8.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*tmask(ji,:,1) 500 !6e 501 sshn(ji,:) = ( -3.5_wp + 7.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*tmask(ji,:,1) 502 !6f 503 !sshn(ji,:) = ( 0.5_wp + 3.75_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*tmask(ji,:,1) 497 504 end do 498 505 ! 506 do ji = mi0(jpiglo/2), mi0(jpiglo) 507 tsn(ji,:,:,jp_sal) = 30._wp 508 tsb(ji,:,:,jp_sal) = 30._wp 509 end do 499 510 ! 500 511 ! ! =========================== … … 511 522 do jj = 1,jpj 512 523 do ji = 1,jpi 513 IF( ht_ 0(ji,jj) + sshn(ji,jj) < rn_wdmin1 ) THEN514 sshn(ji,jj) = tmask(ji,jj,1)*( rn_wdmin1 - ht_ 0(ji,jj) )524 IF( ht_wd(ji,jj) + sshn(ji,jj) < rn_wdmin1 ) THEN 525 sshn(ji,jj) = tmask(ji,jj,1)*( rn_wdmin1 - ht_wd(ji,jj) ) 515 526 ENDIF 516 527 end do
Note: See TracChangeset
for help on using the changeset viewer.