New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 1408 – NEMO

Changeset 1408


Ignore:
Timestamp:
2009-04-17T11:59:49+02:00 (15 years ago)
Author:
rblod
Message:

dev_004_VVL: cosmetic changes and code review

Location:
branches/dev_004_VVL/NEMO/OPA_SRC
Files:
16 edited

Legend:

Unmodified
Added
Removed
  • branches/dev_004_VVL/NEMO/OPA_SRC/DOM/dom_oce.F90

    r1385 r1408  
    127127   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &   !: 
    128128      hur, hvr,          &  !: inverse of u and v-points ocean depth (1/m) 
    129       hu , hv               !: depth at u- and v-points (meters) 
     129      hu , hv,           &  !: depth at u- and v-points (meters) 
     130      hu_0 , hv_0           !: refernce depth at u- and v-points (meters) 
    130131 
    131132   !! z-coordinate with full steps (also used in the other cases as reference z-coordinate) 
  • branches/dev_004_VVL/NEMO/OPA_SRC/DOM/domain.F90

    r1335 r1408  
    44   !! Ocean initialization : domain initialization 
    55   !!============================================================================== 
    6  
     6   !! History :  OPA  !  1990-10  (C. Levy - G. Madec)  Original code 
     7   !!                 !  1991-11  (G. Madec) 
     8   !!                 !  1992-01  (M. Imbard) insert time step initialization 
     9   !!                 !  1996-06  (G. Madec) generalized vertical coordinate  
     10   !!                 !  1997-02  (G. Madec) creation of domwri.F 
     11   !!                 !  2001-05  (E.Durand - G. Madec) insert closed sea 
     12   !!   NEMO     1.0  !  2002-08  (G. Madec)  F90: Free form and module 
     13   !!            2.0  !  2005-11  (V. Garnier) Surface pressure gradient organization 
     14   !!---------------------------------------------------------------------- 
     15    
    716   !!---------------------------------------------------------------------- 
    817   !!   dom_init       : initialize the space and time domain 
     
    1019   !!   dom_ctl        : control print for the ocean domain 
    1120   !!---------------------------------------------------------------------- 
    12    !! * Modules used 
    1321   USE oce             !  
    1422   USE dom_oce         ! ocean space and time domain 
     
    3038   PRIVATE 
    3139 
    32    !! * Routine accessibility 
    33    PUBLIC dom_init       ! called by opa.F90 
     40   PUBLIC   dom_init   ! called by opa.F90 
    3441 
    3542   !! * Substitutions 
    3643#  include "domzgr_substitute.h90" 
    37    !!---------------------------------------------------------------------- 
    38    !!   OPA 9.0 , LOCEAN-IPSL (2005)  
     44   !!------------------------------------------------------------------------- 
     45   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)  
    3946   !! $Id$ 
    40    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    41    !!---------------------------------------------------------------------- 
     47   !! Software is governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     48   !!------------------------------------------------------------------------- 
    4249 
    4350CONTAINS 
     
    5865      !!      - dom_stp: defined the model time step 
    5966      !!      - dom_wri: create the meshmask file if nmsh=1 
    60       !! 
    61       !! History : 
    62       !!        !  90-10  (C. Levy - G. Madec)  Original code 
    63       !!        !  91-11  (G. Madec) 
    64       !!        !  92-01  (M. Imbard) insert time step initialization 
    65       !!        !  96-06  (G. Madec) generalized vertical coordinate  
    66       !!        !  97-02  (G. Madec) creation of domwri.F 
    67       !!        !  01-05  (E.Durand - G. Madec) insert closed sea 
    68       !!   8.5  !  02-08  (G. Madec)  F90: Free form and module 
    69       !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization 
    70       !!---------------------------------------------------------------------- 
    71       !! * Local declarations 
     67      !!---------------------------------------------------------------------- 
    7268      INTEGER ::   jk                ! dummy loop argument 
    7369      INTEGER ::   iconf = 0         ! temporary integers 
     
    9086      CALL dom_msk                        ! Masks 
    9187 
    92       IF( lk_vvl )   CALL dom_vvl_ini     ! Vertical variable mesh 
     88      IF( lk_vvl )   CALL dom_vvl         ! Vertical variable mesh 
    9389 
    9490      ! Local depth or Inverse of the local depth of the water column at u- and v-points 
    9591      ! ------------------------------ 
    9692      ! Ocean depth at U- and V-points 
    97       hu(:,:) = 0. 
    98       hv(:,:) = 0. 
    99  
     93      hu(:,:) = 0.e0 
     94      hv(:,:) = 0.e0 
    10095      DO jk = 1, jpk 
    10196         hu(:,:) = hu(:,:) + fse3u(:,:,jk) * umask(:,:,jk) 
     
    105100      hur(:,:) = fse3u(:,:,1)             ! Lower bound : thickness of the first model level 
    106101      hvr(:,:) = fse3v(:,:,1) 
    107  
    108102      DO jk = 2, jpk                      ! Sum of the vertical scale factors 
    109103         hur(:,:) = hur(:,:) + fse3u(:,:,jk) * umask(:,:,jk) 
    110104         hvr(:,:) = hvr(:,:) + fse3v(:,:,jk) * vmask(:,:,jk) 
    111105      END DO 
    112  
    113106      ! Compute and mask the inverse of the local depth 
    114107      hur(:,:) = 1. / hur(:,:) * umask(:,:,1) 
    115108      hvr(:,:) = 1. / hvr(:,:) * vmask(:,:,1) 
    116109 
    117  
    118110      CALL dom_stp                        ! Time step 
    119111 
     
    121113 
    122114      IF( .NOT.ln_rstart )   CALL dom_ctl    ! Domain control 
    123  
     115      ! 
    124116   END SUBROUTINE dom_init 
    125117 
     
    134126      !!              - namdom namelist 
    135127      !!              - namcla namelist 
    136       !! 
    137       !! History : 
    138       !!   9.0  !  03-08  (G. Madec)  Original code 
    139       !!---------------------------------------------------------------------- 
    140       !! * Modules used 
     128      !!---------------------------------------------------------------------- 
    141129      USE ioipsl 
    142130      NAMELIST/namrun/ no    , cexper, cn_ocerst_in, cn_ocerst_out, ln_rstart, nrstdt,   & 
     
    156144      ENDIF 
    157145 
    158       ! Namelist namrun : parameters of the run 
    159       REWIND( numnam ) 
     146      REWIND( numnam )              ! Namelist namrun : parameters of the run 
    160147      READ  ( numnam, namrun ) 
    161  
    162148      IF(lwp) THEN 
    163149         WRITE(numout,*) '        Namelist namrun' 
     
    228214      ENDIF 
    229215 
    230       ! Namelist namdom : space/time domain (bathymetry, mesh, timestep) 
    231       REWIND( numnam ) 
     216      REWIND( numnam )              ! Namelist namdom : space/time domain (bathymetry, mesh, timestep) 
    232217      READ  ( numnam, namdom ) 
    233218 
     
    252237      ENDIF 
    253238 
    254       ! Default values 
    255239      n_cla = 0 
    256  
    257       ! Namelist cross land advection 
    258       REWIND( numnam ) 
     240      REWIND( numnam )              ! Namelist cross land advection 
    259241      READ  ( numnam, namcla ) 
    260242      IF(lwp) THEN 
     
    264246      ENDIF 
    265247 
    266       IF( nbit_cmp == 1 .AND. n_cla /= 0 ) THEN 
    267          CALL ctl_stop( ' Reproductibility tests (nbit_cmp=1) require n_cla = 0' ) 
    268       END IF 
    269  
     248      IF( nbit_cmp == 1 .AND. n_cla /= 0 )   CALL ctl_stop( ' Reproductibility tests (nbit_cmp=1) require n_cla = 0' ) 
     249      ! 
    270250   END SUBROUTINE dom_nam 
    271251 
     
    278258      !! 
    279259      !! ** Method  :   compute and print extrema of masked scale factors 
    280       !! 
    281       !! History : 
    282       !!   8.5  !  02-08  (G. Madec)    Original code 
    283       !!---------------------------------------------------------------------- 
    284       !! * Local declarations 
     260      !!---------------------------------------------------------------------- 
    285261      INTEGER ::   iimi1, ijmi1, iimi2, ijmi2, iima1, ijma1, iima2, ijma2 
    286262      INTEGER, DIMENSION(2) ::   iloc      !  
    287263      REAL(wp) ::   ze1min, ze1max, ze2min, ze2max 
    288264      !!---------------------------------------------------------------------- 
    289  
    290       ! Extrema of the scale factors 
    291265 
    292266      IF(lwp)WRITE(numout,*) 
     
    325299         WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, iimi2, ijmi2 
    326300      ENDIF 
    327  
     301      ! 
    328302   END SUBROUTINE dom_ctl 
    329303 
  • branches/dev_004_VVL/NEMO/OPA_SRC/DOM/domvvl.F90

    r1389 r1408  
    44   !! Ocean :  
    55   !!====================================================================== 
    6    !! History :   9.0  !  06-06  (B. Levier, L. Marie)  original code 
    7    !!              "   !  07-07  (D. Storkey) Bug fixes and code for BDY option. 
     6   !! History :  2.0  !  2006-06  (B. Levier, L. Marie)  original code 
     7   !!            3.1  !  2009-02  (G. Madec, M. Leclair, R. Benshila)  pure z* coordinate 
    88   !!---------------------------------------------------------------------- 
    9  
     9#if defined key_vvl 
    1010   !!---------------------------------------------------------------------- 
    1111   !!   'key_vvl'                              variable volume 
    1212   !!---------------------------------------------------------------------- 
     13   !!   dom_vvl     : defined coefficients to distribute ssh on each layers 
    1314   !!---------------------------------------------------------------------- 
    14    !!   dom_vvl         : empty routine 
    15    !!   dom_vvl_ini     : defined coefficients to distribute ssh on each layers 
    16    !!---------------------------------------------------------------------- 
    17    !! * Modules used 
    1815   USE oce             ! ocean dynamics and tracers 
    1916   USE dom_oce         ! ocean space and time domain 
    2017   USE sbc_oce         ! surface boundary condition: ocean 
    21    USE dynspg_oce      ! surface pressure gradient variables 
    2218   USE phycst          ! physical constants 
    2319   USE in_out_manager  ! I/O manager 
    2420   USE lib_mpp         ! distributed memory computing library 
    2521   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    26    USE bdy_oce         ! unstructured open boundary conditions 
    2722 
    2823   IMPLICIT NONE 
    2924   PRIVATE 
    3025 
    31    !! * Routine accessibility 
    32    PUBLIC dom_vvl_ini    ! called by dom_init.F90 
    33    PUBLIC dom_vvl        ! called by istate.F90 and step.F90 
     26   PUBLIC dom_vvl        ! called by domain.F90 
    3427 
    35    !! * Module variables 
    36    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   hu_0, hv_0 
    37    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ee_t, ee_u, ee_v, ee_f   !: 
     28   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ee_t, ee_u, ee_v, ee_f   !: ??? 
    3829 
    39    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   &  !: 
    40       mut, muu, muv, muf                            !: 
     30   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   mut, muu, muv, muf   !: ???  
    4131 
    4232   REAL(wp), DIMENSION(jpk) ::   r2dt               ! vertical profile time-step, = 2 rdttra  
     
    4737#  include "vectopt_loop_substitute.h90" 
    4838   !!---------------------------------------------------------------------- 
    49    !!   OPA 9.0 , LOCEAN-IPSL (2005)  
     39   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)  
    5040   !! $Id$ 
    5141   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     
    5444CONTAINS        
    5545 
    56 #if defined key_vvl 
    57  
    58    SUBROUTINE dom_vvl_ini 
     46   SUBROUTINE dom_vvl 
    5947      !!---------------------------------------------------------------------- 
    60       !!                ***  ROUTINE dom_vvl_ini  *** 
     48      !!                ***  ROUTINE dom_vvl  *** 
    6149      !!                    
    6250      !! ** Purpose :  compute coefficients muX at T-U-V-F points to spread 
     
    7058      IF(lwp)   THEN 
    7159         WRITE(numout,*) 
    72          WRITE(numout,*) 'dom_vvl_ini : Variable volume activated' 
    73          WRITE(numout,*) '~~~~~~~~~~~   compute coef. used to spread ssh over each layers' 
     60         WRITE(numout,*) 'dom_vvl : Variable volume activated' 
     61         WRITE(numout,*) '~~~~~~~~  compute coef. used to spread ssh over each layers' 
    7462      ENDIF 
    7563 
     
    127115 
    128116 
    129 !!debug print 
    130 !         ii=50   ;   ij = 50 
    131 !         do jk= 1, jpk  
    132 !         WRITE(numout,*) 'domvvl GM  : h0=', SUM( fse3t_0(ii,ij,1:jk) * tmask(ii,ij,1:jk) ), 'e3t0=', fse3t_0(ii,ij,jk),   & 
    133 !            &            'e3t =', fse3t_0(ii,ij,jk) * ( 1 +  mut(ii,ij,jk) ), 'mut', mut(ii,ij,jk),   & 
    134 !            &            'h =', SUM( fse3t_0(ii,ij,1:jk) * ( 1 +  mut(ii,ij,1:jk) ) * tmask(ii,ij,1:jk) ) 
    135 !         end do 
    136 !!end debug print 
    137  
    138117      ! Reference ocean depth at U- and V-points 
    139118      hu_0(:,:) = 0.e0     
    140119      hv_0(:,:) = 0.e0 
    141120      DO jk = 1, jpk 
    142          hu_0(:,:) = hu_0(:,:) + fse3u(:,:,jk) * umask(:,:,jk) 
    143          hv_0(:,:) = hv_0(:,:) + fse3v(:,:,jk) * vmask(:,:,jk) 
     121         hu_0(:,:) = hu_0(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk) 
     122         hv_0(:,:) = hv_0(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk) 
    144123      END DO 
    145124 
     
    149128            zcoefu = 0.5  * umask(ji,jj,1) / ( e1u(ji,jj) * e2u(ji,jj) ) 
    150129            zcoefv = 0.5  * vmask(ji,jj,1) / ( e1v(ji,jj) * e2v(ji,jj) ) 
    151             zcoeff = 0.25 * fmask(ji,jj,1) 
    152 !!gm bug used of fmask, even if thereafter multiplied by muf  which is correctly masked) 
     130            zcoeff = 0.25 * umask(ji,jj,1) * umask(ji,jj+1,1) 
    153131            sshu_b(ji,jj) = zcoefu * ( e1t(ji  ,jj) * e2t(ji  ,jj) * sshb(ji  ,jj)     & 
    154132               &                     + e1t(ji+1,jj) * e2t(ji+1,jj) * sshb(ji+1,jj) )    
     
    170148      CALL lbc_lnk( sshf_b, 'F', 1. )   ;   CALL lbc_lnk( sshf_n, 'F', 1. ) 
    171149      ! 
    172  
    173  
    174    END SUBROUTINE dom_vvl_ini 
    175  
    176  
    177    SUBROUTINE dom_vvl 
    178       !!---------------------------------------------------------------------- 
    179       !!                ***  ROUTINE dom_vvl  *** 
    180       !!                    
    181       !! ** Purpose :  compute ssh at U-V-F points, T-W scale factors and local 
    182       !!               depths at each time step. 
    183       !!---------------------------------------------------------------------- 
    184       !! * Local declarations 
    185       INTEGER  ::   ji, jj, jk                 ! dummy loop indices 
    186       !!---------------------------------------------------------------------- 
    187  
    188  !     IF( kt == nit000 ) THEN 
    189  !        IF(lwp) WRITE(numout,*) 
    190  !        IF(lwp) WRITE(numout,*) 'dom_vvl : ' 
    191  !        IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
    192  !     ENDIF 
    193  
    194        WRITE(*,*) 'dom_vvl  : empty routine, you should not be here' 
    195  
    196150   END SUBROUTINE dom_vvl 
    197151 
     
    200154   !!   Default option :                                      Empty routine 
    201155   !!---------------------------------------------------------------------- 
    202    SUBROUTINE dom_vvl_ini 
    203    END SUBROUTINE dom_vvl_ini 
     156CONTAINS 
    204157   SUBROUTINE dom_vvl 
    205158   END SUBROUTINE dom_vvl 
  • branches/dev_004_VVL/NEMO/OPA_SRC/DOM/domzgr_substitute.h90

    r1385 r1408  
    7878#   define  fse3uw_a(i,j,k)  (fse3uw_0(i,j,k)*(1+sshu_a(i,j)*muu(i,j,k))) 
    7979#   define  fse3vw_a(i,j,k)  (fse3vw_0(i,j,k)*(1+sshv_a(i,j)*muv(i,j,k))) 
     80 
    8081#else 
    8182! z- or s-coordinate (1D or 3D + no time dependency) use reference in all cases 
  • branches/dev_004_VVL/NEMO/OPA_SRC/DYN/dynnxt.F90

    r1390 r1408  
    1515   !!            2.0  !  2005-11  (V. Garnier) Surface pressure gradient organization 
    1616   !!            2.3  !  2007-07  (D. Storkey) Calls to BDY routines.  
    17    !!            3.1  !  2009-02  (G. Madec)  re-introduce the vvl option 
     17   !!            3.2  !  2009-04  (G. Madec, R.Benshila))  re-introduce the vvl option 
    1818   !!---------------------------------------------------------------------- 
    1919   
     
    2121   !!   dyn_nxt      : update the horizontal velocity from the momentum trend 
    2222   !!---------------------------------------------------------------------- 
    23    !! * Modules used 
    2423   USE oce             ! ocean dynamics and tracers 
    2524   USE dom_oce         ! ocean space and time domain 
     
    4241   PRIVATE 
    4342 
    44    !! * Accessibility 
    45    PUBLIC dyn_nxt                ! routine called by step.F90 
     43   PUBLIC    dyn_nxt   ! routine called by step.F90 
     44 
    4645   !! * Substitutions 
    4746#  include "domzgr_substitute.h90" 
    48    !!---------------------------------------------------------------------- 
    49    !!  OPA 9.0 , LOCEAN-IPSL (2005)  
     47   !!------------------------------------------------------------------------- 
     48   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)  
    5049   !! $Id$  
    51    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    52    !!---------------------------------------------------------------------- 
     50   !! Software is governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)  
     51   !!------------------------------------------------------------------------- 
    5352 
    5453CONTAINS 
     
    8584      REAL(wp) ::   ze3u_b, ze3u_n, ze3u_a   !    -         - 
    8685      REAL(wp) ::   ze3v_b, ze3v_n, ze3v_a   !    -         -  
    87       !!---------------------------------------------------------------------- 
     86      REAL(wp) ::   zuf   , zvf              !    -         -  
     87     !!---------------------------------------------------------------------- 
    8888 
    8989      IF( kt == nit000 ) THEN 
     
    100100 
    101101      ! Lateral boundary conditions on ( ua, va ) 
    102       CALL lbc_lnk( ua, 'U', -1. ) 
    103       CALL lbc_lnk( va, 'V', -1. ) 
     102      CALL lbc_lnk( ua, 'U', -1. )   ;   CALL lbc_lnk( va, 'V', -1. )  
    104103 
    105104      ! Next velocity 
     
    110109      IF( lk_vvl ) THEN                                  ! Varying levels 
    111110         DO jk = 1, jpkm1 
    112             ua(:,:,jk) = (        ub(:,:,jk) * fse3u_b(:,:,jk)     & 
     111            ua(:,:,jk) = (          ub(:,:,jk) * fse3u_b(:,:,jk)     & 
    113112               &           + z2dt * ua(:,:,jk) * fse3u_n(:,:,jk) )   & 
    114113               &         / fse3u_a(:,:,jk) * umask(:,:,jk) 
    115             va(:,:,jk) = (        vb(:,:,jk) * fse3v_b(:,:,jk)     & 
     114            va(:,:,jk) = (          vb(:,:,jk) * fse3v_b(:,:,jk)     & 
    116115               &           + z2dt * va(:,:,jk) * fse3v_n(:,:,jk) )   & 
    117116               &         / fse3v_a(:,:,jk) * vmask(:,:,jk) 
     
    130129 
    131130      IF ( lk_dynspg_exp .OR. lk_dynspg_ts ) THEN 
    132          !Flather boundary condition : 
     131         ! Flather boundary condition : 
    133132         !        - Update sea surface height on each open boundary 
    134133         !                 sshn (= after ssh) for explicit case 
     
    136135         !        - Correct the barotropic velocities 
    137136         CALL obc_dyn_bt( kt ) 
    138  
    139          !Boundary conditions on sshn ( after ssh) 
     137         ! 
     138         ! Boundary conditions on sshn ( after ssh) 
    140139         CALL lbc_lnk( sshn, 'T', 1. ) 
    141  
    142          IF(ln_ctl) THEN         ! print sum trends (used for debugging) 
    143             CALL prt_ctl(tab2d_1=sshn, clinfo1=' ssh      : ', mask1=tmask) 
    144          ENDIF 
    145  
    146          IF ( ln_vol_cst ) CALL obc_vol( kt ) 
    147  
     140         ! 
     141         IF( ln_vol_cst )   CALL obc_vol( kt ) 
     142         ! 
     143         IF(ln_ctl)   CALL prt_ctl( tab2d_1=sshn, clinfo1=' ssh      : ', mask1=tmask ) 
    148144      ENDIF 
    149145 
    150146# elif defined key_bdy  
    151147      ! Update (ua,va) along open boundaries (for exp or ts options). 
    152       IF ( lk_dynspg_exp .or. lk_dynspg_ts ) THEN 
    153    
     148      IF( lk_dynspg_exp .OR. lk_dynspg_ts ) THEN 
     149         ! 
    154150         CALL bdy_dyn_frs( kt ) 
    155  
    156          IF ( ln_bdy_fla ) THEN 
    157  
    158             ua_e(:,:)=0.0 
    159             va_e(:,:)=0.0 
    160  
     151         ! 
     152         IF( ln_bdy_fla ) THEN 
     153            ua_e(:,:) = 0.e0 
     154            va_e(:,:) = 0.e0 
    161155            ! Set these variables for use in bdy_dyn_fla 
    162156            hu_e(:,:) = hu(:,:) 
    163157            hv_e(:,:) = hv(:,:) 
    164  
    165             DO jk = 1, jpkm1 
    166                !! Vertically integrated momentum trends 
     158            DO jk = 1, jpkm1   !! Vertically integrated momentum trends 
    167159               ua_e(:,:) = ua_e(:,:) + fse3u(:,:,jk) * umask(:,:,jk) * ua(:,:,jk) 
    168160               va_e(:,:) = va_e(:,:) + fse3v(:,:,jk) * vmask(:,:,jk) * va(:,:,jk) 
    169161            END DO 
    170  
    171162            DO jk = 1 , jpkm1 
    172163               ua(:,:,jk) = ua(:,:,jk) - ua_e(:,:) * hur(:,:) 
    173164               va(:,:,jk) = va(:,:,jk) - va_e(:,:) * hvr(:,:) 
    174165            END DO 
    175  
    176166            CALL bdy_dta_bt( kt+1, 0) 
    177167            CALL bdy_dyn_fla 
    178  
    179168         ENDIF 
    180  
     169         ! 
    181170         DO jk = 1 , jpkm1 
    182171            ua(:,:,jk) = ua(:,:,jk) + ua_e(:,:) * hur(:,:) 
    183172            va(:,:,jk) = va(:,:,jk) + va_e(:,:) * hvr(:,:) 
    184173         END DO 
    185  
    186       ENDIF 
    187  
     174      ENDIF 
    188175# endif 
     176 
    189177# if defined key_agrif 
    190178      CALL Agrif_dyn( kt ) 
     
    196184      IF( neuler == 0 .AND. kt == nit000 ) THEN 
    197185         DO jk = 1, jpkm1                                  
    198             ub(:,:,jk) = un(:,:,jk)  
    199             vb(:,:,jk) = vn(:,:,jk) 
    200186            un(:,:,jk) = ua(:,:,jk) 
    201187            vn(:,:,jk) = va(:,:,jk) 
     
    230216            END DO 
    231217         ELSE                                             ! Fixed levels 
    232 !RB_vvl : should be done as in tranxt ? 
    233218            DO jk = 1, jpkm1                              ! filter applied on velocities 
    234                ub(:,:,jk) = atfp * ( ub(:,:,jk) + ua(:,:,jk) ) + atfp1 * un(:,:,jk) 
    235                vb(:,:,jk) = atfp * ( vb(:,:,jk) + va(:,:,jk) ) + atfp1 * vn(:,:,jk) 
    236                un(:,:,jk) = ua(:,:,jk) 
    237                vn(:,:,jk) = va(:,:,jk) 
     219               DO jj = 1, jpj 
     220                  DO ji = 1, jpi 
     221                     zuf = atfp * ( ub(ji,jj,jk) + ua(ji,jj,jk) ) + atfp1 * un(ji,jj,jk) 
     222                     zvf = atfp * ( vb(ji,jj,jk) + va(ji,jj,jk) ) + atfp1 * vn(ji,jj,jk) 
     223                     ub(ji,jj,jk) = zuf 
     224                     vb(ji,jj,jk) = zvf 
     225                     un(ji,jj,jk) = ua(ji,jj,jk) 
     226                     vn(ji,jj,jk) = va(ji,jj,jk) 
     227                  END DO 
     228               END DO 
    238229            END DO 
    239230         ENDIF 
    240       ENDIF 
    241  
    242       IF(ln_ctl)   THEN 
    243          CALL prt_ctl(tab3d_1=un, clinfo1=' nxt  - Un: ', mask1=umask, & 
    244             &         tab3d_2=vn, clinfo2=' Vn: ', mask2=vmask) 
    245231      ENDIF 
    246232 
     
    249235#endif       
    250236 
     237      IF(ln_ctl)   CALL prt_ctl( tab3d_1=un, clinfo1=' nxt  - Un: ', mask1=umask,   & 
     238         &                       tab3d_2=vn, clinfo2=' Vn: '       , mask2=vmask ) 
     239      !  
    251240   END SUBROUTINE dyn_nxt 
    252241 
  • branches/dev_004_VVL/NEMO/OPA_SRC/DYN/dynspg_exp.F90

    r1390 r1408  
    44   !! Ocean dynamics:  surface pressure gradient trend 
    55   !!====================================================================== 
     6   !! History :  9.0  !  2005-11  (V. Garnier, G. Madec, L. Bessieres) Original code 
     7   !!---------------------------------------------------------------------- 
    68#if defined key_dynspg_exp   ||   defined key_esopa 
    79   !!---------------------------------------------------------------------- 
     
    1214   !!                      volume case with vector optimization 
    1315   !!---------------------------------------------------------------------- 
    14    !! * Modules used 
    1516   USE oce             ! ocean dynamics and tracers  
    1617   USE dom_oce         ! ocean space and time domain  
     
    3031   PRIVATE 
    3132 
    32    !! * Accessibility 
    33    PUBLIC dyn_spg_exp  ! routine called by step.F90 
     33   PUBLIC   dyn_spg_exp   ! routine called by step.F90 
    3434 
    3535   !! * Substitutions 
     
    3737#  include "vectopt_loop_substitute.h90" 
    3838   !!---------------------------------------------------------------------- 
    39    !!  OPA 9.0 , LOCEAN-IPSL (2005)  
     39   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)  
    4040   !! $Id$ 
    41    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     41   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)  
    4242   !!---------------------------------------------------------------------- 
    43  
    4443 
    4544CONTAINS 
     
    6261      !! 
    6362      !! ** Action : - Update (ua,va) with the surf. pressure gradient trend 
     63      !!--------------------------------------------------------------------- 
     64      INTEGER, INTENT( in )  ::   kt         ! ocean time-step index 
    6465      !! 
    65       !! References : 
    66       !! 
    67       !! History : 
    68       !!   9.0  !  05-11  (V. Garnier, G. Madec, L. Bessieres) Original code 
    69       !! 
    70       !!--------------------------------------------------------------------- 
    71       !! * Arguments 
    72       INTEGER, INTENT( in )  ::   kt         ! ocean time-step index 
    73  
    74       !! * Local declarations 
    7566      INTEGER  ::   ji, jj, jk               ! dummy loop indices 
    7667      !!---------------------------------------------------------------------- 
     
    8677      ENDIF 
    8778 
    88       ! 0. Initialization 
    89       ! ----------------- 
    9079      ! read or estimate sea surface height and vertically integrated velocities 
    9180      IF( lk_obc )   CALL obc_dta_bt( kt, 0 ) 
    9281 
    93       ! 1. Surface pressure gradient (now) 
    94       ! ---------------------------- 
     82      ! Surface pressure gradient (now) 
    9583      DO jj = 2, jpjm1 
    9684         DO ji = fs_2, fs_jpim1   ! vector opt. 
     
    10088      END DO 
    10189 
    102       ! 2. Add the surface pressure trend to the general trend 
    103       ! ------------------------------------------------------ 
     90      ! Add the surface pressure trend to the general trend 
    10491      DO jk = 1, jpkm1 
    10592         DO jj = 2, jpjm1 
  • branches/dev_004_VVL/NEMO/OPA_SRC/DYN/dynspg_flt.F90

    r1390 r1408  
    44   !! Ocean dynamics:  surface pressure gradient trend 
    55   !!====================================================================== 
    6    !! History    8.0  !  98-05  (G. Roullet)  free surface 
    7    !!                 !  98-10  (G. Madec, M. Imbard)  release 8.2 
    8    !!            8.5  !  02-08  (G. Madec)  F90: Free form and module 
    9    !!            " "  !  02-11  (C. Talandier, A-M Treguier) Open boundaries 
    10    !!            9.0  !  04-08  (C. Talandier) New trends organization 
    11    !!            " "  !  05-11  (V. Garnier) Surface pressure gradient organization 
    12    !!            " "  !  06-07  (S. Masson)  distributed restart using iom 
    13    !!            " "  !  05-01  (J.Chanut, A.Sellar) Calls to BDY routines.  
     6   !! History    OPA  !  1998-05  (G. Roullet)  free surface 
     7   !!                 !  1998-10  (G. Madec, M. Imbard)  release 8.2 
     8   !!   NEMO     O.1  !  2002-08  (G. Madec)  F90: Free form and module 
     9   !!             -   !  2002-11  (C. Talandier, A-M Treguier) Open boundaries 
     10   !!            1.0  !  2004-08  (C. Talandier) New trends organization 
     11   !!             -   !  2005-11  (V. Garnier) Surface pressure gradient organization 
     12   !!            2.0  !  2006-07  (S. Masson)  distributed restart using iom 
     13   !!             -   !  2006-08  (J.Chanut, A.Sellar) Calls to BDY routines.  
     14   !!            3.2  !  2009-03  (G. Madec, M. Leclair, R. Benshila) introduce sshwzv module 
    1415   !!---------------------------------------------------------------------- 
    1516#if defined key_dynspg_flt   ||   defined key_esopa   
    1617   !!---------------------------------------------------------------------- 
    1718   !!   'key_dynspg_flt'                              filtered free surface 
    18    !!---------------------------------------------------------------------- 
    1919   !!---------------------------------------------------------------------- 
    2020   !!   dyn_spg_flt  : update the momentum trend with the surface pressure 
     
    5353   PRIVATE 
    5454 
    55    PUBLIC dyn_spg_flt  ! routine called by step.F90 
    56    PUBLIC flt_rst      ! routine called by istate.F90 
     55   PUBLIC   dyn_spg_flt  ! routine called by step.F90 
     56   PUBLIC   flt_rst      ! routine called by istate.F90 
    5757 
    5858   !! * Substitutions 
     
    6060#  include "vectopt_loop_substitute.h90" 
    6161   !!---------------------------------------------------------------------- 
    62    !!   OPA 9.0 , LOCEAN-IPSL (2005)  
     62   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)  
    6363   !! $Id$ 
    6464   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)  
     
    8080      !!         spgv = 1/rau0 d/dy(ps) =  1/e2v dj( sshn + rnu btda ) 
    8181      !!      where sshn is the free surface elevation and btda is the after 
    82       !!      of the free surface elevation 
     82      !!      time derivative of the free surface elevation 
    8383      !!       -1- evaluate the surface presure trend (including the addi- 
    8484      !!      tional force) in three steps: 
     
    104104      !! References : Roullet and Madec 1999, JGR. 
    105105      !!--------------------------------------------------------------------- 
    106       !! * Modules used 
    107       USE oce    , ONLY :   zub   => ta,             &  ! ta used as workspace 
    108                             zvb   => sa                 ! sa   "          " 
    109  
    110       INTEGER, INTENT( in )  ::   kt         ! ocean time-step index 
    111       INTEGER, INTENT( out ) ::   kindic     ! solver convergence flag (<0 if not converge) 
     106      USE oce, ONLY :   zub   => ta   ! ta used as workspace 
     107      USE oce, ONLY :   zvb   => sa   ! ta used as workspace 
     108      !! 
     109      INTEGER, INTENT( in )  ::   kt       ! ocean time-step index 
     110      INTEGER, INTENT( out ) ::   kindic   ! solver convergence flag (<0 if not converge) 
    112111      !!                                    
    113112      INTEGER  ::   ji, jj, jk                          ! dummy loop indices 
    114       REAL(wp) ::   z2dt, z2dtg, zraur, znugdt,      &  ! temporary scalars 
    115          &          znurau, zgcb, zbtd,              &  !   "          " 
    116          &          ztdgu, ztdgv                        !   "          " 
     113      REAL(wp) ::   z2dt, z2dtg, zraur, znugdt          ! temporary scalars 
     114      REAL(wp) ::   znurau, zgcb, zbtd                  !   "          " 
     115      REAL(wp) ::   ztdgu, ztdgv                        !   "          " 
    117116      !!---------------------------------------------------------------------- 
    118117      ! 
     
    130129         ! when using agrif, sshn, gcx have to be read in istate 
    131130         IF (.NOT. lk_agrif) CALL flt_rst( nit000, 'READ' )       ! read or initialize the following fields: 
    132          !                                                        ! gcx, gcxb, sshb, sshn 
     131         !                                                        ! gcx, gcxb 
    133132      ENDIF 
    134133 
     
    354353         DO jj = 2, jpjm1 
    355354            DO ji = fs_2, fs_jpim1   ! vector opt. 
    356                ua(ji,jj,jk) = (ua(ji,jj,jk) + spgu(ji,jj)) * umask(ji,jj,jk) 
    357                va(ji,jj,jk) = (va(ji,jj,jk) + spgv(ji,jj)) * vmask(ji,jj,jk) 
     355               ua(ji,jj,jk) = ( ua(ji,jj,jk) + spgu(ji,jj) ) * umask(ji,jj,jk) 
     356               va(ji,jj,jk) = ( va(ji,jj,jk) + spgv(ji,jj) ) * vmask(ji,jj,jk) 
    358357            END DO 
    359358         END DO 
     
    363362      ! -------------------------------------------------- 
    364363      IF( lrst_oce ) CALL flt_rst( kt, 'WRITE' ) 
    365  
    366       ! print sum trends (used for debugging) 
    367       IF(ln_ctl)     CALL prt_ctl( tab2d_1=sshn, clinfo1=' spg  - ssh: ', mask1=tmask ) 
    368364      ! 
    369365   END SUBROUTINE dyn_spg_flt 
     
    386382           CALL iom_get( numror, jpdom_autoglo, 'gcx' , gcx (1:jpi,1:jpj) ) 
    387383           CALL iom_get( numror, jpdom_autoglo, 'gcxb', gcxb(1:jpi,1:jpj) ) 
    388            IF( neuler == 0 ) THEN 
    389               gcxb(:,:) = gcx (:,:) 
    390            ENDIF 
     384           IF( neuler == 0 )   gcxb(:,:) = gcx (:,:) 
    391385        ELSE 
    392386           gcx (:,:) = 0.e0 
     
    396390! Caution : extra-hallow 
    397391! gcx and gcxb are defined as: DIMENSION(1-jpr2di:jpi+jpr2di,1-jpr2dj:jpj+jpr2dj) 
    398         CALL iom_rstput( kt, nitrst, numrow, 'gcx' , gcx( 1:jpi,1:jpj) ) 
     392        CALL iom_rstput( kt, nitrst, numrow, 'gcx' , gcx (1:jpi,1:jpj) ) 
    399393        CALL iom_rstput( kt, nitrst, numrow, 'gcxb', gcxb(1:jpi,1:jpj) ) 
    400394     ENDIF 
  • branches/dev_004_VVL/NEMO/OPA_SRC/DYN/dynspg_ts.F90

    r1405 r1408  
    11MODULE dynspg_ts 
    22   !!====================================================================== 
    3    !! History :   9.0  !  04-12  (L. Bessieres, G. Madec)  Original code 
    4    !!             " "  !  05-11  (V. Garnier, G. Madec)  optimization 
    5    !!             9.0  !  06-08  (S. Masson)  distributed restart using iom 
    6    !!              "   !  08-01  (R. Benshila)  change averaging method 
    7    !!              "   !  07-07  (D. Storkey) calls to BDY routines 
    8    !!--------------------------------------------------------------------- 
     3   !! History :   1.0  !  04-12  (L. Bessieres, G. Madec)  Original code 
     4   !!              -   !  05-11  (V. Garnier, G. Madec)  optimization 
     5   !!              -   !  06-08  (S. Masson)  distributed restart using iom 
     6   !!             2.0  !  07-07  (D. Storkey) calls to BDY routines 
     7   !!              -   !  08-01  (R. Benshila)  change averaging method 
     8  !!--------------------------------------------------------------------- 
    99#if defined key_dynspg_ts   ||   defined key_esopa 
    1010   !!---------------------------------------------------------------------- 
     
    1616   !!   ts_rst      : read/write the time-splitting restart fields in the ocean restart file 
    1717   !!---------------------------------------------------------------------- 
    18    !! * Modules used 
    1918   USE oce             ! ocean dynamics and tracers 
    2019   USE dom_oce         ! ocean space and time domain 
     
    4645   PUBLIC ts_rst      ! routine called by istate.F90 
    4746 
    48    REAL(wp), DIMENSION(jpi,jpj) ::  ftnw, ftne,   &  ! triad of coriolis parameter 
    49       &                             ftsw, ftse       ! (only used with een vorticity scheme) 
    50  
     47   REAL(wp), DIMENSION(jpi,jpj) ::  ftnw, ftne   ! triad of coriolis parameter 
     48   REAL(wp), DIMENSION(jpi,jpj) ::  ftsw, ftse   ! (only used with een vorticity scheme) 
    5149 
    5250   !! * Substitutions 
    5351#  include "domzgr_substitute.h90" 
    5452#  include "vectopt_loop_substitute.h90" 
    55    !!---------------------------------------------------------------------- 
    56    !!  OPA 9.0 , LOCEAN-IPSL (2005)  
     53   !!------------------------------------------------------------------------- 
     54   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)  
    5755   !! $Id$ 
    58    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    59    !!---------------------------------------------------------------------- 
     56   !! Software is governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)  
     57   !!------------------------------------------------------------------------- 
    6058 
    6159CONTAINS 
     
    9088      !!--------------------------------------------------------------------- 
    9189      INTEGER, INTENT( in )  ::   kt           ! ocean time-step index 
    92  
    93       !! * Local declarations 
     90      !! 
    9491      INTEGER  ::  ji, jj, jk, jit             ! dummy loop indices 
    9592      INTEGER  ::  icycle                      ! temporary scalar 
  • branches/dev_004_VVL/NEMO/OPA_SRC/DYN/dynvor.F90

    r1383 r1408  
    55   !!                 planetary vorticity trends 
    66   !!====================================================================== 
    7    !! History :  1.0  !  89-12  (P. Andrich)  vor_ens: Original code 
    8    !!            5.0  !  91-11  (G. Madec) vor_ene, vor_mix: Original code 
    9    !!            6.0  !  96-01  (G. Madec)  s-coord, suppress work arrays 
    10    !!            8.5  !  02-08  (G. Madec)  F90: Free form and module 
    11    !!            8.5  !  04-02  (G. Madec)  vor_een: Original code 
    12    !!            9.0  !  03-08  (G. Madec)  vor_ctl: Original code 
    13    !!            9.0  !  05-11  (G. Madec)  dyn_vor: Original code (new step architecture) 
    14    !!            9.0  !  06-11  (G. Madec)  flux form advection: add metric term 
     7   !! History :  OPA  !  1989-12  (P. Andrich)  vor_ens: Original code 
     8   !!            5.0  !  1991-11  (G. Madec) vor_ene, vor_mix: Original code 
     9   !!            6.0  !  1996-01  (G. Madec)  s-coord, suppress work arrays 
     10   !!            8.5  !  2002-08  (G. Madec)  F90: Free form and module 
     11   !!   NEMO     1.0  !  2004-02  (G. Madec)  vor_een: Original code 
     12   !!             -   !  2003-08  (G. Madec)  add vor_ctl 
     13   !!             -   !  2005-11  (G. Madec)  add dyn_vor (new step architecture) 
     14   !!            2.0  !  2006-11  (G. Madec)  flux form advection: add metric term 
     15   !!            3.2  !  2009-04  (R. Benshila)  vvl: correction of een scheme 
    1516   !!---------------------------------------------------------------------- 
    1617 
     
    3738   PUBLIC   dyn_vor   ! routine called by step.F90 
    3839 
    39    !!* Namelist nam_dynvor: vorticity term 
     40   !                                             !!* Namelist nam_dynvor: vorticity term 
    4041   LOGICAL, PUBLIC ::   ln_dynvor_ene = .FALSE.   !: energy conserving scheme 
    4142   LOGICAL, PUBLIC ::   ln_dynvor_ens = .TRUE.    !: enstrophy conserving scheme 
     
    5253#  include "vectopt_loop_substitute.h90" 
    5354   !!---------------------------------------------------------------------- 
    54    !!   OPA 9.0 , LOCEAN-IPSL (2006)  
     55   !! NEMO/OPA 3,2 , LOCEAN-IPSL (2009)  
    5556   !! $Id$ 
    5657   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     
    174175      IF(ln_ctl) CALL prt_ctl( tab3d_1=ua, clinfo1=' vor  - Ua: ', mask1=umask, & 
    175176         &                     tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' ) 
    176  
     177      ! 
    177178   END SUBROUTINE dyn_vor 
    178179 
     
    206207      INTEGER , INTENT(in   )                         ::   kt     ! ocean time-step index 
    207208      INTEGER , INTENT(in   )                         ::   kvor   ! =ncor (planetary) ; =ntot (total) ; 
    208          !                                                        ! =nrvm (relative vorticity or metric) 
     209      !                                                           ! =nrvm (relative vorticity or metric) 
    209210      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pua    ! total u-trend 
    210211      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pva    ! total v-trend 
     
    222223      ENDIF 
    223224 
    224       ! Local constant initialization 
    225       zfact2 = 0.5 * 0.5 
     225      zfact2 = 0.5 * 0.5      ! Local constant initialization 
    226226 
    227227!CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz ) 
     
    229229      DO jk = 1, jpkm1                                 ! Horizontal slab 
    230230         !                                             ! =============== 
     231         ! 
    231232         ! Potential vorticity and horizontal fluxes 
    232233         ! ----------------------------------------- 
     
    315316      INTEGER, INTENT(in) ::   kt   ! ocean timestep index 
    316317      !! 
    317       INTEGER ::   ji, jj, jk   ! dummy loop indices 
     318      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    318319      REAL(wp) ::   zfact1, zua, zcua, zx1, zy1   ! temporary scalars 
    319320      REAL(wp) ::   zfact2, zva, zcva, zx2, zy2   !    "         " 
     
    327328      ENDIF 
    328329 
    329       ! Local constant initialization 
    330       zfact1 = 0.5 * 0.25 
     330      zfact1 = 0.5 * 0.25      ! Local constant initialization 
    331331      zfact2 = 0.5 * 0.5 
    332332 
     
    335335      DO jk = 1, jpkm1                                 ! Horizontal slab 
    336336         !                                             ! =============== 
    337  
     337         ! 
    338338         ! Relative and planetary potential vorticity and horizontal fluxes 
    339339         ! ---------------------------------------------------------------- 
     
    438438      ENDIF 
    439439 
    440       ! Local constant initialization 
    441       zfact1 = 0.5 * 0.25 
     440      zfact1 = 0.5 * 0.25      ! Local constant initialization 
    442441 
    443442!CDIR PARALLEL DO PRIVATE( zwx, zwy, zwz ) 
     
    445444      DO jk = 1, jpkm1                                 ! Horizontal slab 
    446445         !                                             ! =============== 
     446         ! 
    447447         ! Potential vorticity and horizontal fluxes 
    448448         ! ----------------------------------------- 
     
    465465                       &       + (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    466466                       &           - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    467                        &       * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )                                               & 
     467                       &       * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )                                                & 
    468468                       &       ) 
    469469               END DO 
    470470            END DO 
    471471         END SELECT 
    472  
     472         ! 
    473473         IF( ln_sco ) THEN 
    474474            DO jj = 1, jpj                      ! caution: don't use (:,:) for this loop  
     
    487487            END DO 
    488488         ENDIF 
    489  
     489         ! 
    490490         ! Compute and add the vorticity term trend 
    491491         ! ---------------------------------------- 
     
    514514      !! 
    515515      !! ** Method  :   Trend evaluated using now fields (centered in time)  
    516       !!      and the Arakawa and Lamb (19XX) flux form formulation : conserves  
     516      !!      and the Arakawa and Lamb (1980) flux form formulation : conserves  
    517517      !!      both the horizontal kinetic energy and the potential enstrophy 
    518       !!      when horizontal divergence is zero. 
    519       !!      The trend of the vorticity term is given by: 
    520       !!       * s-coordinate (ln_sco=T), the e3. are inside the derivatives: 
    521       !!       * z-coordinate (default key), e3t=e3u=e3v, the trend becomes: 
    522       !!      Add this trend to the general momentum trend (ua,va): 
    523       !!          (ua,va) = (ua,va) + ( voru , vorv ) 
     518      !!      when horizontal divergence is zero (see the NEMO documentation) 
     519      !!      Add this trend to the general momentum trend (ua,va). 
    524520      !! 
    525521      !! ** Action : - Update (ua,va) with the now vorticity term trend 
     
    531527      INTEGER , INTENT(in   )                         ::   kt     ! ocean time-step index 
    532528      INTEGER , INTENT(in   )                         ::   kvor   ! =ncor (planetary) ; =ntot (total) ; 
    533          !                                                        ! =nrvm (relative vorticity or metric) 
     529      !                                                           ! =nrvm (relative vorticity or metric) 
    534530      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pua    ! total u-trend 
    535531      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pva    ! total v-trend 
    536532      !! 
    537       INTEGER ::   ji, jj, jk          ! dummy loop indices 
     533      INTEGER  ::   ji, jj, jk         ! dummy loop indices 
    538534      REAL(wp) ::   zfac12, zua, zva   ! temporary scalars 
    539535      REAL(wp), DIMENSION(jpi,jpj) ::   zwx, zwy, zwz            ! temporary 2D workspace 
    540536      REAL(wp), DIMENSION(jpi,jpj) ::   ztnw, ztne, ztsw, ztse   ! temporary 3D workspace 
     537#if defined key_vvl 
     538      REAL(wp), DIMENSION(jpi,jpj,jpk)       ::   ze3f 
     539#else 
    541540      REAL(wp), DIMENSION(jpi,jpj,jpk), SAVE ::   ze3f 
     541#endif 
    542542      !!---------------------------------------------------------------------- 
    543543 
     
    548548      ENDIF 
    549549 
    550       IF( kt == nit000 .OR. lk_vvl ) THEN 
     550      IF( kt == nit000 .OR. lk_vvl ) THEN      ! reciprocal of e3 at F-point (masked averaging of e3t) 
    551551         DO jk = 1, jpk 
    552552            DO jj = 1, jpjm1 
     
    561561      ENDIF 
    562562 
    563       ! Local constant initialization 
    564       zfac12 = 1.e0 / 12.e0 
     563      zfac12 = 1.e0 / 12.e0      ! Local constant initialization 
    565564 
    566565       
     
    573572         ! ----------------------------------------- 
    574573         SELECT CASE( kvor )      ! vorticity considered 
    575          CASE ( 1 )   ;   zwz(:,:) = ff(:,:)      * ze3f(:,:,jk)   ! planetary vorticity (Coriolis) 
    576          CASE ( 2 )   ;   zwz(:,:) = rotn(:,:,jk) * ze3f(:,:,jk)   ! relative  vorticity 
     574         CASE ( 1 )                                                ! planetary vorticity (Coriolis) 
     575            zwz(:,:) = ff(:,:)      * ze3f(:,:,jk) 
     576         CASE ( 2 )                                                ! relative  vorticity 
     577            zwz(:,:) = rotn(:,:,jk) * ze3f(:,:,jk) 
    577578         CASE ( 3 )                                                ! metric term 
    578579            DO jj = 1, jpjm1 
     
    583584               END DO 
    584585            END DO 
    585          CASE ( 4 )   ;   zwz(:,:) = ( rotn(:,:,jk) + ff(:,:) ) * ze3f(:,:,jk) ! total (relative + planetary vorticity) 
     586         CASE ( 4 )                                                ! total (relative + planetary vorticity) 
     587            zwz(:,:) = ( rotn(:,:,jk) + ff(:,:) ) * ze3f(:,:,jk) 
    586588         CASE ( 5 )                                                ! total (coriolis + metric) 
    587589            DO jj = 1, jpjm1 
     
    590592                       &       + (   ( vn(ji+1,jj  ,jk) + vn (ji,jj,jk) ) * ( e2v(ji+1,jj  ) - e2v(ji,jj) )       & 
    591593                       &           - ( un(ji  ,jj+1,jk) + un (ji,jj,jk) ) * ( e1u(ji  ,jj+1) - e1u(ji,jj) )   )   & 
    592                        &       * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )                                               & 
     594                       &       * 0.5 / ( e1f(ji,jj) * e2f(ji,jj) )                                                & 
    593595                       &       ) * ze3f(ji,jj,jk) 
    594596               END DO 
     
    601603         ! Compute and add the vorticity term trend 
    602604         ! ---------------------------------------- 
    603          jj=2 
    604          ztne(1,:) = 0 ; ztnw(1,:) = 0 ; ztse(1,:) = 0 ; ztsw(1,:) = 0 
     605         jj = 2 
     606         ztne(1,:) = 0   ;   ztnw(1,:) = 0   ;   ztse(1,:) = 0   ;  ztsw(1,:) = 0 
    605607         DO ji = 2, jpi    
    606608               ztne(ji,jj) = zwz(ji-1,jj  ) + zwz(ji  ,jj  ) + zwz(ji  ,jj-1) 
     
    638640      !! 
    639641      !! ** Purpose :   Control the consistency between cpp options for 
    640       !!      tracer advection schemes 
     642      !!              tracer advection schemes 
    641643      !!---------------------------------------------------------------------- 
    642644      INTEGER ::   ioptio          ! temporary integer 
  • branches/dev_004_VVL/NEMO/OPA_SRC/DYN/wzvmod.F90

    r1390 r1408  
    11MODULE wzvmod 
     2!! MODULE sshwzv    
    23   !!============================================================================== 
    3    !!                       ***  MODULE  wzvmod  *** 
    4    !! Ocean diagnostic variable : vertical velocity 
     4   !!                       ***  MODULE  sshwzv  *** 
     5   !! Ocean dynamics : sea surface height and vertical velocity 
    56   !!============================================================================== 
    6    !! History :   5.0  !  90-10  (C. Levy, G. Madec)  Original code 
    7    !!             7.0  !  96-01  (G. Madec)  Statement function for e3 
    8    !!             8.5  !  02-07  (G. Madec)  Free form, F90 
    9    !!              "   !  07-07  (D. Storkey) Zero zhdiv at open boundary (BDY)  
    10    !!---------------------------------------------------------------------- 
    11    !!   wzv            : empty routine  
     7   !! History :  3.1  !  2009-02  (G. Madec, M. Leclair)  Original code 
     8   !!---------------------------------------------------------------------- 
     9 
     10   !!---------------------------------------------------------------------- 
    1211   !!   ssh_wzv        : after ssh & now vertical velocity 
    1312   !!   ssh_nxt        : filter ans swap the ssh arrays 
    1413   !!   ssh_rst        : read/write ssh restart fields in the ocean restart file 
    1514   !!---------------------------------------------------------------------- 
    16    !! * Modules used 
    1715   USE oce             ! ocean dynamics and tracers variables 
    1816   USE dom_oce         ! ocean space and time domain variables  
     
    2422   USE prtctl          ! Print control 
    2523   USE phycst 
    26    USE bdy_oce         ! unstructured open boundaries 
    2724   USE lbclnk          ! ocean lateral boundary condition (or mpp link) 
    2825   USE obc_par         ! open boundary cond. parameter 
     
    3229   PRIVATE 
    3330 
    34    !! * Routine accessibility 
    35    PUBLIC wzv          ! routine called by step.F90 and inidtr.F90 
    36    PUBLIC ssh_wzv 
    37    PUBLIC ssh_nxt 
     31   PUBLIC   ssh_wzv    ! called by step.F90 
     32   PUBLIC   ssh_nxt    ! called by step.F90 
    3833 
    3934   !! * Substitutions 
     
    4237 
    4338   !!---------------------------------------------------------------------- 
    44    !!  OPA 9.0 , LOCEAN-IPSL (2005)  
     39   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)  
    4540   !! $Id$ 
    4641   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     
    4843 
    4944CONTAINS 
    50  
    51    SUBROUTINE wzv( kt ) 
    52       !!---------------------------------------------------------------------- 
    53       !!                    ***  ROUTINE wzv  *** 
    54       !! 
    55       !! ** Purpose :   Compute the now vertical velocity after the array swap 
    56       !! 
    57       !! ** Method  :   Using the incompressibility hypothesis, the vertical 
    58       !!      velocity is computed by integrating the horizontal divergence  
    59       !!      from the bottom to the surface. 
    60       !!        The boundary conditions are w=0 at the bottom (no flux) and, 
    61       !!      in regid-lid case, w=0 at the sea surface. 
    62       !! 
    63       !! ** action  :   wn array : the now vertical velocity 
    64       !!---------------------------------------------------------------------- 
    65       INTEGER, INTENT(in) :: kt 
    66        
    67       ! Empty routine 
    68  
    69       WRITE(*,*) 'wzv : you should not be here : error ?'   
    70  
    71    END SUBROUTINE wzv 
    72  
    7345 
    7446   SUBROUTINE ssh_wzv( kt )  
     
    9062      !! if lk_vvl=T:   sshu_a, sshv_a, sshf_a  : after sea surface height 
    9163      !!                          at u-, v-, f-point s 
    92       !!                .._1    : now vertical coordinate arrays 
    9364      !!                hu, hv, hur, hvr : ocean depth and its inverse at u-,v-points 
    9465      !!---------------------------------------------------------------------- 
     
    10374      IF( kt == nit000 ) THEN 
    10475         IF(lwp) WRITE(numout,*) 
    105          IF(lwp) WRITE(numout,*) 'ssh_wzv : vertical velocity from continuity eq. (vvl option)' 
     76         IF(lwp) WRITE(numout,*) 'ssh_wzv : after sea surface height and now vertical velocity ' 
    10677         IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
    10778         ! 
     
    152123 
    153124      !                                                ! Sea surface elevation time stepping 
    154       ssha(:,:) = (  sshb(:,:) - z2dt * ( zraur * emp(:,:)                       + zhdiv(:,:) )  ) * tmask(:,:,1) 
     125      ssha(:,:) = (  sshb(:,:) - z2dt * ( zraur * emp(:,:) + zhdiv(:,:) )  ) * tmask(:,:,1) 
    155126 
    156127#if defined key_obc 
     
    158129      IF ( Agrif_Root() ) THEN  
    159130# endif 
    160          ssha(:,:) = ssha(:,:)*obctmsk(:,:) 
    161          CALL lbc_lnk(ssha,'T',1.)  ! absolutly compulsory !! (jmm) 
     131         ssha(:,:) = ssha(:,:) * obctmsk(:,:) 
     132         CALL lbc_lnk( ssha, 'T', 1. )  ! absolutly compulsory !! (jmm) 
    162133# if defined key_agrif 
    163134      ENDIF 
     
    169140         DO jj = 1, jpjm1 
    170141            DO ji = 1, fs_jpim1   ! Vector Opt. 
    171                sshu_a(ji,jj) = 0.5  * umask(ji,jj,1)  / ( e1u(ji  ,jj) * e2u(ji  ,jj) )                   & 
    172                   &                                   * ( e1t(ji  ,jj) * e2t(ji  ,jj) * ssha(ji  ,jj)     & 
    173                   &                                     + e1t(ji+1,jj) * e2t(ji+1,jj) * ssha(ji+1,jj) ) 
    174                sshv_a(ji,jj) = 0.5  * vmask(ji,jj,1)  / ( e1v(ji,jj  ) * e2v(ji,jj  ) )                   & 
    175                   &                                   * ( e1t(ji,jj  ) * e2t(ji,jj  ) * ssha(ji,jj  )     & 
    176                   &                                     + e1t(ji,jj+1) * e2t(ji,jj+1) * ssha(ji,jj+1) ) 
    177    !!gm bug used of fmask, even if thereafter multiplied by muf  which is correctly masked) 
    178                sshf_a(ji,jj) = 0.25 * fmask(ji,jj,1) * ( ssha(ji  ,jj) + ssha(ji  ,jj+1)                 & 
     142               sshu_a(ji,jj) = 0.5  * umask(ji,jj,1) / ( e1u(ji  ,jj) * e2u(ji  ,jj) )                   & 
     143                  &                                  * ( e1t(ji  ,jj) * e2t(ji  ,jj) * ssha(ji  ,jj)     & 
     144                  &                                    + e1t(ji+1,jj) * e2t(ji+1,jj) * ssha(ji+1,jj) ) 
     145               sshv_a(ji,jj) = 0.5  * vmask(ji,jj,1) / ( e1v(ji,jj  ) * e2v(ji,jj  ) )                   & 
     146                  &                                  * ( e1t(ji,jj  ) * e2t(ji,jj  ) * ssha(ji,jj  )     & 
     147                  &                                    + e1t(ji,jj+1) * e2t(ji,jj+1) * ssha(ji,jj+1) ) 
     148               sshf_a(ji,jj) = 0.25 * umask(ji,jj,1) * umask (ji,jj+1,1)   &   ! Caution : fmask not used 
     149                  &                                  * ( ssha(ji  ,jj) + ssha(ji  ,jj+1)                 & 
    179150                  &                                    + ssha(ji+1,jj) + ssha(ji+1,jj+1) ) 
    180151            END DO 
     
    199170      !                                           !------------------------------! 
    200171      IF( lk_vvl ) THEN                           ! only in vvl case) 
    201       !                                                ! now local depth and scale factors (stored in fse3. arrays) 
     172         !                                             ! now local depth and scale factors (stored in fse3. arrays) 
    202173         DO jk = 1, jpkm1 
    203174            fsdept(:,:,jk) = fsdept_n(:,:,jk)               ! depths 
     
    217188         hv(:,:) = hv_0(:,:) + sshv_n(:,:) 
    218189         !                                             ! masked inverse of the ocean depth (at u- and v-points) 
    219          hur(:,:) = 1. / MAX( hu(:,:), fse3u_0(:,:,1) ) * umask(:,:,1) 
    220          hvr(:,:) = 1. / MAX( hv(:,:), fse3v_0(:,:,1) ) * vmask(:,:,1) 
    221 !!gm to be corrected (the above case does not work properly with 1 ocean level only) 
    222 !        hur(:,:) = 1. / MAX( hu(:,:), 1.e-15 ) * umask(:,:,1) 
    223 !        hvr(:,:) = 1. / MAX( hv(:,:), 1.e-15 ) * vmask(:,:,1) 
    224 !!gm 
    225 !!add end 
     190         hur(:,:) = umask(:,:,1) / ( hu(:,:) + 1.e0 - umask(:,:,1) ) 
     191         hvr(:,:) = vmask(:,:,1) / ( hv(:,:) + 1.e0 - vmask(:,:,1) ) 
    226192 
    227193      ENDIF 
     
    247213      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index 
    248214      !! 
    249       INTEGER  ::   ji, jj, jk           ! dummy loop indices 
     215      INTEGER  ::   ji, jj               ! dummy loop indices 
    250216      !!---------------------------------------------------------------------- 
    251217 
     
    330296         ENDIF 
    331297      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN 
    332          CALL iom_rstput( kt, nitrst, numrow, 'sshb'  , sshb  (:,:) ) 
    333          CALL iom_rstput( kt, nitrst, numrow, 'sshn'  , sshn  (:,:) ) 
     298         CALL iom_rstput( kt, nitrst, numrow, 'sshb'  , sshb(:,:) ) 
     299         CALL iom_rstput( kt, nitrst, numrow, 'sshn'  , sshn(:,:) ) 
    334300      ENDIF 
    335301      ! 
  • branches/dev_004_VVL/NEMO/OPA_SRC/TRA/tranxt.F90

    r1382 r1408  
    1414   !!            2.0  !  2006-02  (L. Debreu, C. Mazauric) Agrif implementation 
    1515   !!            3.0  !  2008-06  (G. Madec)  time stepping always done in trazdf 
     16   !!            3.1  !  2009-02  (G. Madec, R. Benshila)  re-introduce the vvl option 
    1617   !!---------------------------------------------------------------------- 
    1718 
    1819   !!---------------------------------------------------------------------- 
    1920   !!   tra_nxt       : time stepping on temperature and salinity 
     21   !!   tra_nxt_fix   : time stepping on temperature and salinity : fixed    volume case 
     22   !!   tra_nxt_vvl   : time stepping on temperature and salinity : variable volume case 
    2023   !!---------------------------------------------------------------------- 
    2124   USE oce             ! ocean dynamics and tracers variables 
     
    127130#endif       
    128131 
    129       ! trends diagnostics :  Asselin filter trend : (tb filtered - tb)/2dt 
    130       IF( l_trdtra ) THEN      
     132      ! trends computation 
     133      IF( l_trdtra ) THEN              ! trend of the Asselin filter (tb filtered - tb)/dt      
    131134         DO jk = 1, jpkm1 
    132135            zfact = 1.e0 / r2dt_t(jk)              
     
    136139         CALL trd_mod( ztrdt, ztrds, jptra_trd_atf, 'TRA', kt ) 
    137140      END IF 
     141 
    138142      !                        ! control print 
    139143      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tn, clinfo1=' nxt  - Tn: ', mask1=tmask,   & 
     
    182186         !                                           ! ----------------------- ! 
    183187         ! 
    184          IF( neuler == 0 .AND. kt == nit000 ) THEN        ! case of Euler time-stepping at first time-step 
    185             DO jk = 1, jpkm1 
    186                DO jj = 1, jpj 
    187                   DO ji = 1, jpi 
    188                      ztm = 0.25 * ( ta(ji,jj,jk) + 2. * tn(ji,jj,jk) + tb(ji,jj,jk) )      ! mean t 
    189                      zsm = 0.25 * ( sa(ji,jj,jk) + 2. * sn(ji,jj,jk) + sb(ji,jj,jk) ) 
    190                      tb(ji,jj,jk) = tn(ji,jj,jk)                                           ! tb <-- tn 
    191                      sb(ji,jj,jk) = sn(ji,jj,jk) 
     188         IF( neuler == 0 .AND. kt == nit000 ) THEN        ! Euler time-stepping at first time-step 
     189            DO jk = 1, jpkm1                              ! (only swap) 
     190               DO jj = 1, jpj 
     191                  DO ji = 1, jpi 
    192192                     tn(ji,jj,jk) = ta(ji,jj,jk)                                           ! tb <-- tn 
    193193                     sn(ji,jj,jk) = sa(ji,jj,jk) 
    194                      ta(ji,jj,jk) = ztm                                                    ! ta <-- mean t 
    195                      sa(ji,jj,jk) = zsm 
    196194                  END DO 
    197195               END DO 
     
    219217         !                                           ! ----------------------- ! 
    220218         ! 
    221          IF( neuler == 0 .AND. kt == nit000 ) THEN        ! case of Euler time-stepping at first time-step 
    222             DO jk = 1, jpkm1 
    223                DO jj = 1, jpj 
    224                   DO ji = 1, jpi 
    225                      tb(ji,jj,jk) = tn(ji,jj,jk)                                           ! tb <-- tn 
    226                      sb(ji,jj,jk) = sn(ji,jj,jk) 
     219         IF( neuler == 0 .AND. kt == nit000 ) THEN        ! Euler time-stepping at first time-step 
     220            DO jk = 1, jpkm1 
     221               DO jj = 1, jpj 
     222                  DO ji = 1, jpi 
    227223                     tn(ji,jj,jk) = ta(ji,jj,jk)                                           ! tn <-- ta 
    228224                     sn(ji,jj,jk) = sa(ji,jj,jk) 
     
    234230               DO jj = 1, jpj 
    235231                  DO ji = 1, jpi 
    236 !RBvvl for reproducibility 
    237 !                     ztf = atfp * ( ta(ji,jj,jk) - 2.* tn(ji,jj,jk) + tb(ji,jj,jk) )       ! Asselin filter on t  
    238 !                     zsf = atfp * ( sa(ji,jj,jk) - 2.* sn(ji,jj,jk) + sb(ji,jj,jk) ) 
    239 !                     tb(ji,jj,jk) = tn(ji,jj,jk) + ztf                                     ! tb <-- filtered tn  
    240 !                     sb(ji,jj,jk) = sn(ji,jj,jk) + zsf 
    241                      tb(ji,jj,jk) = atfp * ( tb(ji,jj,jk) + ta(ji,jj,jk) ) + atfp1 * tn(ji,jj,jk) 
    242                      sb(ji,jj,jk) = atfp * ( sb(ji,jj,jk) + sa(ji,jj,jk) ) + atfp1 * sn(ji,jj,jk) 
     232                     ztf = atfp * ( ta(ji,jj,jk) - 2.* tn(ji,jj,jk) + tb(ji,jj,jk) )       ! Asselin filter on t  
     233                     zsf = atfp * ( sa(ji,jj,jk) - 2.* sn(ji,jj,jk) + sb(ji,jj,jk) ) 
     234                     tb(ji,jj,jk) = tn(ji,jj,jk) + ztf                                     ! tb <-- filtered tn  
     235                     sb(ji,jj,jk) = sn(ji,jj,jk) + zsf 
    243236                     tn(ji,jj,jk) = ta(ji,jj,jk)                                           ! tn <-- ta 
    244237                     sn(ji,jj,jk) = sa(ji,jj,jk) 
     
    295288         !                                           ! ----------------------- ! 
    296289         ! 
    297          IF( neuler == 0 .AND. kt == nit000 ) THEN 
    298             DO jk = 1, jpkm1 
    299                DO jj = 1, jpj 
    300                   DO ji = 1, jpi 
    301                      ze3t_b = fse3t_b(ji,jj,jk) 
    302                      ze3t_n = fse3t_n(ji,jj,jk) 
    303                      ze3t_a = fse3t_a(ji,jj,jk) 
    304                      !                                         ! tracer content at Before, now and after 
    305                      ztcb = tb(ji,jj,jk) *  ze3t_b   ;   zscb = sb(ji,jj,jk) * ze3t_b 
    306                      ztcn = tn(ji,jj,jk) *  ze3t_n   ;   zscn = sn(ji,jj,jk) * ze3t_n 
    307                      ztca = ta(ji,jj,jk) *  ze3t_a   ;   zsca = sa(ji,jj,jk) * ze3t_a 
    308                      ! 
    309                      !                                         ! mean thickness and tracer 
    310                      ze3mr= 1.e0  / ( ze3t_a + 2.* ze3t_n + ze3t_b ) 
    311                      ztm  = ze3mr * ( ztca   + 2.* ztcn   + ztcb   ) 
    312                      zsm  = ze3mr * ( zsca   + 2.* zscn   + zscb   ) 
    313 !!gm mean e3t have to be saved and used in dynhpg  or it can be recomputed in dynhpg !! 
    314 !!gm                 e3t_m(ji,jj,jk) = 0.25 / ze3mr 
    315                      !                                         ! swap of arrays 
    316                      tb(ji,jj,jk) = tn(ji,jj,jk)                    ! tb <-- tn 
    317                      sb(ji,jj,jk) = sn(ji,jj,jk) 
     290         IF( neuler == 0 .AND. kt == nit000 ) THEN        ! Euler time-stepping at first time-step 
     291            DO jk = 1, jpkm1                              ! (only swap) 
     292               DO jj = 1, jpj 
     293                  DO ji = 1, jpi 
    318294                     tn(ji,jj,jk) = ta(ji,jj,jk)                    ! tn <-- ta 
    319295                     sn(ji,jj,jk) = sa(ji,jj,jk) 
    320                      ta(ji,jj,jk) = ztm                             ! ta <-- mean t 
    321                      sa(ji,jj,jk) = zsm 
    322296                  END DO 
    323297               END DO 
     
    369343               DO jj = 1, jpj                             ! ONLY swap 
    370344                  DO ji = 1, jpi 
    371                      tb(ji,jj,jk) = tn(ji,jj,jk)                                 ! tb <-- tn 
    372                      sb(ji,jj,jk) = sn(ji,jj,jk) 
    373345                     tn(ji,jj,jk) = ta(ji,jj,jk)                                 ! tn <-- ta 
    374346                     sn(ji,jj,jk) = sa(ji,jj,jk) 
     
    394366                     zsc_f  = atfp * ( zsca   - 2.* zscn   + zscb   )  
    395367                     ! 
    396 !!gm tmask useless below 
    397368                     !                                         ! filtered tracer including the correction  
    398                      ze3fr = tmask(ji,jj,jk) / ( ze3t_n + ze3t_f ) 
     369                     ze3fr = 1.e0 / ( ze3t_n + ze3t_f ) 
    399370                     ztf   =  ( ztcn  + ztc_f ) * ze3fr 
    400371                     zsf   =  ( zscn  + zsc_f ) * ze3fr 
  • branches/dev_004_VVL/NEMO/OPA_SRC/TRA/trazdf_exp.F90

    r1362 r1408  
    4242#  include "vectopt_loop_substitute.h90" 
    4343   !!---------------------------------------------------------------------- 
    44    !! NEMO/OPA  3.0 , LOCEAN-IPSL (2008)  
     44   !! NEMO/OPA  3.2 , LOCEAN-IPSL (2009)  
    4545   !! $Id$  
    4646   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     
    7878      INTEGER  ::   ji, jj, jk, jl            ! dummy loop indices 
    7979      REAL(wp) ::   zlavmr, zave3r, ze3tr     ! temporary scalars 
    80       REAL(wp) ::   zta, zsa, ze3tb, zcoef    ! temporary scalars 
     80      REAL(wp) ::   zta, zsa, ze3tb           ! temporary scalars 
    8181      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zwx, zwy, zwz, zww   ! 3D workspace 
    8282      !!--------------------------------------------------------------------- 
     
    130130            DO jj = 2, jpjm1  
    131131               DO ji = fs_2, fs_jpim1   ! vector opt. 
    132                   ze3tb = fse3t_b(ji,jj,jk)                                            ! before e3t 
     132                  ze3tb = fse3t_b(ji,jj,jk) / fse3t(ji,jj,jk)                          ! before e3t 
    133133                  zta   = zwx(ji,jj,jk) - tb(ji,jj,jk) + p2dt(jk) * ta(ji,jj,jk)       ! total trends * 2*rdt 
    134134                  zsa   = zwz(ji,jj,jk) - sb(ji,jj,jk) + p2dt(jk) * sa(ji,jj,jk)      
    135                   zcoef = 1.e0 / fse3t_n(ji,jj,jk) * tmask(ji,jj,jk) 
    136 !RBvvl : which e3t ? 
    137                   ta(ji,jj,jk) = (  ze3tb * tb(ji,jj,jk) + fse3t(ji,jj,jk) * zta  ) * zcoef 
    138                   sa(ji,jj,jk) = (  ze3tb * sb(ji,jj,jk) + fse3t(ji,jj,jk) * zsa  ) * zcoef 
     135                  ta(ji,jj,jk) = (  ze3tb * tb(ji,jj,jk) + zta  ) * tmask(ji,jj,jk) 
     136                  sa(ji,jj,jk) = (  ze3tb * sb(ji,jj,jk) + zsa  ) * tmask(ji,jj,jk) 
    139137               END DO 
    140138            END DO 
     
    144142            DO jj = 2, jpjm1  
    145143               DO ji = fs_2, fs_jpim1   ! vector opt. 
    146                   ta(ji,jj,jk) = ( zwx(ji,jj,jk) + p2dt(jk) * ta(ji,jj,jk) ) *tmask(ji,jj,jk) 
    147                   sa(ji,jj,jk) = ( zwz(ji,jj,jk) + p2dt(jk) * sa(ji,jj,jk) ) *tmask(ji,jj,jk) 
     144                  ta(ji,jj,jk) = ( zwx(ji,jj,jk) + p2dt(jk) * ta(ji,jj,jk) ) * tmask(ji,jj,jk) 
     145                  sa(ji,jj,jk) = ( zwz(ji,jj,jk) + p2dt(jk) * sa(ji,jj,jk) ) * tmask(ji,jj,jk) 
    148146               END DO 
    149147            END DO 
  • branches/dev_004_VVL/NEMO/OPA_SRC/TRA/trazdf_imp.F90

    r1362 r1408  
    11MODULE trazdf_imp 
    2    !!============================================================================== 
     2   !!====================================================================== 
    33   !!                 ***  MODULE  trazdf_imp  *** 
    44   !! Ocean active tracers:  vertical component of the tracer mixing trend 
    5    !!============================================================================== 
    6    !! History : 
    7    !!   6.0  !  90-10  (B. Blanke)  Original code 
    8    !!   7.0  !  91-11 (G. Madec) 
    9    !!        !  92-06 (M. Imbard) correction on tracer trend loops 
    10    !!        !  96-01 (G. Madec) statement function for e3 
    11    !!        !  97-05 (G. Madec) vertical component of isopycnal 
    12    !!        !  97-07 (G. Madec) geopotential diffusion in s-coord 
    13    !!        !  00-08 (G. Madec) double diffusive mixing 
    14    !!   8.5  !  02-08 (G. Madec)  F90: Free form and module 
    15    !!   9.0  !  06-11 (G. Madec)  New step reorganisation 
     5   !!====================================================================== 
     6   !! History :  OPA  !  1990-10  (B. Blanke)  Original code 
     7   !!            7.0  !  1991-11  (G. Madec) 
     8   !!                 !  1992-06  (M. Imbard) correction on tracer trend loops 
     9   !!                 !  1996-01  (G. Madec) statement function for e3 
     10   !!                 !  1997-05  (G. Madec) vertical component of isopycnal 
     11   !!                 !  1997-07  (G. Madec) geopotential diffusion in s-coord 
     12   !!                 !  2000-08  (G. Madec) double diffusive mixing 
     13   !!   NEMO     1.0  !  2002-08  (G. Madec) F90: Free form and module 
     14   !!            2.0  !  2006-11  (G. Madec) New step reorganisation 
     15   !!            3.2  !  2009-03  (G. Madec)  heat and salt content trends 
     16   !!---------------------------------------------------------------------- 
     17   
    1618   !!---------------------------------------------------------------------- 
    1719   !!   tra_zdf_imp : Update the tracer trend with the diagonal vertical   
    1820   !!                 part of the mixing tensor. 
    1921   !!---------------------------------------------------------------------- 
    20    !! * Modules used 
    2122   USE oce             ! ocean dynamics and tracers variables 
    2223   USE dom_oce         ! ocean space and time domain variables  
     
    3637   PRIVATE 
    3738 
    38    !! * Routine accessibility 
    39    PUBLIC tra_zdf_imp   !  routine called by step.F90 
     39   PUBLIC   tra_zdf_imp   !  routine called by step.F90 
    4040 
    4141   !! * Substitutions 
     
    4545#  include "vectopt_loop_substitute.h90" 
    4646   !!---------------------------------------------------------------------- 
    47    !!---------------------------------------------------------------------- 
    48    !!  OPA 9.0 , LOCEAN-IPSL (2005)  
     47   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)  
    4948   !! $Id$ 
    5049   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     
    9089      !! 
    9190      !!--------------------------------------------------------------------- 
    92       !! * Modules used 
    93       USE oce    , ONLY :   zwd   => ua,  &  ! ua used as workspace 
    94                             zws   => va      ! va   "          " 
    95       !! * Arguments 
    96       INTEGER, INTENT( in ) ::   kt          ! ocean time-step index 
    97       REAL(wp), DIMENSION(jpk), INTENT( in ) ::   & 
    98          p2dt                                ! vertical profile of tracer time-step 
    99  
    100       !! * Local declarations 
    101       INTEGER  ::   ji, jj, jk               ! dummy loop indices 
    102       REAL(wp) ::   zavi, zrhs, znvvl,     & ! temporary scalars 
    103          ze3tb, ze3tn, ze3ta                 ! variable vertical scale factors 
    104       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   & 
    105          zwi, zwt, zavsi                     ! workspace arrays 
     91      USE oce    , ONLY :   zwd   => ua   ! ua used as workspace 
     92      USE oce    , ONLY :   zws   => va   ! va  -          - 
     93      !! 
     94      INTEGER                 , INTENT(in) ::   kt     ! ocean time-step index 
     95      REAL(wp), DIMENSION(jpk), INTENT(in) ::   p2dt   ! vertical profile of tracer time-step 
     96      !! 
     97      INTEGER  ::   ji, jj, jk            ! dummy loop indices 
     98      REAL(wp) ::   zavi, zrhs, znvvl     ! temporary scalars 
     99      REAL(wp) ::   ze3tb, ze3tn, ze3ta   ! variable vertical scale factors 
     100      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zwi, zwt, zavsi   ! workspace arrays 
    106101      !!--------------------------------------------------------------------- 
    107102 
     
    169164         DO jj = 2, jpjm1 
    170165            DO ji = fs_2, fs_jpim1   ! vector opt. 
    171                ze3ta = ( 1. - znvvl ) + znvvl*fse3t_a(ji,jj,jk)                        ! after scale factor at T-point 
    172                ze3tn = ( 1. - znvvl )*fse3t_n(ji,jj,jk) + znvvl                        ! now   scale factor at T-point 
     166               ze3ta =  ( 1. - znvvl )   &                            ! after scale factor at T-point 
     167                  &   +        znvvl    * fse3t_a(ji,jj,jk)  
     168               ze3tn =    znvvl          &                            ! now   scale factor at T-point 
     169                  &   + ( 1. - znvvl ) * fse3t_n(ji,jj,jk) 
    173170               zwi(ji,jj,jk) = - p2dt(jk) * zwt(ji,jj,jk  ) / ( ze3tn * fse3w(ji,jj,jk  ) ) 
    174171               zws(ji,jj,jk) = - p2dt(jk) * zwt(ji,jj,jk+1) / ( ze3tn * fse3w(ji,jj,jk+1) ) 
     
    181178      DO jj = 2, jpjm1 
    182179         DO ji = fs_2, fs_jpim1   ! vector opt. 
    183             ze3ta = ( 1. - znvvl ) + znvvl*fse3t_a(ji,jj,1)                             ! after scale factor at T-point 
     180            ze3ta = ( 1. - znvvl )        &                           ! after scale factor at T-point 
     181               &   +       znvvl      * fse3t_a(ji,jj,1)  
    184182            zwi(ji,jj,1) = 0.e0 
    185183            zwd(ji,jj,1) = ze3ta - zws(ji,jj,1) 
     
    225223      DO jj = 2, jpjm1 
    226224         DO ji = fs_2, fs_jpim1 
    227             ze3tb = ( 1. - znvvl ) + znvvl*fse3t_b(ji,jj,1) 
    228 !RBvvl which  
    229             ze3tn = ( 1. - znvvl ) + znvvl*fse3t(ji,jj,1) 
     225            ze3tb = ( 1. - znvvl ) + znvvl * fse3t_b(ji,jj,1) 
     226            ze3tn = ( 1. - znvvl ) + znvvl * fse3t(ji,jj,1) 
    230227            ta(ji,jj,1) = ze3tb * tb(ji,jj,1) + p2dt(1) * ze3tn * ta(ji,jj,1) 
    231228         END DO 
     
    234231         DO jj = 2, jpjm1 
    235232            DO ji = fs_2, fs_jpim1 
    236                ze3tb = ( 1. - znvvl ) + znvvl*fse3t_b(ji,jj,jk) 
    237 !RB_vvl 
    238                ze3tn = ( 1. - znvvl ) + znvvl*fse3t (ji,jj,jk) 
     233               ze3tb = ( 1. - znvvl ) + znvvl * fse3t_b(ji,jj,jk) 
     234               ze3tn = ( 1. - znvvl ) + znvvl * fse3t  (ji,jj,jk) 
    239235               zrhs = ze3tb * tb(ji,jj,jk) + p2dt(jk) * ze3tn * ta(ji,jj,jk)   ! zrhs=right hand side  
    240                ta(ji,jj,jk) = zrhs - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) *ta(ji,jj,jk-1) 
     236               ta(ji,jj,jk) = zrhs - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) * ta(ji,jj,jk-1) 
    241237            END DO 
    242238         END DO 
     
    269265         DO jj = 2, jpjm1 
    270266            DO ji = fs_2, fs_jpim1   ! vector opt. 
    271                ze3ta = ( 1. - znvvl ) + znvvl*fse3t_a(ji,jj,jk)                            ! after scale factor at T-point 
    272                ze3tn = ( 1. - znvvl )*fse3t_n(ji,jj,jk) + znvvl                            ! now   scale factor at T-point 
     267               ze3ta =  ( 1. - znvvl )                        &         ! after scale factor at T-point 
     268                  &   +        znvvl   * fse3t_a(ji,jj,jk)            
     269               ze3tn =    znvvl                               &         ! now   scale factor at T-point 
     270                  &   + ( 1. - znvvl ) * fse3t_n(ji,jj,jk) 
    273271               zwi(ji,jj,jk) = - p2dt(jk) * zavsi(ji,jj,jk  ) / ( ze3tn * fse3w(ji,jj,jk  ) ) 
    274272               zws(ji,jj,jk) = - p2dt(jk) * zavsi(ji,jj,jk+1) / ( ze3tn * fse3w(ji,jj,jk+1) ) 
     
    281279      DO jj = 2, jpjm1 
    282280         DO ji = fs_2, fs_jpim1   ! vector opt. 
    283             ze3ta = ( 1. - znvvl ) + znvvl*fse3t_a(ji,jj,1) 
     281            ze3ta = ( 1. - znvvl ) + znvvl * fse3t_a(ji,jj,1) 
    284282            zwi(ji,jj,1) = 0.e0 
    285283            zwd(ji,jj,1) = ze3ta - zws(ji,jj,1) 
     
    324322      DO jj = 2, jpjm1 
    325323         DO ji = fs_2, fs_jpim1 
    326             ze3tb = ( 1. - znvvl ) + znvvl*fse3t_b(ji,jj,1)                      ! before scale factor at T-point 
    327             ze3tn = ( 1. - znvvl ) + znvvl*fse3t(ji,jj,1)                        ! now    scale factor at T-point 
     324            ze3tb = ( 1. - znvvl )   &                                           ! before scale factor at T-point 
     325               &   +  znvvl       * fse3t_b(ji,jj,1) 
     326            ze3tn = ( 1. - znvvl ) + znvvl * fse3t  (ji,jj,1)                    ! now    scale factor at T-point 
    328327            sa(ji,jj,1) = ze3tb * sb(ji,jj,1) + p2dt(1) * ze3tn * sa(ji,jj,1) 
    329328         END DO 
     
    332331         DO jj = 2, jpjm1 
    333332            DO ji = fs_2, fs_jpim1 
    334                ze3tb = ( 1. - znvvl ) + znvvl*fse3t_b(ji,jj,jk)                  ! before scale factor at T-point 
    335                ze3tn = ( 1. - znvvl ) + znvvl*fse3t(ji,jj,jk)                    ! now    scale factor at T-point 
     333               ze3tb = ( 1. - znvvl )   &                                        ! before scale factor at T-point 
     334                  &   +  znvvl       * fse3t_b(ji,jj,jk) 
     335               ze3tn = ( 1. - znvvl ) + znvvl * fse3t  (ji,jj,jk)                ! now    scale factor at T-point 
    336336               zrhs = ze3tb * sb(ji,jj,jk) + p2dt(jk) * ze3tn * sa(ji,jj,jk)     ! zrhs=right hand side 
    337337               sa(ji,jj,jk) = zrhs - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) *sa(ji,jj,jk-1) 
     
    355355         END DO 
    356356      END DO 
    357  
     357      ! 
    358358   END SUBROUTINE tra_zdf_imp 
    359359 
  • branches/dev_004_VVL/NEMO/OPA_SRC/ZDF/zdfevd.F90

    r1388 r1408  
    55   !!                of vertical eddy mixing coefficient 
    66   !!====================================================================== 
     7   !! History :  OPA  !  1997-06  (G. Madec, A. Lazar)  Original code 
     8   !!   NEMO     1.0  !  2002-06  (G. Madec)  F90: Free form and module 
     9   !!             -   !  2005-06  (C. Ethe) KPP parameterization 
     10   !!            3.2  !  2009-03  (M. Leclair, G. Madec, R. Benshila) test on both before & after 
     11   !!---------------------------------------------------------------------- 
    712 
    813   !!---------------------------------------------------------------------- 
    9    !!   zdf_evd      : update momentum and tracer Kz at the location of 
    10    !!                  statically unstable portion of the water column 
    11    !!                  (called if ln_zdfevd=T) 
     14   !!   zdf_evd      : increase the momentum and tracer Kz at the location of 
     15   !!                  statically unstable portion of the water column (ln_zdfevd=T) 
    1216   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
    1417   USE oce             ! ocean dynamics and tracers variables 
    1518   USE dom_oce         ! ocean space and time domain variables 
     
    2225   PRIVATE 
    2326 
    24    !! * Routine accessibility 
    25    PUBLIC zdf_evd      ! called by step.F90 
     27   PUBLIC   zdf_evd    ! called by step.F90 
    2628 
    2729   !! * Substitutions 
    2830#  include "domzgr_substitute.h90" 
    2931   !!---------------------------------------------------------------------- 
    30    !!   OPA 9.0 , LOCEAN-IPSL (2005)  
    31    !! $Id$  
    32    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     32   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)  
     33   !! $Id:$ 
     34   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)  
    3335   !!---------------------------------------------------------------------- 
    3436 
     
    4850      !! ** Action  :   Update avt, avmu, avmv in statically instable cases 
    4951      !!                and avt_evd which is avt due to convection 
    50       !! References : 
    51       !!      Lazar, A., these de l'universite Paris VI, France, 1997 
    52       !! History : 
    53       !!   7.0  !  97-06  (G. Madec, A. Lazar)  Original code 
    54       !!   8.5  !  02-06  (G. Madec)  F90: Free form and module 
    55       !!   9.0  !  05-06  (C. Ethe) KPP parameterization 
     52      !! References :   Lazar, A., these de l'universite Paris VI, France, 1997 
    5653      !!---------------------------------------------------------------------- 
    57       !! * Arguments 
    5854      INTEGER, INTENT( in ) ::   kt         ! ocean time-step indexocean time step 
    59  
    60       !! * Local declarations 
     55      !! 
    6156      INTEGER ::   ji, jj, jk               ! dummy loop indices 
    6257      !!---------------------------------------------------------------------- 
     
    8075            !                                             ! =============== 
    8176#if defined key_vectopt_loop 
    82 !!!         WHERE( rn2(:,:,jk) <= -1.e-12 ) avt(:,:,jk) = tmask(:,:,jk) * avevd   ! agissant sur T SEUL! 
    83             jj = 1                     ! big loop forced 
    84             DO ji = jpi+2, jpij    
    85 # if defined key_zdfkpp 
    86 !! no implicit mixing in the boundary layer with KPP 
    87                IF( ( MIN( rn2(ji,jj,jk), rn2b(ji,jj,jk) ) <= -1.e-12 ) .AND. ( fsdepw(ji,jj,jk) > hkpp(ji,jj) ) ) THEN 
    88 # else 
    89                IF(   MIN( rn2(ji,jj,jk), rn2b(ji,jj,jk) ) <= -1.e-12 ) THEN 
    90 # endif 
    91                   avt (ji  ,jj  ,jk) = avevd * tmask(ji  ,jj  ,jk) 
    92                   avmu(ji  ,jj  ,jk) = avevd * umask(ji  ,jj  ,jk) 
    93                   avmu(ji-1,jj  ,jk) = avevd * umask(ji-1,jj  ,jk) 
    94                   avmv(ji  ,jj  ,jk) = avevd * vmask(ji  ,jj  ,jk) 
    95                   avmv(ji  ,jj-1,jk) = avevd * vmask(ji  ,jj-1,jk) 
    96                ENDIF 
    97             END DO 
     77            DO jj = 1, 1                     ! big loop forced 
     78               DO ji = jpi+2, jpij    
    9879#else 
    9980            DO jj = 2, jpj             ! no vector opt. 
    10081               DO ji = 2, jpi 
    101 # if defined key_zdfkpp 
     82#endif 
     83#if defined key_zdfkpp 
    10284!! no implicit mixing in the boundary layer with KPP 
    103                IF( ( MIN( rn2(ji,jj,jk),  rn2b(ji,jj,jk) ) <= -1.e-12 ) .AND. ( fsdepw(ji,jj,jk) > hkpp(ji,jj) ) ) THEN 
    104 # else 
    105                IF(   MIN( rn2(ji,jj,jk),  rn2b(ji,jj,jk) ) <= -1.e-12 ) THEN 
    106 # endif 
     85                  IF( ( MIN( rn2(ji,jj,jk),  rn2b(ji,jj,jk) ) <= -1.e-12 ) .AND. ( fsdepw(ji,jj,jk) > hkpp(ji,jj) ) ) THEN 
     86#else 
     87                  IF(   MIN( rn2(ji,jj,jk),  rn2b(ji,jj,jk) ) <= -1.e-12 ) THEN 
     88#endif 
    10789                     avt (ji  ,jj  ,jk) = avevd * tmask(ji  ,jj  ,jk) 
    10890                     avmu(ji  ,jj  ,jk) = avevd * umask(ji  ,jj  ,jk) 
     
    11395               END DO 
    11496            END DO 
    115 #endif 
    11697            !                                             ! =============== 
    11798         END DO                                           !   End of slab 
     
    130111!!!         WHERE( rn2(:,:,jk) <= -1.e-12 ) avt(:,:,jk) = tmask(:,:,jk) * avevd   ! agissant sur T SEUL!  
    131112#if defined key_vectopt_loop 
    132             jj = 1                     ! big loop forced 
    133             DO ji = 1, jpij    
    134 # if defined key_zdfkpp 
    135 !! no implicit mixing in the boundary layer with KPP 
    136                IF( ( MIN( rn2(ji,jj,jk), rn2b(ji,jj,jk) ) <= -1.e-12 ) .AND. ( fsdepw(ji,jj,jk) > hkpp(ji,jj) ) ) &               
    137 # else 
    138                IF(   MIN( rn2(ji,jj,jk), rn2b(ji,jj,jk) ) <= -1.e-12 ) & 
    139 # endif 
    140                   avt(ji,jj,jk) = avevd * tmask(ji,jj,jk) 
    141             END DO 
     113            DO jj = 1, 1                     ! big loop forced 
     114               DO ji = 1, jpij    
    142115#else 
    143116            DO jj = 1, jpj             ! loop over the whole domain (no lbc_lnk call) 
    144117               DO ji = 1, jpi 
    145 # if defined key_zdfkpp 
     118#endif 
     119#if defined key_zdfkpp 
    146120!! no implicit mixing in the boundary layer with KPP 
    147121                  IF( ( MIN( rn2(ji,jj,jk), rn2b(ji,jj,jk) ) <= -1.e-12 ) .AND. ( fsdepw(ji,jj,jk) > hkpp(ji,jj) ) ) &           
    148 # else 
     122#else 
    149123                  IF(   MIN( rn2(ji,jj,jk), rn2b(ji,jj,jk) ) <= -1.e-12 ) & 
    150 # endif 
     124#endif 
    151125                     avt(ji,jj,jk) = avevd * tmask(ji,jj,jk) 
    152126               END DO 
    153127            END DO 
    154 #endif 
    155128            !                                             ! =============== 
    156129         END DO                                           !   End of slab 
     
    159132 
    160133      ! update of avt_evd and avmu_evd 
    161       avt_evd  (:,:,:) = avt (:,:,:)  - avt_evd (:,:,:)  
    162       avmu_evd (:,:,:) = avmu(:,:,:)  - avmu_evd (:,:,:)  
     134      avt_evd (:,:,:) = avt (:,:,:)  - avt_evd (:,:,:)  
     135      avmu_evd(:,:,:) = avmu(:,:,:)  - avmu_evd(:,:,:)  
    163136 
    164137   END SUBROUTINE zdf_evd 
  • branches/dev_004_VVL/NEMO/OPA_SRC/oce.F90

    r1388 r1408  
    44   !! Ocean        :  dynamics and active tracers defined in memory  
    55   !!====================================================================== 
    6    !! History : 
    7    !!   8.5  !  02-11  (G. Madec)  F90: Free form and module 
    8    !!   9.0  !  05-11  (V. Garnier) Surface pressure gradient organization 
     6   !! History :  0.1  !  2002-11  (G. Madec)  F90: Free form and module 
     7   !!            1.0  !  2005-11  (V. Garnier) Surface pressure gradient organization 
     8   !!            3.1  !  2009-02  (G. Madec, M. Leclair)  pure z* coordinate 
    99   !!---------------------------------------------------------------------- 
    10    !!  OPA 9.0 , LOCEAN-IPSL (2005)  
    11    !! $Id$  
    12    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    13    !!---------------------------------------------------------------------- 
    14    !! * Modules used 
    1510   USE par_oce      ! ocean parameters 
    1611 
     
    2015   !! Physics and algorithm flags 
    2116   !! --------------------------- 
    22    LOGICAL, PUBLIC ::   l_traldf_rot    = .FALSE.  !: rotated laplacian operator for lateral diffusion  
     17   LOGICAL, PUBLIC ::   l_traldf_rot    = .FALSE.  !: rotated laplacian operator for lateral diffusion 
    2318   LOGICAL, PUBLIC ::   ln_dynhpg_imp   = .FALSE.  !: semi-implicite hpg flag 
    2419   INTEGER, PUBLIC ::   nn_dynhpg_rst   = 0        !: add dynhpg implicit variables in restart ot not 
    2520 
    26    !! dynamics and tracer fields 
    27    !! -------------------------- 
    28    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   &   !: 
    29       ! before !  now      !  after  !      ! the after trends becomes the fields 
    30       ! fields !  fields   !  trends !      ! only in dyn(tra)_zdf and dyn(tra)_nxt 
    31       ub       ,  un       ,  ua     ,   &  !: i-horizontal velocity (m/s) 
    32       vb       ,  vn       ,  va     ,   &  !: j-horizontal velocity (m/s) 
    33                   wn       ,             &  !: vertical velocity (m/s) 
    34       rotb     ,  rotn     ,             &  !: relative vorticity (1/s) 
    35       hdivb    ,  hdivn    ,             &  !: horizontal divergence (1/s) 
    36       tb       ,  tn       ,  ta     ,   &  !: potential temperature (celcius) 
    37       sb       ,  sn       ,  sa            !: salinity (psu) 
    38    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   &   !: 
    39       rhd ,                              &  !: in situ density anomalie rhd=(rho-rau0)/rau0 (no units) 
    40       rhop,                              &  !: potential volumic mass (kg/m3) 
    41       rn2,                               &  !: now    brunt-vaisala frequency (1/s2) 
    42       rn2b                                  !: before brunt-vaisala frequency (1/s2) 
     21   !! dynamics and tracer fields             !  before  !  now     !  after   ! the after trends becomes the fields 
     22   !! --------------------------             !  fields  !  fields  !  trends  ! only after tra_zdf and dyn_spg 
     23   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   ub     ,  un      ,  ua      !: i-horizontal velocity      [m/s] 
     24   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   vb     ,  vn      ,  va      !: j-horizontal velocity      [m/s] 
     25   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::             wn                 !: vertical velocity          [m/s] 
     26   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   rotb   ,  rotn               !: relative vorticity         [s-1] 
     27   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   hdivb  ,  hdivn              !: horizontal divergence      [s-1] 
     28   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   tb     ,  tn      ,  ta      !: potential temperature      [Celcius] 
     29   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   sb     ,  sn      ,  sa      !: salinity                   [psu] 
     30   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   rn2b   ,  rn2                !: brunt-vaisala frequency**2 [s-2] 
     31   ! 
     32   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   rhd    !: in situ density anomalie rhd=(rho-rau0)/rau0     [no units] 
     33   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   rhop   !: potential volumic mass                           [kg/m3] 
    4334 
    44  
    45       !! advection scheme choice 
    46       !! ----------------------- 
    47       CHARACTER(len=3), PUBLIC  ::   l_adv   !: 'ce2' centre scheme used 
    48          !                                   !: 'tvd' TVD scheme used 
    49          !                                   !: 'mus' MUSCL scheme used 
    50          !                                   !: 'mu2' MUSCL2 scheme used 
     35   !! advection scheme choice 
     36   !! ----------------------- 
     37   CHARACTER(len=3), PUBLIC  ::   l_adv   !: flag for the advection scheme used (= 'ce2', 'tvd', 'mus' or ...) 
    5138 
    5239   !! surface pressure gradient 
    5340   !! ------------------------- 
    54    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &   !: 
    55       spgu, spgv             !: horizontal surface pressure gradient 
     41   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   spgu, spgv      !: horizontal surface pressure gradient 
    5642 
    5743   !! interpolated gradient (only used in zps case) 
    5844   !! --------------------- 
    59    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &   !: 
    60       gtu, gsu, gru,      &  !: t-, s- and rd horizontal gradient at u- and  
    61       gtv, gsv, grv          !: v-points at bottom ocean level  
     45   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   gtu, gsu, gru   !: horizontal gradient of T, S and rd at bottom u-point 
     46   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   gtv, gsv, grv   !: horizontal gradient of T, S and rd at bottom v-point  
    6247 
    6348   !! free surface                       !  before  !  now     !  after   ! 
     
    6853   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sshf_b ,  sshf_n  ,  sshf_a  !: sea surface height at f-point [m] 
    6954 
    70    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sshbb         !: before before & after sea surface height at t-point 
    71    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sshn_f        !: filtered now sea surface height at t-point 
    72  
    7355#if defined key_dynspg_rl   ||   defined key_esopa 
    7456   !! rigid-lid formulation 
    7557   !! --------------------- 
    76    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &   !: 
    77       bsfb, bsfn,         &  !: before, now barotropic streamfunction (m3/s) 
    78       bsfd                   !: now trend of barotropic streamfunction (m3/s2) 
     58   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   bsfb, bsfn   !: before, now barotropic streamfunction (m3/s) 
     59   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   bsfd         !: now trend of barotropic streamfunction (m3/s2) 
    7960#endif 
    8061   !!---------------------------------------------------------------------- 
     62   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2008)  
     63   !! $Id:$  
     64   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     65   !!====================================================================== 
    8166END MODULE oce 
  • branches/dev_004_VVL/NEMO/OPA_SRC/step.F90

    r1388 r1408  
    44   !! Time-stepping    : manager of the ocean, tracer and ice time stepping 
    55   !!====================================================================== 
    6    !! History :        !  91-03  (G. Madec)  Original code 
    7    !!                  !  92-06  (M. Imbard)  add a first output record 
    8    !!                  !  96-04  (G. Madec)  introduction of dynspg 
    9    !!                  !  96-04  (M.A. Foujols)  introduction of passive tracer 
    10    !!             8.0  !  97-06  (G. Madec)  new architecture of call 
    11    !!             8.2  !  97-06  (G. Madec, M. Imbard, G. Roullet)  free surface 
    12    !!             8.2  !  99-02  (G. Madec, N. Grima)  hpg implicit 
    13    !!             8.2  !  00-07  (J-M Molines, M. Imbard)  Open Bondary Conditions 
    14    !!             9.0  !  02-06  (G. Madec)  free form, suppress macro-tasking 
    15    !!             " "  !  04-08  (C. Talandier) New trends organization 
    16    !!             " "  !  05-01  (C. Ethe) Add the KPP closure scheme 
    17    !!             " "  !  05-11  (V. Garnier) Surface pressure gradient organization 
    18    !!             " "  !  05-11  (G. Madec)  Reorganisation of tra and dyn calls 
    19    !!             " "  !  06-01  (L. Debreu, C. Mazauric)  Agrif implementation 
    20    !!             " "  !  06-07  (S. Masson)  restart using iom 
    21    !!             " "  !  06-08  (G. Madec)  surface module  
    22    !!             " "  !  07-07  (J. Chanut, A. Sellar) Unstructured open boundaries (BDY) 
     6   !! History :  OPA  !  1991-03  (G. Madec)  Original code 
     7   !!                 !  1991-11  (G. Madec) 
     8   !!                 !  1992-06  (M. Imbard)  add a first output record 
     9   !!                 !  1996-04  (G. Madec)  introduction of dynspg 
     10   !!                 !  1996-04  (M.A. Foujols)  introduction of passive tracer 
     11   !!            8.0  !  1997-06  (G. Madec)  new architecture of call 
     12   !!            8.2  !  1997-06  (G. Madec, M. Imbard, G. Roullet)  free surface 
     13   !!             -   !  1999-02  (G. Madec, N. Grima)  hpg implicit 
     14   !!             -   !  2000-07  (J-M Molines, M. Imbard)  Open Bondary Conditions 
     15   !!   NEMO     1.0  !  2002-06  (G. Madec)  free form, suppress macro-tasking 
     16   !!             -   !  2004-08  (C. Talandier) New trends organization 
     17   !!                 !  2005-01  (C. Ethe) Add the KPP closure scheme 
     18   !!                 !  2005-11  (G. Madec)  Reorganisation of tra and dyn calls 
     19   !!                 !  2006-01  (L. Debreu, C. Mazauric)  Agrif implementation 
     20   !!                 !  2006-07  (S. Masson)  restart using iom 
     21   !!            3.2  !  2009-02  (G. Madec, R. Benshila)  reintroduicing z*-coordinate 
    2322   !!---------------------------------------------------------------------- 
    2423 
     
    116115   USE restart         ! ocean restart                    (rst_wri routine) 
    117116   USE prtctl          ! Print control                    (prt_ctl routine) 
    118    USE domvvl          ! variable volume                  (dom_vvl routine) 
    119117 
    120118#if defined key_agrif 
     
    140138#if defined key_agrif 
    141139   SUBROUTINE stp( ) 
     140      INTEGER             ::   kstp   ! ocean time-step index 
    142141#else 
    143142   SUBROUTINE stp( kstp ) 
     143      INTEGER, INTENT(in) ::   kstp   ! ocean time-step index 
    144144#endif 
    145145      !!---------------------------------------------------------------------- 
     
    159159      !!              -8- Outputs and diagnostics 
    160160      !!---------------------------------------------------------------------- 
    161       !! * Arguments 
    162 #if defined key_agrif    
    163       INTEGER             ::   kstp   ! ocean time-step index 
    164 #else 
    165       INTEGER, INTENT(in) ::   kstp   ! ocean time-step index 
    166 #endif       
    167161      INTEGER ::   jk       ! dummy loop indice 
    168162      INTEGER ::   indic    ! error indicator if < 0 
     
    202196 
    203197      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    204       ! Computation of diagnostic variables 
    205       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    206       ! N.B. ua, va, ta, sa arrays are used as workspace in this section 
    207       !----------------------------------------------------------------------- 
     198      !                Ocean dynamics : ssh, wn, hdiv, rot                   ! 
     199      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    208200                       CALL div_cur( kstp )                 ! Horizontal divergence & Relative vorticity 
    209201      IF( n_cla == 1 ) CALL div_cla( kstp )                 ! Cross Land Advection (Update Hor. divergence) 
     
    213205      ! Ocean physics update 
    214206      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     207      ! N.B. ua, va, ta, sa arrays are used as workspace in this section 
     208      !----------------------------------------------------------------------- 
    215209#if defined key_zdftke2 
    216210      IF ( ln_dynhpg_imp ) THEN 
    217211      !----------------------------------------------------------------------- 
    218212      !  LATERAL PHYSICS 
    219       !----------------------------------------------------------------------- 
    220       ! N.B. ua, va, ta, sa arrays are used as workspace in this section 
    221213      !----------------------------------------------------------------------- 
    222214                               CALL zdf_mxl( kstp )                 ! mixed layer depth 
     
    229221      !----------------------------------------------------------------------- 
    230222      !  VERTICAL PHYSICS 
    231       !----------------------------------------------------------------------- 
    232       ! N.B. ua, va, ta, sa arrays are used as workspace in this section 
    233223      !----------------------------------------------------------------------- 
    234224      IF( neuler == 0 .AND. kstp == nit000 ) THEN 
     
    273263      !----------------------------------------------------------------------- 
    274264      !  LATERAL PHYSICS 
    275       !----------------------------------------------------------------------- 
    276       ! N.B. ua, va, ta, sa arrays are used as workspace in this section 
    277265      !----------------------------------------------------------------------- 
    278266      IF( lk_ldfslp     )   CALL ldf_slp( kstp, rhd, rn2b )      ! before slope of the lateral mixing 
Note: See TracChangeset for help on using the changeset viewer.