- Timestamp:
- 2016-04-07T16:32:24+02:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5518_GC3p0_package/NEMOGCM/TOOLS/SIREN/src/vgrid.f90
r5037 r6440 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) … … 118 119 !> 119 120 !> @author G. Madec 120 !> - 03,08- G. Madec:F90: Free form and module121 !> @date Marsh,2008 - F90: Free form and module 121 122 ! 122 123 !> @note Reference : Marti, Madec & Delecluse, 1992, JGR, 97, No8, 12,763-12,766. … … 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 !> @brief This subroutine 294 !> 295 !> @todo add subroutine description 296 !> 297 !> @param[inout] dd_bathy 298 !> @param[in] dd_gdepw 299 !> @param[in] dd_hmin 300 !> @param[in] dd_fill 301 !------------------------------------------------------------------- 302 SUBROUTINE vgrid_zgr_bat( dd_bathy, dd_gdepw, dd_hmin, dd_fill ) 303 IMPLICIT NONE 304 ! Argument 305 REAL(dp), DIMENSION(:,:), INTENT(INOUT) :: dd_bathy 306 REAL(dp), DIMENSION(:) , INTENT(IN ) :: dd_gdepw 307 REAL(dp) , INTENT(IN ) :: dd_hmin 308 REAL(dp) , INTENT(IN ), OPTIONAL :: dd_fill 309 310 ! local 311 INTEGER(i4) :: il_jpk 312 313 REAL(dp) :: dl_hmin 314 REAL(dp) :: dl_fill 315 316 ! loop indices 317 INTEGER(i4) :: jk 318 !---------------------------------------------------------------- 319 il_jpk = SIZE(dd_gdepw(:)) 320 321 dl_fill=0._dp 322 IF( PRESENT(dd_fill) ) dl_fill=dd_fill 323 324 IF( dd_hmin < 0._dp ) THEN 325 jk = - INT( dd_hmin ) ! from a nb of level 326 ELSE 327 jk = MINLOC( dd_gdepw, mask = dd_gdepw > dd_hmin, dim = 1 ) ! from a depth 328 ENDIF 329 330 dl_hmin = dd_gdepw(jk+1) ! minimum depth = ik+1 w-levels 331 WHERE( dd_bathy(:,:) <= 0._wp .OR. dd_bathy(:,:) == dl_fill ) 332 dd_bathy(:,:) = dl_fill ! min=0 over the lands 333 ELSE WHERE 334 dd_bathy(:,:) = MAX( dl_hmin , dd_bathy(:,:) ) ! min=dl_hmin over the oceans 335 END WHERE 336 WRITE(*,*) 'Minimum ocean depth: ', dl_hmin, ' minimum number of ocean levels : ', jk 337 338 END SUBROUTINE vgrid_zgr_bat 271 339 !------------------------------------------------------------------- 272 340 !> @brief This subroutine set the depth and vertical scale factor in partial step … … 311 379 !> - gdept, gdepw and e3 are positives 312 380 !> - gdept_ps, gdepw_ps and e3_ps are positives 313 ! 381 !> 314 382 !> @author A. Bozec, G. Madec 315 !> - 02-09 (A. Bozec, G. Madec) F90: Free form and module 316 !> - 02-09 (A. de Miranda) rigid-lid + islands 383 !> @date February, 2009 - F90: Free form and module 384 !> @date February, 2009 385 !> - A. de Miranda : rigid-lid + islands 317 386 !> 318 387 !> @note Reference : Pacanowsky & Gnanadesikan 1997, Mon. Wea. Rev., 126, 3248-3270. … … 325 394 !> @param[in] dd_e3zps_min 326 395 !> @param[in] dd_e3zps_rat 396 !> @param[in] dd_fill 327 397 !------------------------------------------------------------------- 328 398 SUBROUTINE vgrid_zgr_zps( id_mbathy, dd_bathy, id_jpkmax, & 329 & dd_gdepw, dd_e3t, & 330 & dd_e3zps_min, dd_e3zps_rat ) 399 & dd_gdepw, dd_e3t, & 400 & dd_e3zps_min, dd_e3zps_rat, & 401 & dd_fill ) 331 402 IMPLICIT NONE 332 403 ! Argument … … 336 407 REAL(dp) , DIMENSION(:) , INTENT(IN ) :: dd_gdepw 337 408 REAL(dp) , DIMENSION(:) , INTENT(IN ) :: dd_e3t 338 REAL(dp) :: dd_e3zps_min 339 REAL(dp) :: dd_e3zps_rat 409 REAL(dp) , INTENT(IN ) :: dd_e3zps_min 410 REAL(dp) , INTENT(IN ) :: dd_e3zps_rat 411 REAL(dp) , INTENT(IN ), OPTIONAL :: dd_fill 340 412 341 413 ! local variable 342 414 REAL(dp) :: dl_zmax ! Maximum depth 343 REAL(dp) :: dl_zmin ! Minimum depth415 !REAL(dp) :: dl_zmin ! Minimum depth 344 416 REAL(dp) :: dl_zdepth ! Ajusted ocean depth to avoid too small e3t 417 REAL(dp) :: dl_fill 345 418 346 419 INTEGER(i4) :: il_jpk … … 359 432 il_jpjglo=SIZE(id_mbathy(:,:),DIM=2) 360 433 434 dl_fill=0._dp 435 IF( PRESENT(dd_fill) ) dl_fill=dd_fill 436 361 437 ! Initialization of constant 362 dl_zmax = dd_gdepw(il_jpk) + dd_e3t(il_jpk) 363 dl_zmin = dd_gdepw(4) 438 dl_zmax = dd_gdepw(il_jpk) + dd_e3t(il_jpk) ! maximum depth (i.e. the last ocean level thickness <= 2*e3t_1d(jpkm1) ) 439 440 ! bounded value of bathy (min already set at the end of zgr_bat) 441 WHERE( dd_bathy(:,:) /= dl_fill ) 442 dd_bathy(:,:) = MIN( dl_zmax , dd_bathy(:,:) ) 443 END WHERE 364 444 365 445 ! bathymetry in level (from bathy_meter) … … 372 452 DO jj = 1, il_jpjglo 373 453 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 ) 454 IF( dd_bathy(ji,jj) <= 0._dp )THEN 455 id_mbathy(ji,jj) = INT(dd_bathy(ji,jj),i4) 456 ELSEIF( dd_bathy(ji,jj) == dl_fill )THEN 457 id_mbathy(ji,jj) = 0_i4 389 458 ENDIF 390 459 END DO … … 401 470 DO jj = 1, il_jpjglo 402 471 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 472 IF( dd_bathy(ji,jj) /= dl_fill )THEN 473 IF( 0. < dd_bathy(ji,jj) .AND. & 474 & dd_bathy(ji,jj) <= dl_zdepth ) id_mbathy(ji,jj) = jk-1 475 ENDIF 404 476 END DO 405 477 END DO … … 432 504 !> ** Action : - update mbathy: level bathymetry (in level index) 433 505 !> - update bathy : meter bathymetry (in meters) 434 506 !> 435 507 !> @author G.Madec 436 !> - 03-08Original code508 !> @date Marsh, 2008 - Original code 437 509 ! 438 510 !> @param[in] id_mbathy … … 543 615 !> 544 616 !> @author J.Paul 545 !> - November, 2013- Initial Version617 !> @date November, 2013 - Initial Version 546 618 ! 547 619 !> @param[in] td_bathy Bathymetry file structure … … 567 639 REAL(dp) , DIMENSION(:) , ALLOCATABLE :: dl_e3w 568 640 REAL(dp) , DIMENSION(:) , ALLOCATABLE :: dl_e3t 641 REAL(dp) , DIMENSION(:) , ALLOCATABLE :: dl_e3w_1d 642 REAL(dp) , DIMENSION(:) , ALLOCATABLE :: dl_e3t_1d 569 643 570 644 INTEGER(i4) :: il_status … … 710 784 ALLOCATE( dl_gdepw(in_nlevel), dl_gdept(in_nlevel) ) 711 785 ALLOCATE( dl_e3w(in_nlevel), dl_e3t(in_nlevel) ) 786 ALLOCATE( dl_e3w_1d(in_nlevel), dl_e3t_1d(in_nlevel) ) 712 787 CALL vgrid_zgr_z( dl_gdepw(:), dl_gdept(:), dl_e3w(:), dl_e3t(:), & 788 & dl_e3w_1d, dl_e3t_1d, & 713 789 & dn_ppkth, dn_ppkth2, dn_ppacr, dn_ppacr2, & 714 790 & dn_ppdzmin, dn_pphmax, dn_pp_to_be_computed, &
Note: See TracChangeset
for help on using the changeset viewer.