Changeset 5609 for trunk/NEMOGCM/TOOLS/SIREN/src/vgrid.f90
- Timestamp:
- 2015-07-17T17:42:15+02:00 (9 years ago)
- File:
-
- 1 edited
-
trunk/NEMOGCM/TOOLS/SIREN/src/vgrid.f90 (modified) (15 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/TOOLS/SIREN/src/vgrid.f90
r5037 r5609 70 70 !> @date Spetember, 2014 71 71 !> - add header 72 !> @date June, 2015 - update subroutine with NEMO 3.6 72 73 !> 73 74 !> @note Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) … … 139 140 !------------------------------------------------------------------- 140 141 SUBROUTINE vgrid_zgr_z( dd_gdepw, dd_gdept, dd_e3w, dd_e3t, & 142 & dd_e3w_1d, dd_e3t_1d, & 141 143 & dd_ppkth, dd_ppkth2, dd_ppacr, dd_ppacr2, & 142 144 & dd_ppdzmin, dd_pphmax, dd_pp_to_be_computed, & … … 148 150 REAL(dp), DIMENSION(:), INTENT(INOUT) :: dd_e3w 149 151 REAL(dp), DIMENSION(:), INTENT(INOUT) :: dd_e3t 152 REAL(dp), DIMENSION(:), INTENT(INOUT) :: dd_e3w_1d 153 REAL(dp), DIMENSION(:), INTENT(INOUT) :: dd_e3t_1d 150 154 151 155 REAL(dp) , INTENT(IN ) :: dd_ppkth … … 226 230 DO jk = 1, il_jpk 227 231 dl_zw = REAL(jk,dp) 228 dl_zt = REAL(jk,dp) + 0.5 232 dl_zt = REAL(jk,dp) + 0.5_dp 229 233 dd_gdepw(jk) = ( dl_zw - 1.0 ) * dl_za1 230 234 dd_gdept(jk) = ( dl_zt - 1.0 ) * dl_za1 … … 237 241 DO jk = 1, il_jpk 238 242 dl_zw = REAL( jk,dp) 239 dl_zt = REAL( jk,dp) + 0.5 243 dl_zt = REAL( jk,dp) + 0.5_dp 240 244 dd_gdepw(jk) = ( dl_zsur + dl_za0 * dl_zw + & 241 245 & dl_za1 * dl_zacr * LOG( COSH( (dl_zw-dl_zkth)/dl_zacr ) ) + & … … 255 259 ENDIF 256 260 261 ! need to be like this to compute the pressure gradient with ISF. 262 ! If not, level beneath the ISF are not aligned (sum(e3t) /= depth) 263 ! define e3t_0 and e3w_0 as the differences between gdept and gdepw respectively 264 DO jk = 1, il_jpk-1 265 dd_e3t_1d(jk) = dd_gdepw(jk+1)-dd_gdepw(jk) 266 END DO 267 dd_e3t_1d(il_jpk) = dd_e3t_1d(il_jpk-1) ! we don't care because this level is masked in NEMO 268 269 DO jk = 2, il_jpk 270 dd_e3w_1d(jk) = dd_gdept(jk) - dd_gdept(jk-1) 271 END DO 272 dd_e3w_1d(1 ) = 2._dp * (dd_gdept(1) - dd_gdepw(1)) 273 257 274 ! Control and print 258 275 ! ================== … … 260 277 DO jk = 1, il_jpk 261 278 IF( dd_e3w(jk) <= 0. .OR. dd_e3t(jk) <= 0. )then 262 CALL logger_debug("VGRID ZGR Z: e3w or e3t =<0 ")279 CALL logger_debug("VGRID ZGR Z: e3w or e3t <= 0 ") 263 280 ENDIF 281 282 IF( dd_e3w_1d(jk) <= 0. .OR. dd_e3t_1d(jk) <= 0. )then 283 CALL logger_debug("VGRID ZGR Z: e3w_1d or e3t_1d <= 0 ") 284 ENDIF 264 285 265 286 IF( dd_gdepw(jk) < 0. .OR. dd_gdept(jk) < 0. )then … … 269 290 270 291 END SUBROUTINE vgrid_zgr_z 292 !------------------------------------------------------------------- 293 !------------------------------------------------------------------- 294 SUBROUTINE vgrid_zgr_bat( dd_bathy, dd_gdepw, dd_hmin, dd_fill ) 295 IMPLICIT NONE 296 ! Argument 297 REAL(dp), DIMENSION(:,:), INTENT(INOUT) :: dd_bathy 298 REAL(dp), DIMENSION(:) , INTENT(IN ) :: dd_gdepw 299 REAL(dp) , INTENT(IN ) :: dd_hmin 300 REAL(dp) , INTENT(IN ), OPTIONAL :: dd_fill 301 302 ! local 303 INTEGER(i4) :: il_jpk 304 305 REAL(dp) :: dl_hmin 306 REAL(dp) :: dl_fill 307 308 ! loop indices 309 INTEGER(i4) :: jk 310 !---------------------------------------------------------------- 311 il_jpk = SIZE(dd_gdepw(:)) 312 313 dl_fill=0._dp 314 IF( PRESENT(dd_fill) ) dl_fill=dd_fill 315 316 IF( dd_hmin < 0._dp ) THEN 317 jk = - INT( dd_hmin ) ! from a nb of level 318 ELSE 319 jk = MINLOC( dd_gdepw, mask = dd_gdepw > dd_hmin, dim = 1 ) ! from a depth 320 ENDIF 321 322 dl_hmin = dd_gdepw(jk+1) ! minimum depth = ik+1 w-levels 323 WHERE( dd_bathy(:,:) <= 0._wp .OR. dd_bathy(:,:) == dl_fill ) 324 dd_bathy(:,:) = dl_fill ! min=0 over the lands 325 ELSE WHERE 326 dd_bathy(:,:) = MAX( dl_hmin , dd_bathy(:,:) ) ! min=dl_hmin over the oceans 327 END WHERE 328 WRITE(*,*) 'Minimum ocean depth: ', dl_hmin, ' minimum number of ocean levels : ', jk 329 330 END SUBROUTINE vgrid_zgr_bat 271 331 !------------------------------------------------------------------- 272 332 !> @brief This subroutine set the depth and vertical scale factor in partial step … … 327 387 !------------------------------------------------------------------- 328 388 SUBROUTINE vgrid_zgr_zps( id_mbathy, dd_bathy, id_jpkmax, & 329 & dd_gdepw, dd_e3t, & 330 & dd_e3zps_min, dd_e3zps_rat ) 389 & dd_gdepw, dd_e3t, & 390 & dd_e3zps_min, dd_e3zps_rat, & 391 & dd_fill ) 331 392 IMPLICIT NONE 332 393 ! Argument … … 336 397 REAL(dp) , DIMENSION(:) , INTENT(IN ) :: dd_gdepw 337 398 REAL(dp) , DIMENSION(:) , INTENT(IN ) :: dd_e3t 338 REAL(dp) :: dd_e3zps_min 339 REAL(dp) :: dd_e3zps_rat 399 REAL(dp) , INTENT(IN ) :: dd_e3zps_min 400 REAL(dp) , INTENT(IN ) :: dd_e3zps_rat 401 REAL(dp) , INTENT(IN ), OPTIONAL :: dd_fill 340 402 341 403 ! local variable 342 404 REAL(dp) :: dl_zmax ! Maximum depth 343 REAL(dp) :: dl_zmin ! Minimum depth405 !REAL(dp) :: dl_zmin ! Minimum depth 344 406 REAL(dp) :: dl_zdepth ! Ajusted ocean depth to avoid too small e3t 407 REAL(dp) :: dl_fill 345 408 346 409 INTEGER(i4) :: il_jpk … … 359 422 il_jpjglo=SIZE(id_mbathy(:,:),DIM=2) 360 423 424 dl_fill=0._dp 425 IF( PRESENT(dd_fill) ) dl_fill=dd_fill 426 361 427 ! Initialization of constant 362 dl_zmax = dd_gdepw(il_jpk) + dd_e3t(il_jpk) 363 dl_zmin = dd_gdepw(4) 428 dl_zmax = dd_gdepw(il_jpk) + dd_e3t(il_jpk) ! maximum depth (i.e. the last ocean level thickness <= 2*e3t_1d(jpkm1) ) 429 430 ! bounded value of bathy (min already set at the end of zgr_bat) 431 WHERE( dd_bathy(:,:) /= dl_fill ) 432 dd_bathy(:,:) = MIN( dl_zmax , dd_bathy(:,:) ) 433 END WHERE 364 434 365 435 ! bathymetry in level (from bathy_meter) … … 372 442 DO jj = 1, il_jpjglo 373 443 DO ji= 1, il_jpiglo 374 IF( dd_bathy(ji,jj) <= 0. ) id_mbathy(ji,jj) = INT(dd_bathy(ji,jj),i4) 375 END DO 376 END DO 377 378 ! bounded value of bathy 379 ! minimum depth == 3 levels 380 ! maximum depth == gdepw(jpk)+e3t(jpk) 381 ! i.e. the last ocean level thickness cannot exceed e3t(jpkm1)+e3t(jpk) 382 DO jj = 1, il_jpjglo 383 DO ji= 1, il_jpiglo 384 IF( dd_bathy(ji,jj) <= 0. ) THEN 385 dd_bathy(ji,jj) = 0.e0 386 ELSE 387 dd_bathy(ji,jj) = MAX( dd_bathy(ji,jj), dl_zmin ) 388 dd_bathy(ji,jj) = MIN( dd_bathy(ji,jj), dl_zmax ) 444 IF( dd_bathy(ji,jj) <= 0._dp )THEN 445 id_mbathy(ji,jj) = INT(dd_bathy(ji,jj),i4) 446 ELSEIF( dd_bathy(ji,jj) == dl_fill )THEN 447 id_mbathy(ji,jj) = 0_i4 389 448 ENDIF 390 449 END DO … … 401 460 DO jj = 1, il_jpjglo 402 461 DO ji = 1, il_jpiglo 403 IF( 0. < dd_bathy(ji,jj) .AND. dd_bathy(ji,jj) <= dl_zdepth ) id_mbathy(ji,jj) = jk-1 462 IF( dd_bathy(ji,jj) /= dl_fill )THEN 463 IF( 0. < dd_bathy(ji,jj) .AND. & 464 & dd_bathy(ji,jj) <= dl_zdepth ) id_mbathy(ji,jj) = jk-1 465 ENDIF 404 466 END DO 405 467 END DO … … 567 629 REAL(dp) , DIMENSION(:) , ALLOCATABLE :: dl_e3w 568 630 REAL(dp) , DIMENSION(:) , ALLOCATABLE :: dl_e3t 631 REAL(dp) , DIMENSION(:) , ALLOCATABLE :: dl_e3w_1d 632 REAL(dp) , DIMENSION(:) , ALLOCATABLE :: dl_e3t_1d 569 633 570 634 INTEGER(i4) :: il_status … … 710 774 ALLOCATE( dl_gdepw(in_nlevel), dl_gdept(in_nlevel) ) 711 775 ALLOCATE( dl_e3w(in_nlevel), dl_e3t(in_nlevel) ) 776 ALLOCATE( dl_e3w_1d(in_nlevel), dl_e3t_1d(in_nlevel) ) 712 777 CALL vgrid_zgr_z( dl_gdepw(:), dl_gdept(:), dl_e3w(:), dl_e3t(:), & 778 & dl_e3w_1d, dl_e3t_1d, & 713 779 & dn_ppkth, dn_ppkth2, dn_ppacr, dn_ppacr2, & 714 780 & dn_ppdzmin, dn_pphmax, dn_pp_to_be_computed, &
Note: See TracChangeset
for help on using the changeset viewer.
