Changeset 6973


Ignore:
Timestamp:
2016-10-03T14:09:01+02:00 (5 years ago)
Author:
flavoni
Message:

fix small bug

Location:
branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90

    r6904 r6973  
    106106      CALL dom_msk( ik_top, ik_bot )   ! Masks 
    107107      ! 
    108       IF( ln_sco )  CALL dom_stiff     ! Maximum stiffness ratio/hydrostatic consistency 
    109       ! 
    110108      DO jj = 1, jpj                   ! depth of the iceshelves 
    111109         DO ji = 1, jpj 
  • branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OPA_SRC/DOM/domwri.F90

    r6923 r6973  
    171171      CALL iom_rstput( 0, 0, inum, 'e2f', e2f, ktype = jp_r8 ) 
    172172       
    173       CALL iom_rstput( 0, 0, inum, 'ff_f', ff_f, ktype = jp_r8 )           !    ! coriolis factor 
     173      CALL iom_rstput( 0, 0, inum, 'ff_f', ff_f, ktype = jp_r8 )       !    ! coriolis factor 
    174174      CALL iom_rstput( 0, 0, inum, 'ff_t', ff_t, ktype = jp_r8 ) 
    175175       
     
    180180      CALL iom_rstput( 0, 0, inum, 'misf', zprt, ktype = jp_i4 )       !    ! nb of ocean T-points 
    181181      zprt(:,:) = ssmask(:,:) * REAL( risfdep(:,:) , wp ) 
    182       CALL iom_rstput( 0, 0, inum, 'isfdraft', zprt, ktype = jp_r8 )       !    ! nb of ocean T-points 
     182      CALL iom_rstput( 0, 0, inum, 'isfdraft', zprt, ktype = jp_r8 )   !    ! nb of ocean T-points 
    183183             
    184       CALL iom_rstput( 0, 0, inum, 'e3t_0', e3t_0 )         !    ! scale factors 
    185       CALL iom_rstput( 0, 0, inum, 'e3u_0', e3u_0 ) 
    186       CALL iom_rstput( 0, 0, inum, 'e3v_0', e3v_0 ) 
    187       CALL iom_rstput( 0, 0, inum, 'e3w_0', e3w_0 ) 
    188       ! 
    189       CALL dom_stiff( zprt ) 
    190       CALL iom_rstput( 0, 0, inum, 'stiffness', zprt )      !    ! Max. grid stiffness ratio 
    191       ! 
    192       CALL iom_rstput( 0, 0, inum, 'gdept_1d' , gdept_1d )  !    ! stretched system 
    193       CALL iom_rstput( 0, 0, inum, 'gdepw_1d' , gdepw_1d ) 
    194       CALL iom_rstput( 0, 0, inum, 'gdept_0', gdept_0, ktype = jp_r8 ) 
    195       CALL iom_rstput( 0, 0, inum, 'gdepw_0', gdepw_0, ktype = jp_r8 ) 
    196        
     184      !                                                         ! vertical mesh 
     185      CALL iom_rstput( 0, 0, inum, 'e3t_0', e3t_0, ktype = jp_r8  )    !    ! scale factors 
     186      CALL iom_rstput( 0, 0, inum, 'e3u_0', e3u_0, ktype = jp_r8  ) 
     187      CALL iom_rstput( 0, 0, inum, 'e3v_0', e3v_0, ktype = jp_r8  ) 
     188      CALL iom_rstput( 0, 0, inum, 'e3w_0', e3w_0, ktype = jp_r8  ) 
     189      ! 
     190      CALL iom_rstput( 0, 0, inum, 'gdept_1d' , gdept_1d , ktype = jp_r8 )  ! stretched system 
     191      CALL iom_rstput( 0, 0, inum, 'gdepw_1d' , gdepw_1d , ktype = jp_r8 ) 
     192      CALL iom_rstput( 0, 0, inum, 'gdept_0'  , gdept_0  , ktype = jp_r8 ) 
     193      CALL iom_rstput( 0, 0, inum, 'gdepw_0'  , gdepw_0  , ktype = jp_r8 ) 
     194      ! 
     195      IF( ln_sco ) THEN                                         ! s-coordinate stiffness 
     196         CALL dom_stiff( zprt ) 
     197         CALL iom_rstput( 0, 0, inum, 'stiffness', zprt )      !    ! Max. grid stiffness ratio 
     198      ENDIF 
     199      ! 
    197200      !                                     ! ============================ 
    198201      CALL iom_close( inum )                !        close the files  
  • branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/TOOLS/DOMAINcfg/src/dom_oce.f90

    r6951 r6973  
    222222   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hifv  , hiff       !: interface depth between stretching at v--f 
    223223   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hift  , hifu       !: and quasi-uniform spacing             t--u points (m) 
    224    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   rx1                !: Maximum grid stiffness ratio 
    225224!!gm end 
    226225 
     
    362361         &      scosrf(jpi,jpj) , scobot(jpi,jpj) ,     & 
    363362         &      hifv  (jpi,jpj) , hiff  (jpi,jpj) ,     & 
    364          &      hift  (jpi,jpj) , hifu  (jpi,jpj) , rx1(jpi,jpj) , STAT=ierr(8) ) 
     363         &      hift  (jpi,jpj) , hifu  (jpi,jpj) , STAT=ierr(8) ) 
    365364 
    366365      ALLOCATE( mbathy(jpi,jpj) , bathy  (jpi,jpj) ,                                       & 
  • branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/TOOLS/DOMAINcfg/src/domain.f90

    r6951 r6973  
    2020   !!   dom_nam        : read and contral domain namelists 
    2121   !!   dom_ctl        : control print for the ocean domain 
    22    !!   dom_stiff      : diagnose maximum grid stiffness/hydrostatic consistency (s-coordinate) 
    2322   !!---------------------------------------------------------------------- 
    2423   USE oce             ! ocean variables 
     
    8887                     CALL dom_zgr               ! Vertical mesh and bathymetry 
    8988                     CALL dom_msk               ! Masks 
    90       IF( ln_sco )   CALL dom_stiff             ! Maximum stiffness ratio/hydrostatic consistency 
    9189      ! 
    9290      ht_0(:,:) = e3t_0(:,:,1) * tmask(:,:,1)   ! Reference ocean thickness 
     
    370368 
    371369 
    372    SUBROUTINE dom_stiff 
    373       !!---------------------------------------------------------------------- 
    374       !!                  ***  ROUTINE dom_stiff  *** 
    375       !!                      
    376       !! ** Purpose :   Diagnose maximum grid stiffness/hydrostatic consistency 
    377       !! 
    378       !! ** Method  :   Compute Haney (1991) hydrostatic condition ratio 
    379       !!                Save the maximum in the vertical direction 
    380       !!                (this number is only relevant in s-coordinates) 
    381       !! 
    382       !!                Haney, R. L., 1991: On the pressure gradient force 
    383       !!                over steep topography in sigma coordinate ocean models.  
    384       !!                J. Phys. Oceanogr., 21, 610???619. 
    385       !!---------------------------------------------------------------------- 
    386       INTEGER  ::   ji, jj, jk  
    387       REAL(wp) ::   zrxmax 
    388       REAL(wp), DIMENSION(4) ::   zr1 
    389       !!---------------------------------------------------------------------- 
    390       rx1(:,:) = 0._wp 
    391       zrxmax   = 0._wp 
    392       zr1(:)   = 0._wp 
    393       ! 
    394       DO ji = 2, jpim1 
    395          DO jj = 2, jpjm1 
    396             DO jk = 1, jpkm1 
    397                zr1(1) = ABS(  ( gdepw_0(ji  ,jj,jk  )-gdepw_0(ji-1,jj,jk  )               &  
    398                     &          +gdepw_0(ji  ,jj,jk+1)-gdepw_0(ji-1,jj,jk+1) )             & 
    399                     &       / ( gdepw_0(ji  ,jj,jk  )+gdepw_0(ji-1,jj,jk  )               & 
    400                     &          -gdepw_0(ji  ,jj,jk+1)-gdepw_0(ji-1,jj,jk+1) + rsmall )  ) * umask(ji-1,jj,jk) 
    401                zr1(2) = ABS(  ( gdepw_0(ji+1,jj,jk  )-gdepw_0(ji  ,jj,jk  )               & 
    402                     &          +gdepw_0(ji+1,jj,jk+1)-gdepw_0(ji  ,jj,jk+1) )             & 
    403                     &       / ( gdepw_0(ji+1,jj,jk  )+gdepw_0(ji  ,jj,jk  )               & 
    404                     &          -gdepw_0(ji+1,jj,jk+1)-gdepw_0(ji  ,jj,jk+1) + rsmall )  ) * umask(ji  ,jj,jk) 
    405                zr1(3) = ABS(  ( gdepw_0(ji,jj+1,jk  )-gdepw_0(ji,jj  ,jk  )               & 
    406                     &          +gdepw_0(ji,jj+1,jk+1)-gdepw_0(ji,jj  ,jk+1) )             & 
    407                     &       / ( gdepw_0(ji,jj+1,jk  )+gdepw_0(ji,jj  ,jk  )               & 
    408                     &          -gdepw_0(ji,jj+1,jk+1)-gdepw_0(ji,jj  ,jk+1) + rsmall )  ) * vmask(ji,jj  ,jk) 
    409                zr1(4) = ABS(  ( gdepw_0(ji,jj  ,jk  )-gdepw_0(ji,jj-1,jk  )               & 
    410                     &          +gdepw_0(ji,jj  ,jk+1)-gdepw_0(ji,jj-1,jk+1) )             & 
    411                     &       / ( gdepw_0(ji,jj  ,jk  )+gdepw_0(ji,jj-1,jk  )               & 
    412                     &          -gdepw_0(ji,jj  ,jk+1)-gdepw_0(ji,jj-1,jk+1) + rsmall )  ) * vmask(ji,jj-1,jk) 
    413                zrxmax = MAXVAL( zr1(1:4) ) 
    414                rx1(ji,jj) = MAX( rx1(ji,jj) , zrxmax ) 
    415             END DO 
    416          END DO 
    417       END DO 
    418       CALL lbc_lnk( rx1, 'T', 1. ) 
    419       ! 
    420       zrxmax = MAXVAL( rx1 ) 
    421       ! 
    422       IF( lk_mpp )   CALL mpp_max( zrxmax ) ! max over the global domain 
    423       ! 
    424       IF(lwp) THEN 
    425          WRITE(numout,*) 
    426          WRITE(numout,*) 'dom_stiff : maximum grid stiffness ratio: ', zrxmax 
    427          WRITE(numout,*) '~~~~~~~~~' 
    428       ENDIF 
    429       ! 
    430    END SUBROUTINE dom_stiff 
    431  
    432  
    433370   SUBROUTINE cfg_write 
    434371      !!---------------------------------------------------------------------- 
     
    537474      CALL iom_rstput( 0, 0, inum, 'top    level' , REAL( mikt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points (ISF) 
    538475      ! 
    539       IF( ln_sco ) THEN             ! s-coordinate: store grid stiffness ratio (Not required anyway) 
    540         !SF  CALL dom_stiff( z2d ) !commented because at compilation error:  The 
    541                                    !number of actual arguments cannot be greater than the number of dummy 
    542                                    !arguments.   [DOM_STIFF] 
    543                                    !CALL dom_stiff( z2d ) 
    544                                    ! --------------^ compilation aborted for 
    545                                    !/workgpfs/rech/gzi/rgzi011/commit-simplif2-TOOLS/NEMOGCM/CONFIG/TEST/BLD/ppsrc/nemo/domain.f90 
     476      IF( ln_sco ) THEN             ! s-coordinate: store grid stiffness ratio  (Not required anyway) 
     477         CALL dom_stiff( z2d ) 
    546478         CALL iom_rstput( 0, 0, inum, 'stiffness', z2d )        !    ! Max. grid stiffness ratio 
    547479      ENDIF 
  • branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/TOOLS/DOMAINcfg/src/domwri.f90

    r6951 r6973  
    1313   !!   dom_wri        : create and write mesh and mask file(s) 
    1414   !!   dom_uniq       : identify unique point of a grid (TUVF) 
     15   !!   dom_stiff      : diagnose maximum grid stiffness/hydrostatic consistency (s-coordinate) 
    1516   !!---------------------------------------------------------------------- 
    1617   USE dom_oce         ! ocean space and time domain 
     
    2728   PUBLIC   dom_wri              ! routine called by inidom.F90 
    2829   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 
    3632   !!---------------------------------------------------------------------- 
    3733   !! NEMO/OPA 3.7 , NEMO Consortium (2014) 
    3834   !! $Id: vectopt_loop_substitute.h90 4990 2014-12-15 16:42:49Z timgraham $  
    3935   !! 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) 
    4536   !!---------------------------------------------------------------------- 
    4637CONTAINS 
     
    286277         CALL iom_rstput( 0, 0, inum4, 'e3v_0', e3v_0 ) 
    287278         CALL iom_rstput( 0, 0, inum4, 'e3w_0', e3w_0 ) 
    288          CALL iom_rstput( 0, 0, inum4, 'rx1', rx1 )             !    ! Max. grid stiffness ratio 
    289279         ! 
    290280         CALL iom_rstput( 0, 0, inum4, 'gdept_1d' , gdept_1d )  !    ! stretched system 
     
    292282         CALL iom_rstput( 0, 0, inum4, 'gdept_0', gdept_0, ktype = jp_r8 )      
    293283         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 
    294286      ENDIF 
    295287       
     
    416408   END SUBROUTINE dom_uniq 
    417409 
     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 
    418478   !!====================================================================== 
    419479END MODULE domwri 
Note: See TracChangeset for help on using the changeset viewer.