- Timestamp:
- 2016-10-03T14:09:01+02:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/TOOLS/DOMAINcfg/src/domwri.f90
r6951 r6973 13 13 !! dom_wri : create and write mesh and mask file(s) 14 14 !! dom_uniq : identify unique point of a grid (TUVF) 15 !! dom_stiff : diagnose maximum grid stiffness/hydrostatic consistency (s-coordinate) 15 16 !!---------------------------------------------------------------------- 16 17 USE dom_oce ! ocean space and time domain … … 27 28 PUBLIC dom_wri ! routine called by inidom.F90 28 29 PUBLIC dom_wri_coordinate ! routine called by domhgr.F90 29 !! * Substitutions 30 !!---------------------------------------------------------------------- 31 !! *** vectopt_loop_substitute *** 32 !!---------------------------------------------------------------------- 33 !! ** purpose : substitute the inner loop start/end indices with CPP macro 34 !! allow unrolling of do-loop (useful with vector processors) 35 !!---------------------------------------------------------------------- 30 PUBLIC dom_stiff ! routine called by inidom.F90 31 36 32 !!---------------------------------------------------------------------- 37 33 !! NEMO/OPA 3.7 , NEMO Consortium (2014) 38 34 !! $Id: vectopt_loop_substitute.h90 4990 2014-12-15 16:42:49Z timgraham $ 39 35 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 40 !!----------------------------------------------------------------------41 !!----------------------------------------------------------------------42 !! NEMO/OPA 3.3 , NEMO Consortium (2010)43 !! $Id: domwri.F90 6689 2016-06-13 14:24:43Z flavoni $44 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)45 36 !!---------------------------------------------------------------------- 46 37 CONTAINS … … 286 277 CALL iom_rstput( 0, 0, inum4, 'e3v_0', e3v_0 ) 287 278 CALL iom_rstput( 0, 0, inum4, 'e3w_0', e3w_0 ) 288 CALL iom_rstput( 0, 0, inum4, 'rx1', rx1 ) ! ! Max. grid stiffness ratio289 279 ! 290 280 CALL iom_rstput( 0, 0, inum4, 'gdept_1d' , gdept_1d ) ! ! stretched system … … 292 282 CALL iom_rstput( 0, 0, inum4, 'gdept_0', gdept_0, ktype = jp_r8 ) 293 283 CALL iom_rstput( 0, 0, inum4, 'gdepw_0', gdepw_0, ktype = jp_r8 ) 284 CALL dom_stiff( zprt ) 285 CALL iom_rstput( 0, 0, inum, 'stiffness', zprt ) ! ! Max. grid stiffness ratio 294 286 ENDIF 295 287 … … 416 408 END SUBROUTINE dom_uniq 417 409 410 411 SUBROUTINE dom_stiff( px1 ) 412 !!---------------------------------------------------------------------- 413 !! *** ROUTINE dom_stiff *** 414 !! 415 !! ** Purpose : Diagnose maximum grid stiffness/hydrostatic consistency 416 !! 417 !! ** Method : Compute Haney (1991) hydrostatic condition ratio 418 !! Save the maximum in the vertical direction 419 !! (this number is only relevant in s-coordinates) 420 !! 421 !! Haney, 1991, J. Phys. Oceanogr., 21, 610-619. 422 !!---------------------------------------------------------------------- 423 REAL(wp), DIMENSION(:,:), INTENT(out), OPTIONAL :: px1 ! stiffness 424 ! 425 INTEGER :: ji, jj, jk 426 REAL(wp) :: zrxmax 427 REAL(wp), DIMENSION(4) :: zr1 428 REAL(wp), DIMENSION(jpi,jpj) :: zx1 429 !!---------------------------------------------------------------------- 430 zx1(:,:) = 0._wp 431 zrxmax = 0._wp 432 zr1(:) = 0._wp 433 ! 434 DO ji = 2, jpim1 435 DO jj = 2, jpjm1 436 DO jk = 1, jpkm1 437 !!gm remark: dk(gdepw) = e3t ===>>> possible simplification of the following calculation.... 438 !! especially since it is gde3w which is used to compute the pressure gradient 439 !! furthermore, I think gdept_0 should be used below instead of w point in the numerator 440 !! so that the ratio is computed at the same point (i.e. uw and vw) .... 441 zr1(1) = ABS( ( gdepw_0(ji ,jj,jk )-gdepw_0(ji-1,jj,jk ) & 442 & +gdepw_0(ji ,jj,jk+1)-gdepw_0(ji-1,jj,jk+1) ) & 443 & / ( gdepw_0(ji ,jj,jk )+gdepw_0(ji-1,jj,jk ) & 444 & -gdepw_0(ji ,jj,jk+1)-gdepw_0(ji-1,jj,jk+1) + rsmall ) ) * umask(ji-1,jj,jk) 445 zr1(2) = ABS( ( gdepw_0(ji+1,jj,jk )-gdepw_0(ji ,jj,jk ) & 446 & +gdepw_0(ji+1,jj,jk+1)-gdepw_0(ji ,jj,jk+1) ) & 447 & / ( gdepw_0(ji+1,jj,jk )+gdepw_0(ji ,jj,jk ) & 448 & -gdepw_0(ji+1,jj,jk+1)-gdepw_0(ji ,jj,jk+1) + rsmall ) ) * umask(ji ,jj,jk) 449 zr1(3) = ABS( ( gdepw_0(ji,jj+1,jk )-gdepw_0(ji,jj ,jk ) & 450 & +gdepw_0(ji,jj+1,jk+1)-gdepw_0(ji,jj ,jk+1) ) & 451 & / ( gdepw_0(ji,jj+1,jk )+gdepw_0(ji,jj ,jk ) & 452 & -gdepw_0(ji,jj+1,jk+1)-gdepw_0(ji,jj ,jk+1) + rsmall ) ) * vmask(ji,jj ,jk) 453 zr1(4) = ABS( ( gdepw_0(ji,jj ,jk )-gdepw_0(ji,jj-1,jk ) & 454 & +gdepw_0(ji,jj ,jk+1)-gdepw_0(ji,jj-1,jk+1) ) & 455 & / ( gdepw_0(ji,jj ,jk )+gdepw_0(ji,jj-1,jk ) & 456 & -gdepw_0(ji,jj ,jk+1)-gdepw_0(ji,jj-1,jk+1) + rsmall ) ) * vmask(ji,jj-1,jk) 457 zrxmax = MAXVAL( zr1(1:4) ) 458 zx1(ji,jj) = MAX( zx1(ji,jj) , zrxmax ) 459 END DO 460 END DO 461 END DO 462 CALL lbc_lnk( zx1, 'T', 1. ) 463 ! 464 IF( PRESENT( px1 ) ) px1 = zx1 465 ! 466 zrxmax = MAXVAL( zx1 ) 467 ! 468 IF( lk_mpp ) CALL mpp_max( zrxmax ) ! max over the global domain 469 ! 470 IF(lwp) THEN 471 WRITE(numout,*) 472 WRITE(numout,*) 'dom_stiff : maximum grid stiffness ratio: ', zrxmax 473 WRITE(numout,*) '~~~~~~~~~' 474 ENDIF 475 ! 476 END SUBROUTINE dom_stiff 477 418 478 !!====================================================================== 419 479 END MODULE domwri
Note: See TracChangeset
for help on using the changeset viewer.
