- Timestamp:
- 2020-06-03T16:26:23+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
utils/tools_dev_r12970_AGRIF_CMEMS/DOMAINcfg/src/domzgr.F90
r12414 r13024 17 17 !! 3.4 ! 2012-08 (J. Siddorn) added Siddorn and Furner stretching function 18 18 !! 3.4 ! 2012-12 (R. Bourdalle-Badie and G. Reffray) modify C1D case 19 !! 3.6 ! 2014-11 (P. Mathiot and C. Harris) add ice shelf capability 19 !! 3.6 ! 2014-11 (P. Mathiot and C. Harris) add ice shelf capabilitye 20 20 !! 3.? ! 2015-11 (H. Liu) Modifications for Wetting/Drying 21 21 !!---------------------------------------------------------------------- … … 36 36 !!--------------------------------------------------------------------- 37 37 USE dom_oce ! ocean domain 38 USE depth_e3 ! depth <=> e3 38 39 ! 39 40 USE in_out_manager ! I/O manager … … 92 93 CONTAINS 93 94 94 SUBROUTINE dom_zgr 95 SUBROUTINE dom_zgr( k_top, k_bot ) 95 96 !!---------------------------------------------------------------------- 96 97 !! *** ROUTINE dom_zgr *** … … 109 110 !! ** Action : define gdep., e3., mbathy and bathy 110 111 !!---------------------------------------------------------------------- 112 INTEGER, DIMENSION(:,:), INTENT(out) :: k_top, k_bot ! ocean first and last level indices 113 ! 111 114 INTEGER :: ioptio, ibat ! local integer 112 115 INTEGER :: ios 113 116 ! 117 INTEGER :: jk 118 REAL(wp) :: zrefdep ! depth of the reference level (~10m) 119 120 114 121 NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco, ln_isfcav, ln_linssh 115 122 !!---------------------------------------------------------------------- … … 124 131 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzgr in configuration namelist', lwp ) 125 132 IF(lwm) WRITE ( numond, namzgr ) 133 134 IF(ln_read_cfg) THEN 135 IF(lwp) WRITE(numout,*) 136 IF(lwp) WRITE(numout,*) ' ==>>> Read vertical mesh in ', TRIM( cn_domcfg ), ' file' 137 ! 138 CALL zgr_read ( ln_zco , ln_zps , ln_sco, ln_isfcav, & 139 & gdept_1d, gdepw_1d, e3t_1d, e3w_1d , & ! 1D gridpoints depth 140 & gdept_0 , gdepw_0 , & ! gridpoints depth 141 & e3t_0 , e3u_0 , e3v_0 , e3f_0 , & ! vertical scale factors 142 & e3w_0 , e3uw_0 , e3vw_0 , & ! vertical scale factors 143 & k_top , k_bot ) ! 1st & last ocean level 144 ! 145 !!gm to be remove when removing the OLD definition of e3 scale factors so that gde3w disappears 146 ! ! Compute gde3w_0 (vertical sum of e3w) 147 ! gde3w_0(:,:,1) = 0.5_wp * e3w_0(:,:,1) 148 ! DO jk = 2, jpk 149 ! gde3w_0(:,:,jk) = gde3w_0(:,:,jk-1) + e3w_0(:,:,jk) 150 ! END DO 151 ! 152 153 ! ! top/bottom ocean level indices for t-, u- and v-points (f-point also for top) 154 CALL zgr_top_bot( k_top, k_bot ) ! with a minimum value set to 1 155 156 ! ! deepest/shallowest W level Above/Below ~10m 157 !!gm BUG in s-coordinate this does not work! 158 zrefdep = 10._wp - 0.1_wp * MINVAL( e3w_1d ) ! ref. depth with tolerance (10% of minimum layer thickness) 159 nlb10 = MINLOC( gdepw_1d, mask = gdepw_1d > zrefdep, dim = 1 ) ! shallowest W level Below ~10m 160 nla10 = nlb10 - 1 ! deepest W level Above ~10m 161 !!gm end bug 162 163 ENDIF 126 164 127 165 IF(lwp) THEN ! Control print … … 150 188 IF ( ln_sco .AND. ln_isfcav ) ioptio = ioptio + 1 151 189 IF( ioptio > 0 ) CALL ctl_stop( ' Cavity not tested/compatible with full step (zco) and sigma (ln_sco) ' ) 190 191 IF(.NOT.ln_read_cfg) THEN 152 192 ! 153 193 ! Build the vertical coordinate system … … 163 203 ! ----------------------------------- 164 204 CALL zgr_bot_level ! deepest ocean level for t-, u- and v-points 205 k_bot = mbkt 165 206 CALL zgr_top_level ! shallowest ocean level for T-, U-, V- points 207 k_top = mikt 208 209 ENDIF 166 210 ! 167 211 IF( nprint == 1 .AND. lwp ) THEN 168 WRITE(numout,*) ' MIN val mbathy ', MINVAL( mbathy(:,:) ), ' MAX ', MAXVAL( mbathy(:,:) ) 212 WRITE(numout,*) ' MIN val k_top ', MINVAL( k_top(:,:) ), ' MAX ', MAXVAL( k_top(:,:) ) 213 WRITE(numout,*) ' MIN val k_bot ', MINVAL( k_bot(:,:) ), ' MAX ', MAXVAL( k_bot(:,:) ) 169 214 WRITE(numout,*) ' MIN val depth t ', MINVAL( gdept_0(:,:,:) ), & 170 215 & ' w ', MINVAL( gdepw_0(:,:,:) ) … … 181 226 & ' w ', MAXVAL( e3w_0(:,:,:) ) 182 227 ENDIF 228 183 229 ! 184 230 END SUBROUTINE dom_zgr 185 231 232 SUBROUTINE zgr_read( ld_zco , ld_zps , ld_sco , ld_isfcav, & ! type of vertical coordinate 233 & pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d , & ! 1D reference vertical coordinate 234 & pdept , pdepw , & ! 3D t & w-points depth 235 & pe3t , pe3u , pe3v , pe3f , & ! vertical scale factors 236 & pe3w , pe3uw , pe3vw , & ! - - - 237 & k_top , k_bot ) ! top & bottom ocean level 238 !!--------------------------------------------------------------------- 239 !! *** ROUTINE zgr_read *** 240 !! 241 !! ** Purpose : Read the vertical information in the domain configuration file 242 !! 243 !!---------------------------------------------------------------------- 244 LOGICAL , INTENT(out) :: ld_zco, ld_zps, ld_sco ! vertical coordinate flags 245 LOGICAL , INTENT(out) :: ld_isfcav ! under iceshelf cavity flag 246 REAL(wp), DIMENSION(:) , INTENT(out) :: pdept_1d, pdepw_1d ! 1D grid-point depth [m] 247 REAL(wp), DIMENSION(:) , INTENT(out) :: pe3t_1d , pe3w_1d ! 1D vertical scale factors [m] 248 REAL(wp), DIMENSION(:,:,:), INTENT(out) :: pdept, pdepw ! grid-point depth [m] 249 REAL(wp), DIMENSION(:,:,:), INTENT(out) :: pe3t , pe3u , pe3v , pe3f ! vertical scale factors [m] 250 REAL(wp), DIMENSION(:,:,:), INTENT(out) :: pe3w , pe3uw, pe3vw ! - - - 251 INTEGER , DIMENSION(:,:) , INTENT(out) :: k_top , k_bot ! first & last ocean level 252 ! 253 INTEGER :: jk ! dummy loop index 254 INTEGER :: inum ! local logical unit 255 REAL(WP) :: z_zco, z_zps, z_sco, z_cav 256 REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace 257 !!---------------------------------------------------------------------- 258 ! 259 IF(lwp) THEN 260 WRITE(numout,*) 261 WRITE(numout,*) ' zgr_read : read the vertical coordinates in ', TRIM( cn_domcfg ), ' file' 262 WRITE(numout,*) ' ~~~~~~~~' 263 ENDIF 264 ! 265 CALL iom_open( cn_domcfg, inum ) 266 ! 267 ! !* type of vertical coordinate 268 CALL iom_get( inum, 'ln_zco' , z_zco ) 269 CALL iom_get( inum, 'ln_zps' , z_zps ) 270 CALL iom_get( inum, 'ln_sco' , z_sco ) 271 IF( z_zco == 0._wp ) THEN ; ld_zco = .false. ; ELSE ; ld_zco = .true. ; ENDIF 272 IF( z_zps == 0._wp ) THEN ; ld_zps = .false. ; ELSE ; ld_zps = .true. ; ENDIF 273 IF( z_sco == 0._wp ) THEN ; ld_sco = .false. ; ELSE ; ld_sco = .true. ; ENDIF 274 ! 275 ! !* ocean cavities under iceshelves 276 CALL iom_get( inum, 'ln_isfcav', z_cav ) 277 IF( z_cav == 0._wp ) THEN ; ld_isfcav = .false. ; ELSE ; ld_isfcav = .true. ; ENDIF 278 ! 279 ! !* vertical scale factors 280 CALL iom_get( inum, jpdom_unknown, 'e3t_1d' , pe3t_1d ) ! 1D reference coordinate 281 CALL iom_get( inum, jpdom_unknown, 'e3w_1d' , pe3w_1d ) 282 ! 283 CALL iom_get( inum, jpdom_data, 'e3t_0' , pe3t , lrowattr=ln_use_jattr ) ! 3D coordinate 284 CALL iom_get( inum, jpdom_data, 'e3u_0' , pe3u , lrowattr=ln_use_jattr ) 285 CALL iom_get( inum, jpdom_data, 'e3v_0' , pe3v , lrowattr=ln_use_jattr ) 286 CALL iom_get( inum, jpdom_data, 'e3f_0' , pe3f , lrowattr=ln_use_jattr ) 287 CALL iom_get( inum, jpdom_data, 'e3w_0' , pe3w , lrowattr=ln_use_jattr ) 288 CALL iom_get( inum, jpdom_data, 'e3uw_0' , pe3uw , lrowattr=ln_use_jattr ) 289 CALL iom_get( inum, jpdom_data, 'e3vw_0' , pe3vw , lrowattr=ln_use_jattr ) 290 ! 291 ! !* depths 292 ! !- old depth definition (obsolescent feature) 293 IF( iom_varid( inum, 'gdept_1d', ldstop = .FALSE. ) > 0 .AND. & 294 & iom_varid( inum, 'gdepw_1d', ldstop = .FALSE. ) > 0 .AND. & 295 & iom_varid( inum, 'gdept_0' , ldstop = .FALSE. ) > 0 .AND. & 296 & iom_varid( inum, 'gdepw_0' , ldstop = .FALSE. ) > 0 ) THEN 297 CALL ctl_warn( 'zgr_read : old definition of depths and scale factors used ', & 298 & ' depths at t- and w-points read in the domain configuration file') 299 CALL iom_get( inum, jpdom_unknown, 'gdept_1d', pdept_1d ) 300 CALL iom_get( inum, jpdom_unknown, 'gdepw_1d', pdepw_1d ) 301 CALL iom_get( inum, jpdom_data , 'gdept_0' , pdept , lrowattr=ln_use_jattr ) 302 CALL iom_get( inum, jpdom_data , 'gdepw_0' , pdepw , lrowattr=ln_use_jattr ) 303 ! 304 ELSE !- depths computed from e3. scale factors 305 CALL e3_to_depth( pe3t_1d, pe3w_1d, pdept_1d, pdepw_1d ) ! 1D reference depth 306 CALL e3_to_depth( pe3t , pe3w , pdept , pdepw ) ! 3D depths 307 IF(lwp) THEN 308 WRITE(numout,*) 309 WRITE(numout,*) ' Reference 1D z-coordinate depth and scale factors:' 310 WRITE(numout, "(9x,' level gdept_1d gdepw_1d e3t_1d e3w_1d ')" ) 311 WRITE(numout, "(10x, i4, 4f9.2)" ) ( jk, pdept_1d(jk), pdepw_1d(jk), pe3t_1d(jk), pe3w_1d(jk), jk = 1, jpk ) 312 ENDIF 313 ENDIF 314 ! 315 ! !* ocean top and bottom level 316 CALL iom_get( inum, jpdom_data, 'top_level' , z2d , lrowattr=ln_use_jattr ) ! 1st wet T-points (ISF) 317 k_top(:,:) = NINT( z2d(:,:) ) 318 CALL iom_get( inum, jpdom_data, 'bottom_level' , z2d , lrowattr=ln_use_jattr ) ! last wet T-points 319 k_bot(:,:) = NINT( z2d(:,:) ) 320 ! 321 ! reference depth for negative bathy (wetting and drying only) 322 ! IF( ll_wd ) CALL iom_get( inum, 'rn_wd_ref_depth' , ssh_ref ) 323 ! 324 CALL iom_close( inum ) 325 ! 326 END SUBROUTINE zgr_read 327 328 SUBROUTINE zgr_top_bot( k_top, k_bot ) 329 !!---------------------------------------------------------------------- 330 !! *** ROUTINE zgr_top_bot *** 331 !! 332 !! ** Purpose : defines the vertical index of ocean bottom (mbk. arrays) 333 !! 334 !! ** Method : computes from k_top and k_bot with a minimum value of 1 over land 335 !! 336 !! ** Action : mikt, miku, mikv : vertical indices of the shallowest 337 !! ocean level at t-, u- & v-points 338 !! (min value = 1) 339 !! ** Action : mbkt, mbku, mbkv : vertical indices of the deeptest 340 !! ocean level at t-, u- & v-points 341 !! (min value = 1 over land) 342 !!---------------------------------------------------------------------- 343 INTEGER , DIMENSION(:,:), INTENT(in) :: k_top, k_bot ! top & bottom ocean level indices 344 ! 345 INTEGER :: ji, jj ! dummy loop indices 346 REAL(wp), DIMENSION(jpi,jpj) :: zk ! workspace 347 !!---------------------------------------------------------------------- 348 ! 349 IF(lwp) WRITE(numout,*) 350 IF(lwp) WRITE(numout,*) ' zgr_top_bot : ocean top and bottom k-index of T-, U-, V- and W-levels ' 351 IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~' 352 ! 353 mikt(:,:) = MAX( k_top(:,:) , 1 ) ! top ocean k-index of T-level (=1 over land) 354 ! 355 mbkt(:,:) = MAX( k_bot(:,:) , 1 ) ! bottom ocean k-index of T-level (=1 over land) 356 357 ! ! N.B. top k-index of W-level = mikt 358 ! ! bottom k-index of W-level = mbkt+1 359 DO jj = 1, jpjm1 360 DO ji = 1, jpim1 361 miku(ji,jj) = MAX( mikt(ji+1,jj ) , mikt(ji,jj) ) 362 mikv(ji,jj) = MAX( mikt(ji ,jj+1) , mikt(ji,jj) ) 363 mikf(ji,jj) = MAX( mikt(ji ,jj+1) , mikt(ji,jj), mikt(ji+1,jj ), mikt(ji+1,jj+1) ) 364 ! 365 mbku(ji,jj) = MIN( mbkt(ji+1,jj ) , mbkt(ji,jj) ) 366 mbkv(ji,jj) = MIN( mbkt(ji ,jj+1) , mbkt(ji,jj) ) 367 END DO 368 END DO 369 ! converte into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk 370 zk(:,:) = REAL( miku(:,:), wp ) ; CALL lbc_lnk( 'domzgr', zk, 'U', 1. ) ; miku(:,:) = MAX( NINT( zk(:,:) ), 1 ) 371 zk(:,:) = REAL( mikv(:,:), wp ) ; CALL lbc_lnk( 'domzgr', zk, 'V', 1. ) ; mikv(:,:) = MAX( NINT( zk(:,:) ), 1 ) 372 zk(:,:) = REAL( mikf(:,:), wp ) ; CALL lbc_lnk( 'domzgr', zk, 'F', 1. ) ; mikf(:,:) = MAX( NINT( zk(:,:) ), 1 ) 373 ! 374 zk(:,:) = REAL( mbku(:,:), wp ) ; CALL lbc_lnk( 'domzgr', zk, 'U', 1. ) ; mbku(:,:) = MAX( NINT( zk(:,:) ), 1 ) 375 zk(:,:) = REAL( mbkv(:,:), wp ) ; CALL lbc_lnk( 'domzgr', zk, 'V', 1. ) ; mbkv(:,:) = MAX( NINT( zk(:,:) ), 1 ) 376 ! 377 END SUBROUTINE zgr_top_bot 186 378 187 379 SUBROUTINE zgr_z … … 665 857 ENDIF 666 858 667 zbathy(:,:) = FLOAT( mbathy(:,:) ) 668 CALL lbc_lnk( 'domzgr',zbathy, 'T', 1._wp ) 669 mbathy(:,:) = INT( zbathy(:,:) ) 670 859 IF( lk_mpp ) THEN 860 zbathy(:,:) = FLOAT( mbathy(:,:) ) 861 CALL lbc_lnk( 'toto',zbathy, 'T', 1._wp ) 862 mbathy(:,:) = INT( zbathy(:,:) ) 863 ENDIF 671 864 ! ! East-west cyclic boundary conditions 672 865 IF( jperio == 0 ) THEN … … 1061 1254 END SUBROUTINE zgr_zps 1062 1255 1256 1063 1257 SUBROUTINE zgr_sco 1064 1258 !!---------------------------------------------------------------------- … … 1262 1456 END DO 1263 1457 ! apply lateral boundary condition CAUTION: keep the value when the lbc field is zero 1264 CALL lbc_lnk( ' domzgr',zenv, 'T', 1._wp, 'no0' )1458 CALL lbc_lnk( 'toto',zenv, 'T', 1._wp, 'no0' ) 1265 1459 ! ! ================ ! 1266 1460 END DO ! End loop !
Note: See TracChangeset
for help on using the changeset viewer.
