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 5836 for trunk/NEMOGCM/NEMO/OPA_SRC/step.F90 – NEMO

Ignore:
Timestamp:
2015-10-26T15:49:40+01:00 (8 years ago)
Author:
cetlod
Message:

merge the simplification branch onto the trunk, see ticket #1612

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/step.F90

    r5656 r5836  
    2424   !!             -   !  2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase + merge TRC-TRA 
    2525   !!            3.4  !  2011-04  (G. Madec, C. Ethe) Merge of dtatem and dtasal 
    26    !!                 !  2012-07  (J. Simeon, G. Madec, C. Ethe) Online coarsening of outputs 
    27    !!            3.7  !  2014-04  (F. Roquet, G. Madec) New equations of state 
     26   !!            3.6  !  2012-07  (J. Simeon, G. Madec. C. Ethe)  Online coarsening of outputs 
     27   !!            3.6  !  2014-04  (F. Roquet, G. Madec) New equations of state 
     28   !!            3.7  !  2014-10  (G. Madec)  LDF simplication  
     29   !!             -   !  2014-12  (G. Madec) remove KPP scheme 
    2830   !!---------------------------------------------------------------------- 
    2931 
     
    3739   PRIVATE 
    3840 
    39    PUBLIC   stp   ! called by opa.F90 
     41   PUBLIC   stp   ! called by nemogcm.F90 
    4042 
    4143   !! * Substitutions 
     
    4345!!gm   #  include "zdfddm_substitute.h90" 
    4446   !!---------------------------------------------------------------------- 
    45    !! NEMO/OPA 3.7 , NEMO Consortium (2014) 
     47   !! NEMO/OPA 3.7 , NEMO Consortium (2015) 
    4648   !! $Id$ 
    4749   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     
    6971      !!              -5- Compute the momentum trends 
    7072      !!              -6- Update the horizontal velocity 
    71       !!              -7- Compute the diagnostics variables (rd,N2, div,cur,w) 
     73      !!              -7- Compute the diagnostics variables (rd,N2, hdiv,w) 
    7274      !!              -8- Outputs and diagnostics 
    7375      !!---------------------------------------------------------------------- 
     
    7678      INTEGER ::   kcall    ! optional integer argument (dom_vvl_sf_nxt) 
    7779      !! --------------------------------------------------------------------- 
    78  
    7980#if defined key_agrif 
    8081      kstp = nit000 + Agrif_Nb_Step() 
    81       IF ( lk_agrif_debug ) THEN 
    82          IF ( Agrif_Root() .and. lwp) Write(*,*) '---' 
    83          IF (lwp) Write(*,*) 'Grid Number',Agrif_Fixed(),' time step ',kstp, 'int tstep',Agrif_NbStepint() 
    84       ENDIF 
    85  
    86       IF ( kstp == (nit000 + 1) ) lk_agrif_fstep = .FALSE. 
    87  
     82      IF( lk_agrif_debug ) THEN 
     83         IF( Agrif_Root() .and. lwp)   WRITE(*,*) '---' 
     84         IF(lwp)   WRITE(*,*) 'Grid Number', Agrif_Fixed(),' time step ', kstp, 'int tstep', Agrif_NbStepint() 
     85      ENDIF 
     86      IF( kstp == nit000 + 1 )   lk_agrif_fstep = .FALSE. 
    8887# if defined key_iomput 
    8988      IF( Agrif_Nbstepint() == 0 )   CALL iom_swap( cxios_context ) 
    9089# endif 
    9190#endif 
    92                              indic = 0           ! reset to no error condition 
    93       IF( kstp == nit000 ) THEN 
    94          ! must be done after nemo_init for AGRIF+XIOS+OASIS 
    95                       CALL iom_init(      cxios_context          )  ! iom_put initialization 
    96          IF( ln_crs ) CALL iom_init( TRIM(cxios_context)//"_crs" )  ! initialize context for coarse grid 
    97       ENDIF 
    98  
     91      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     92      ! update I/O and calendar  
     93      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     94                             indic = 0                ! reset to no error condition 
     95                              
     96      IF( kstp == nit000 ) THEN                       ! initialize IOM context (must be done after nemo_init for AGRIF+XIOS+OASIS) 
     97                             CALL iom_init(      cxios_context          )  ! for model grid (including passible AGRIF zoom) 
     98         IF( ln_crs      )   CALL iom_init( TRIM(cxios_context)//"_crs" )  ! for coarse grid 
     99      ENDIF 
    99100      IF( kstp /= nit000 )   CALL day( kstp )         ! Calendar (day was already called at nit000 in day_init) 
    100                              CALL iom_setkt( kstp - nit000 + 1,      cxios_context          )   ! tell iom we are at time step kstp 
    101       IF( ln_crs     )       CALL iom_setkt( kstp - nit000 + 1, TRIM(cxios_context)//"_crs" )   ! tell iom we are at time step kstp 
    102  
    103       !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    104       ! Update data, open boundaries, surface boundary condition (including sea-ice) 
    105       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    106       IF( lk_tide    )   CALL sbc_tide( kstp ) 
    107       IF( lk_bdy     )  THEN 
    108          IF( ln_apr_dyn) CALL sbc_apr( kstp )   ! bdy_dta needs ssh_ib  
    109                          CALL bdy_dta ( kstp, time_offset=+1 )   ! update dynamic & tracer data at open boundaries 
    110       ENDIF 
    111                          CALL sbc    ( kstp )         ! Sea Boundary Condition (including sea-ice) 
    112                                                       ! clem: moved here for bdy ice purpose 
     101                             CALL iom_setkt( kstp - nit000 + 1,      cxios_context          )   ! tell IOM we are at time step kstp 
     102      IF( ln_crs         )   CALL iom_setkt( kstp - nit000 + 1, TRIM(cxios_context)//"_crs" )   ! tell IOM we are at time step kstp 
     103 
     104      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     105      ! Update external forcing (tides, open boundaries, and surface boundary condition (including sea-ice) 
     106      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     107      IF( lk_tide    )   CALL sbc_tide( kstp )                   ! update tide potential 
     108      IF( ln_apr_dyn )   CALL sbc_apr ( kstp )                   ! atmospheric pressure (NB: call before bdy_dta which needs ssh_ib)  
     109      IF( lk_bdy     )   CALL bdy_dta ( kstp, time_offset=+1 )   ! update dynamic & tracer data at open boundaries 
     110                         CALL sbc     ( kstp )                   ! Sea Boundary Condition (including sea-ice) 
     111 
    113112      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    114113      ! Update stochastic parameters and random T/S fluctuations 
    115114      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    116                         CALL sto_par( kstp )          ! Stochastic parameters 
     115                         CALL sto_par( kstp )          ! Stochastic parameters 
    117116 
    118117      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     
    124123                         CALL bn2    ( tsb, rab_b, rn2b ) ! before Brunt-Vaisala frequency 
    125124                         CALL bn2    ( tsn, rab_n, rn2  ) ! now    Brunt-Vaisala frequency 
     125 
    126126      ! 
    127127      !  VERTICAL PHYSICS 
     
    131131      IF( lk_zdftke  )   CALL zdf_tke( kstp )            ! TKE closure scheme for Kz 
    132132      IF( lk_zdfgls  )   CALL zdf_gls( kstp )            ! GLS closure scheme for Kz 
    133       IF( lk_zdfkpp  )   CALL zdf_kpp( kstp )            ! KPP closure scheme for Kz 
    134133      IF( lk_zdfcst  ) THEN                              ! Constant Kz (reset avt, avm[uv] to the background value) 
    135134         avt (:,:,:) = rn_avt0 * wmask (:,:,:) 
     
    137136         avmv(:,:,:) = rn_avm0 * wvmask(:,:,:) 
    138137      ENDIF 
     138 
    139139      IF( ln_rnf_mouth ) THEN                         ! increase diffusivity at rivers mouths 
    140          DO jk = 2, nkrnf   ;   avt(:,:,jk) = avt(:,:,jk) + 2.e0 * rn_avt_rnf * rnfmsk(:,:) * tmask(:,:,jk)   ;   END DO 
     140         DO jk = 2, nkrnf   ;   avt(:,:,jk) = avt(:,:,jk) + 2._wp * rn_avt_rnf * rnfmsk(:,:) * tmask(:,:,jk)   ;   END DO 
    141141      ENDIF 
    142142      IF( ln_zdfevd  )   CALL zdf_evd( kstp )         ! enhanced vertical eddy diffusivity 
     
    144144      IF( lk_zdftmx  )   CALL zdf_tmx( kstp )         ! tidal vertical mixing 
    145145 
    146       IF( lk_zdfddm .AND. .NOT. lk_zdfkpp )   & 
    147          &               CALL zdf_ddm( kstp )         ! double diffusive mixing 
     146      IF( lk_zdfddm  )   CALL zdf_ddm( kstp )         ! double diffusive mixing 
    148147 
    149148                         CALL zdf_mxl( kstp )         ! mixed layer depth 
     
    155154      !  LATERAL  PHYSICS 
    156155      ! 
    157       IF( lk_ldfslp ) THEN                            ! slope of lateral mixing 
     156      IF( l_ldfslp ) THEN                             ! slope of lateral mixing 
     157!!gm : why this here ???? 
    158158         IF(ln_sto_eos ) CALL sto_pts( tsn )          ! Random T/S fluctuations 
     159!!gm 
    159160                         CALL eos( tsb, rhd, gdept_0(:,:,:) )               ! before in situ density 
     161 
    160162         IF( ln_zps .AND. .NOT. ln_isfcav)                               & 
    161163            &            CALL zps_hde    ( kstp, jpts, tsb, gtsu, gtsv,  &  ! Partial steps: before horizontal gradient 
    162164            &                                          rhd, gru , grv    )  ! of t, s, rd at the last ocean level 
     165 
    163166         IF( ln_zps .AND.       ln_isfcav)                               & 
    164             &            CALL zps_hde_isf( kstp, jpts, tsb, gtsu, gtsv, &    ! Partial steps for top cell (ISF) 
     167            &            CALL zps_hde_isf( kstp, jpts, tsb, gtsu, gtsv, gtui, gtvi,  &    ! Partial steps for top cell (ISF) 
    165168            &                                          rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv ,   & 
    166             &                                   gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) ! of t, s, rd at the first ocean level 
    167          IF( ln_traldf_grif ) THEN                           ! before slope for Griffies operator 
    168                          CALL ldf_slp_grif( kstp ) 
    169          ELSE 
    170                          CALL ldf_slp( kstp, rhd, rn2b )     ! before slope for Madec operator 
     169            &                                   grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) ! of t, s, rd at the first ocean level 
     170 
     171         IF( ln_traldf_triad ) THEN  
     172                         CALL ldf_slp_triad( kstp )                       ! before slope for triad operator 
     173         ELSE      
     174                         CALL ldf_slp     ( kstp, rhd, rn2b )             ! before slope for standard operator 
    171175         ENDIF 
    172176      ENDIF 
    173 #if defined key_traldf_c2d 
    174       IF( lk_traldf_eiv )   CALL ldf_eiv( kstp )      ! eddy induced velocity coefficient 
    175 #endif 
    176 #if defined key_traldf_c3d && defined key_traldf_smag 
    177                           CALL ldf_tra_smag( kstp )      ! eddy induced velocity coefficient 
    178 #  endif 
    179 #if defined key_dynldf_c3d && defined key_dynldf_smag 
    180                           CALL ldf_dyn_smag( kstp )      ! eddy induced velocity coefficient 
    181 #  endif 
    182  
    183       !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    184       !  Ocean dynamics : hdiv, rot, ssh, e3, wn 
    185       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    186                          CALL ssh_nxt       ( kstp )  ! after ssh (includes call to div_cur) 
     177      !                                               ! eddy diffusivity coeff. and/or eiv coeff. 
     178      IF( l_ldftra_time .OR. l_ldfeiv_time )   CALL ldf_tra( kstp )  
     179 
     180      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     181      !  Ocean dynamics : hdiv, ssh, e3, wn 
     182      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     183                         CALL ssh_nxt       ( kstp )  ! after ssh (includes call to div_hor) 
    187184      IF( lk_vvl     )   CALL dom_vvl_sf_nxt( kstp )  ! after vertical scale factors  
    188185                         CALL wzv           ( kstp )  ! now cross-level velocity  
     
    192189          ! Note that the computation of vertical velocity above, hence "after" sea level 
    193190          ! is necessary to compute momentum advection for the rhs of barotropic loop: 
     191!!gm : why also here ???? 
    194192            IF(ln_sto_eos ) CALL sto_pts( tsn )                             ! Random T/S fluctuations 
     193!!gm 
    195194                            CALL eos    ( tsn, rhd, rhop, fsdept_n(:,:,:) ) ! now in situ density for hpg computation 
    196             IF( ln_zps .AND. .NOT. ln_isfcav)                               & 
    197                &            CALL zps_hde    ( kstp, jpts, tsn, gtsu, gtsv,  &    ! Partial steps: before horizontal gradient 
    198                &                                          rhd, gru , grv    )  ! of t, s, rd at the last ocean level 
    199             IF( ln_zps .AND.       ln_isfcav)                               & 
    200                &            CALL zps_hde_isf( kstp, jpts, tsn, gtsu, gtsv,  &    ! Partial steps for top cell (ISF) 
     195                             
     196            IF( ln_zps .AND. .NOT. ln_isfcav)   &                           ! Partial steps: bottom before horizontal gradient 
     197               &            CALL zps_hde    ( kstp, jpts, tsn, gtsu, gtsv,  &  ! of t, s, rd at the last ocean level 
     198               &                                          rhd, gru , grv    ) 
     199            IF( ln_zps .AND.       ln_isfcav)   &                           ! Partial steps: top & bottom before horizontal gradient 
     200               &            CALL zps_hde_isf( kstp, jpts, tsn, gtsu, gtsv, gtui, gtvi,   &  
    201201               &                                          rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv ,   & 
    202                &                                   gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) ! of t, s, rd at the last ocean level 
    203  
    204                                   ua(:,:,:) = 0.e0             ! set dynamics trends to zero 
    205                                   va(:,:,:) = 0.e0 
     202               &                                               grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) 
     203 
     204                                  ua(:,:,:) = 0._wp            ! set dynamics trends to zero 
     205                                  va(:,:,:) = 0._wp 
    206206          IF(  lk_asminc .AND. ln_asmiau .AND. & 
    207207             & ln_dyninc       )  CALL dyn_asm_inc  ( kstp )   ! apply dynamics assimilation increment 
    208           IF( ln_neptsimp )       CALL dyn_nept_cor ( kstp )   ! subtract Neptune velocities (simplified) 
    209208          IF( lk_bdy           )  CALL bdy_dyn3d_dmp( kstp )   ! bdy damping trends 
    210209                                  CALL dyn_adv      ( kstp )   ! advection (vector or flux form) 
    211210                                  CALL dyn_vor      ( kstp )   ! vorticity term including Coriolis 
    212211                                  CALL dyn_ldf      ( kstp )   ! lateral mixing 
    213           IF( ln_neptsimp )       CALL dyn_nept_cor ( kstp )   ! add Neptune velocities (simplified) 
    214212#if defined key_agrif 
    215213          IF(.NOT. Agrif_Root())  CALL Agrif_Sponge_dyn        ! momentum sponge 
     
    221219                                  va_sv(:,:,:) = va(:,:,:) 
    222220 
    223                                   CALL div_cur( kstp )         ! Horizontal divergence & Relative vorticity (2nd call in time-split case) 
     221                                  CALL div_hor( kstp )         ! Horizontal divergence (2nd call in time-split case) 
    224222          IF( lk_vvl     )        CALL dom_vvl_sf_nxt( kstp, kcall=2 )  ! after vertical scale factors (update depth average component) 
    225223                                  CALL wzv           ( kstp )  ! now cross-level velocity  
     
    229227      ! diagnostics and outputs             (ua, va, tsa used as workspace) 
    230228      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    231       IF( lk_floats  )      CALL flo_stp( kstp )         ! drifting Floats 
    232       IF( lk_diahth  )      CALL dia_hth( kstp )         ! Thermocline depth (20 degres isotherm depth) 
    233       IF( .NOT. ln_cpl )    CALL dia_fwb( kstp )         ! Fresh water budget diagnostics 
    234       IF( lk_diadct  )      CALL dia_dct( kstp )         ! Transports 
    235       IF( lk_diaar5  )      CALL dia_ar5( kstp )         ! ar5 diag 
    236       IF( lk_diaharm )      CALL dia_harm( kstp )        ! Tidal harmonic analysis 
    237                             CALL dia_wri( kstp )         ! ocean model: outputs 
    238       ! 
    239       IF( ln_crs     )      CALL crs_fld( kstp )         ! ocean model: online field coarsening & output 
     229      IF( lk_floats  )   CALL flo_stp( kstp )         ! drifting Floats 
     230      IF( lk_diahth  )   CALL dia_hth( kstp )         ! Thermocline depth (20 degres isotherm depth) 
     231      IF(.NOT.ln_cpl )   CALL dia_fwb( kstp )         ! Fresh water budget diagnostics 
     232      IF( lk_diadct  )   CALL dia_dct( kstp )         ! Transports 
     233      IF( lk_diaar5  )   CALL dia_ar5( kstp )         ! ar5 diag 
     234      IF( lk_diaharm )   CALL dia_harm( kstp )        ! Tidal harmonic analysis 
     235                         CALL dia_wri( kstp )         ! ocean model: outputs 
     236      ! 
     237      IF( ln_crs     )   CALL crs_fld( kstp )         ! ocean model: online field coarsening & output 
    240238 
    241239#if defined key_top 
     
    246244#endif 
    247245 
    248  
    249246      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    250247      ! Active tracers                              (ua, va used as workspace) 
    251248      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    252                              tsa(:,:,:,:) = 0.e0            ! set tracer trends to zero 
     249                             tsa(:,:,:,:) = 0._wp           ! set tracer trends to zero 
    253250 
    254251      IF(  lk_asminc .AND. ln_asmiau .AND. & 
     
    261258      IF( lk_bdy         )   CALL bdy_tra_dmp( kstp )       ! bdy damping trends 
    262259                             CALL tra_adv    ( kstp )       ! horizontal & vertical advection 
    263       IF( lk_zdfkpp      )   CALL tra_kpp    ( kstp )       ! KPP non-local tracer fluxes 
    264260                             CALL tra_ldf    ( kstp )       ! lateral mixing 
    265261 
     262!!gm : why CALL to dia_ptr has been moved here??? (use trends info?) 
    266263      IF( ln_diaptr      )   CALL dia_ptr                   ! Poleward adv/ldf TRansports diagnostics 
     264!!gm 
    267265 
    268266#if defined key_agrif 
     
    274272         IF( ln_zdfnpc   )   CALL tra_npc( kstp )                ! update after fields by non-penetrative convection 
    275273                             CALL tra_nxt( kstp )                ! tracer fields at next time step 
     274!!gm : why again a call to sto_pts ??? 
    276275            IF( ln_sto_eos ) CALL sto_pts( tsn )                 ! Random T/S fluctuations 
     276!!gm 
    277277                             CALL eos    ( tsa, rhd, rhop, fsdept_n(:,:,:) )  ! Time-filtered in situ density for hpg computation 
    278278            IF( ln_zps .AND. .NOT. ln_isfcav)                                & 
     
    280280               &                                           rhd, gru , grv    )  ! of t, s, rd at the last ocean level 
    281281            IF( ln_zps .AND.       ln_isfcav)                                & 
    282                &             CALL zps_hde_isf( kstp, jpts, tsa, gtsu, gtsv,  &    ! Partial steps for top cell (ISF) 
     282               &             CALL zps_hde_isf( kstp, jpts, tsa, gtsu, gtsv, gtui, gtvi,  &    ! Partial steps for top/bottom cells 
    283283               &                                           rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv ,   & 
    284                &                                    gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) ! of t, s, rd at the last ocean level 
     284               &                                                grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) 
    285285      ELSE                                                  ! centered hpg  (eos then time stepping) 
    286286         IF ( .NOT. lk_dynspg_ts ) THEN                     ! eos already called in time-split case 
     287!!gm : why again a call to sto_pts ??? 
    287288            IF( ln_sto_eos ) CALL sto_pts( tsn )    ! Random T/S fluctuations 
     289!!gm 
    288290                             CALL eos    ( tsn, rhd, rhop, fsdept_n(:,:,:) )  ! now in situ density for hpg computation 
    289291         IF( ln_zps .AND. .NOT. ln_isfcav)                                   & 
    290                &             CALL zps_hde    ( kstp, jpts, tsn, gtsu, gtsv,  &    ! Partial steps: before horizontal gradient 
    291                &                                           rhd, gru , grv    )  ! of t, s, rd at the last ocean level 
     292               &             CALL zps_hde    ( kstp, jpts, tsn, gtsu, gtsv,  &    ! Partial steps: bottom before horizontal gradient 
     293               &                                           rhd, gru , grv    )    ! of t, s, rd at the last ocean level 
    292294         IF( ln_zps .AND.       ln_isfcav)                                   &  
    293                &             CALL zps_hde_isf( kstp, jpts, tsn, gtsu, gtsv,  &    ! Partial steps for top cell (ISF) 
     295               &             CALL zps_hde_isf( kstp, jpts, tsn, gtsu, gtsv, gtui, gtvi,   &    ! Partial steps for top/bottom cells 
    294296               &                                           rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv ,   & 
    295                &                                    gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) ! of t, s, rd at the last ocean level 
     297               &                                    grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) ! of t, s, rd at the last ocean level 
    296298         ENDIF 
    297299         IF( ln_zdfnpc   )   CALL tra_npc( kstp )                ! update after fields by non-penetrative convection 
     
    307309                               ua(:,:,:) = ua_sv(:,:,:) 
    308310                               va(:,:,:) = va_sv(:,:,:) 
    309                                                              ! Revert now divergence and rotational to previously computed ones  
    310                                                              !(needed because of the time swap in div_cur, at the beginning of each time step) 
    311                                hdivn(:,:,:) = hdivb(:,:,:) 
    312                                rotn(:,:,:)  = rotb(:,:,:)  
    313311 
    314312                               CALL dyn_bfr( kstp )         ! bottom friction 
    315313                               CALL dyn_zdf( kstp )         ! vertical diffusion 
    316314      ELSE 
    317                                ua(:,:,:) = 0.e0             ! set dynamics trends to zero 
    318                                va(:,:,:) = 0.e0 
     315                               ua(:,:,:) = 0._wp            ! set dynamics trends to zero 
     316                               va(:,:,:) = 0._wp 
    319317 
    320318        IF(  lk_asminc .AND. ln_asmiau .AND. & 
    321319           & ln_dyninc      )  CALL dyn_asm_inc( kstp )     ! apply dynamics assimilation increment 
    322320        IF( ln_bkgwri )        CALL asm_bkg_wri( kstp )     ! output background fields 
    323         IF( ln_neptsimp )      CALL dyn_nept_cor( kstp )    ! subtract Neptune velocities (simplified) 
    324321        IF( lk_bdy          )  CALL bdy_dyn3d_dmp(kstp )    ! bdy damping trends 
    325322                               CALL dyn_adv( kstp )         ! advection (vector or flux form) 
    326323                               CALL dyn_vor( kstp )         ! vorticity term including Coriolis 
    327324                               CALL dyn_ldf( kstp )         ! lateral mixing 
    328         IF( ln_neptsimp )      CALL dyn_nept_cor( kstp )    ! add Neptune velocities (simplified) 
    329325#if defined key_agrif 
    330326        IF(.NOT. Agrif_Root()) CALL Agrif_Sponge_dyn        ! momemtum sponge 
     
    340336      IF( lk_vvl           )   CALL dom_vvl_sf_swp( kstp )  ! swap of vertical scale factors 
    341337      ! 
     338 
     339!!gm : This does not only concern the dynamics ==>>> add a new title 
     340!!gm2: why ouput restart before AGRIF update? 
    342341      IF( lrst_oce         )   CALL rst_write( kstp )       ! write output ocean restart file 
    343342 
     
    367366                 CALL iom_close( numror )     ! close input  ocean restart file 
    368367         IF(lwm) CALL FLUSH    ( numond )     ! flush output namelist oce 
    369          IF( lwm.AND.numoni /= -1 ) CALL FLUSH    ( numoni )     ! flush output namelist ice 
     368         IF(lwm.AND.numoni /= -1 )   & 
     369            &    CALL FLUSH    ( numoni )     ! flush output namelist ice (if exist) 
    370370      ENDIF 
    371371 
     
    373373      ! Coupled mode 
    374374      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     375!!gm why lk_oasis and not lk_cpl ???? 
    375376      IF( lk_oasis         )   CALL sbc_cpl_snd( kstp )     ! coupled mode : field exchanges 
    376377      ! 
     
    383384      ! 
    384385      IF( nn_timing == 1 .AND.  kstp == nit000  )   CALL timing_reset 
    385       !      
    386386      ! 
    387387   END SUBROUTINE stp 
Note: See TracChangeset for help on using the changeset viewer.