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 11394 – NEMO

Changeset 11394


Ignore:
Timestamp:
2019-08-02T15:14:02+02:00 (5 years ago)
Author:
mattmartin
Message:

First implementation of STOPACK in the GO6 package branch.

Location:
branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC
Files:
18 edited
1 copied

Legend:

Unmodified
Added
Removed

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf.F90

    r6486 r11394  
    3030   USE wrk_nemo        ! Memory Allocation 
    3131   USE timing          ! Timing 
     32   USE stopack 
    3233 
    3334   IMPLICIT NONE 
    … …  
    3839 
    3940   INTEGER ::   nldf = -2   ! type of lateral diffusion used defined from ln_dynldf_... namlist logicals) 
     41 
     42#if defined key_dynldf_c3d 
     43   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ahm10,ahm20,ahm30,ahm40 
     44#endif 
     45#if defined key_dynldf_c2d 
     46   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:  ) :: ahm10,ahm20,ahm30,ahm40 
     47#endif 
    4048 
    4149   !! * Substitutions 
    … …  
    6775         ztrdv(:,:,:) = va(:,:,:)  
    6876      ENDIF 
     77 
     78#if defined key_dynldf_c3d 
     79         IF( kt .eq. nit000 .and. (nn_spp_ahm1+nn_spp_ahm2) .gt. 0 ) THEN 
     80             ALLOCATE ( ahm10(jpi,jpj,jpk), ahm20(jpi,jpj,jpk) ) 
     81             ALLOCATE ( ahm30(jpi,jpj,jpk), ahm40(jpi,jpj,jpk) ) 
     82             ahm10 = ahm1 
     83             ahm20 = ahm2 
     84             ahm30 = ahm3 
     85             ahm40 = ahm4 
     86         ENDIF 
     87#endif 
     88#if defined key_dynldf_c2d 
     89         IF( kt .eq. nit000 .and. (nn_spp_ahm1+nn_spp_ahm2) .gt. 0 ) THEN 
     90             ALLOCATE ( ahm10(jpi,jpj), ahm20(jpi,jpj) ) 
     91             ALLOCATE ( ahm30(jpi,jpj), ahm40(jpi,jpj) ) 
     92             ahm10 = ahm1 
     93             ahm20 = ahm2 
     94             ahm30 = ahm3 
     95             ahm40 = ahm4 
     96         ENDIF 
     97#endif 
     98 
     99#if defined key_traldf_c3d || defined key_traldf_c2d 
     100         IF( nn_spp_ahm1 .GT. 0) THEN 
     101           IF( ln_dynldf_lap ) THEN 
     102             ahm1 = ahm10 
     103             CALL spp_ahm(kt,ahm1,nn_spp_ahm1,rn_ahm1_sd,jk_spp_ahm1) 
     104           ENDIF 
     105           IF( ln_dynldf_bilap ) THEN 
     106             ahm3 = ahm30 
     107             CALL spp_ahm(kt,ahm3,nn_spp_ahm1,rn_ahm1_sd,jk_spp_ahm3) 
     108           ENDIF 
     109         ENDIF 
     110         IF( nn_spp_ahm2 .GT. 0) THEN 
     111           IF( ln_dynldf_lap ) THEN 
     112             ahm2 = ahm20 
     113             CALL spp_ahm(kt,ahm2,nn_spp_ahm2,rn_ahm2_sd,jk_spp_ahm2) 
     114           ENDIF 
     115           IF( ln_dynldf_bilap ) THEN 
     116             ahm4 = ahm40 
     117             CALL spp_ahm(kt,ahm4,nn_spp_ahm2,rn_ahm2_sd,jk_spp_ahm4) 
     118           ENDIF 
     119         ENDIF 
     120#endif 
    69121 
    70122      SELECT CASE ( nldf )                       ! compute lateral mixing trend and add it to the general trend 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/SBC/albedo.F90

    r6498 r11394  
    2222   USE wrk_nemo       ! work arrays 
    2323   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
     24   USE stopack 
    2425 
    2526   IMPLICIT NONE 
    … …  
    154155               END DO 
    155156            END DO 
     157             
     158            IF ( nn_spp_icealb > 0 ) CALL spp_gen( 1, pa_ice_cs(:,:,jl), nn_spp_icealb, rn_icealb_sd, jk_spp_alb, jl ) 
     159                         
    156160         END DO 
    157161 
    … …  
    209213              END DO 
    210214            END DO 
     215             
     216            IF ( nn_spp_icealb > 0 ) CALL spp_gen( 1, pa_ice_os(:,:,jl), nn_spp_icealb, rn_icealb_sd, jk_spp_alb, jl ) 
     217             
    211218         END DO 
    212219         ! Effect of the clouds (2d order polynomial) 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90

    r9288 r11394  
    5151   USE par_ice_2 
    5252#endif 
     53   USE stopack 
    5354 
    5455   IMPLICIT NONE 
    … …  
    8990   REAL(wp) ::   rn_efac     ! multiplication factor for evaporation (clem) 
    9091   REAL(wp) ::   rn_vfac     ! multiplication factor for ice/ocean velocity in the calculation of wind stress (clem) 
     92   REAL(wp), ALLOCATABLE, SAVE ::   rn_vfac0(:,:) ! multiplication factor for ice/ocean velocity in the calculation of wind stress (clem) 
    9193   REAL(wp) ::   rn_zqt      ! z(q,t) : height of humidity and temperature measurements 
    9294   REAL(wp) ::   rn_zu       ! z(u)   : height of wind measurements 
    … …  
    198200         sfx(:,:) = 0._wp                          ! salt flux; zero unless ice is present (computed in limsbc(_2).F90) 
    199201         ! 
     202         ALLOCATE ( rn_vfac0(jpi,jpj) ) 
     203         rn_vfac0(:,:) = rn_vfac 
     204         ! 
     205      ENDIF 
     206 
     207      IF( nn_spp_relw > 0 ) THEN 
     208         rn_vfac0 = rn_vfac 
     209         CALL spp_gen(kt, rn_vfac0, nn_spp_relw, rn_relw_sd, jk_spp_relw ) 
    200210      ENDIF 
    201211 
    … …  
    291301      DO jj = 2, jpjm1 
    292302         DO ji = fs_2, fs_jpim1   ! vect. opt. 
    293             zwnd_i(ji,jj) = (  sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( pu(ji-1,jj  ) + pu(ji,jj) )  ) 
    294             zwnd_j(ji,jj) = (  sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( pv(ji  ,jj-1) + pv(ji,jj) )  ) 
     303            zwnd_i(ji,jj) = (  sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac0(ji,jj) * 0.5 * ( pu(ji-1,jj  ) + pu(ji,jj) )  ) 
     304            zwnd_j(ji,jj) = (  sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac0(ji,jj) * 0.5 * ( pv(ji  ,jj-1) + pv(ji,jj) )  ) 
    295305         END DO 
    296306      END DO 
    … …  
    467477            DO ji = 2, jpim1   ! B grid : NO vector opt 
    468478               ! ... scalar wind at I-point (fld being at T-point) 
    469                zwndi_f = 0.25 * (  sf(jp_wndi)%fnow(ji-1,jj  ,1) + sf(jp_wndi)%fnow(ji  ,jj  ,1)   & 
    470                   &              + sf(jp_wndi)%fnow(ji-1,jj-1,1) + sf(jp_wndi)%fnow(ji  ,jj-1,1)  ) - rn_vfac * u_ice(ji,jj) 
    471                zwndj_f = 0.25 * (  sf(jp_wndj)%fnow(ji-1,jj  ,1) + sf(jp_wndj)%fnow(ji  ,jj  ,1)   & 
    472                   &              + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji  ,jj-1,1)  ) - rn_vfac * v_ice(ji,jj) 
     479               zwndi_f = 0.25 * (  sf(jp_wndi)%fnow(ji-1,jj  ,1) + sf(jp_wndi)%fnow(ji  ,jj  ,1)    & 
     480                  &              + sf(jp_wndi)%fnow(ji-1,jj-1,1) + sf(jp_wndi)%fnow(ji  ,jj-1,1)  ) & 
     481                  &      - rn_vfac0(ji,jj) * u_ice(ji,jj) 
     482               zwndj_f = 0.25 * (  sf(jp_wndj)%fnow(ji-1,jj  ,1) + sf(jp_wndj)%fnow(ji  ,jj  ,1)    & 
     483                  &              + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji  ,jj-1,1)  ) & 
     484                  &      - rn_vfac0(ji,jj) * v_ice(ji,jj) 
    473485               zwnorm_f = zcoef_wnorm * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f ) 
    474486               ! ... ice stress at I-point 
    … …  
    476488               vtau_ice(ji,jj) = zwnorm_f * zwndj_f 
    477489               ! ... scalar wind at T-point (fld being at T-point) 
    478                zwndi_t = sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.25 * (  u_ice(ji,jj+1) + u_ice(ji+1,jj+1)   & 
    479                   &                                                    + u_ice(ji,jj  ) + u_ice(ji+1,jj  )  ) 
    480                zwndj_t = sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.25 * (  v_ice(ji,jj+1) + v_ice(ji+1,jj+1)   & 
    481                   &                                                    + v_ice(ji,jj  ) + v_ice(ji+1,jj  )  ) 
     490               zwndi_t = sf(jp_wndi)%fnow(ji,jj,1)                                         & 
     491                  &      - rn_vfac0(ji,jj) * 0.25 * (  u_ice(ji,jj+1) + u_ice(ji+1,jj+1)   & 
     492                  &                                  + u_ice(ji,jj  ) + u_ice(ji+1,jj  )  ) 
     493               zwndj_t = sf(jp_wndj)%fnow(ji,jj,1)                                         & 
     494                  &      - rn_vfac0(ji,jj) * 0.25 * (  v_ice(ji,jj+1) + v_ice(ji+1,jj+1)   & 
     495                  &                                  + v_ice(ji,jj  ) + v_ice(ji+1,jj  )  ) 
    482496               wndm_ice(ji,jj)  = SQRT( zwndi_t * zwndi_t + zwndj_t * zwndj_t ) * tmask(ji,jj,1) 
    483497            END DO 
    … …  
    490504         DO jj = 2, jpj 
    491505            DO ji = fs_2, jpi   ! vect. opt. 
    492                zwndi_t = (  sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( u_ice(ji-1,jj  ) + u_ice(ji,jj) )  ) 
    493                zwndj_t = (  sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( v_ice(ji  ,jj-1) + v_ice(ji,jj) )  ) 
     506               zwndi_t = (  sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac0(ji,jj) * 0.5 * ( u_ice(ji-1,jj  ) + u_ice(ji,jj) )  ) 
     507               zwndj_t = (  sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac0(ji,jj) * 0.5 * ( v_ice(ji  ,jj-1) + v_ice(ji,jj) )  ) 
    494508               wndm_ice(ji,jj) = SQRT( zwndi_t * zwndi_t + zwndj_t * zwndj_t ) * tmask(ji,jj,1) 
    495509            END DO 
    … …  
    497511         DO jj = 2, jpjm1 
    498512            DO ji = fs_2, fs_jpim1   ! vect. opt. 
    499                utau_ice(ji,jj) = zcoef_wnorm2 * ( wndm_ice(ji+1,jj  ) + wndm_ice(ji,jj) )                          & 
    500                   &          * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) - rn_vfac * u_ice(ji,jj) ) 
    501                vtau_ice(ji,jj) = zcoef_wnorm2 * ( wndm_ice(ji,jj+1  ) + wndm_ice(ji,jj) )                          & 
    502                   &          * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) - rn_vfac * v_ice(ji,jj) ) 
     513               utau_ice(ji,jj) = zcoef_wnorm2 * ( wndm_ice(ji+1,jj  ) + wndm_ice(ji,jj) )             & 
     514                  &              * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) & 
     515                  &                  - rn_vfac0(ji,jj) * u_ice(ji,jj) ) 
     516               vtau_ice(ji,jj) = zcoef_wnorm2 * ( wndm_ice(ji,jj+1  ) + wndm_ice(ji,jj) )             & 
     517                  &              * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) & 
     518                  &                  - rn_vfac0(ji,jj) * v_ice(ji,jj) ) 
    503519            END DO 
    504520         END DO 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/SBC/sbcrnf.F90

    r8302 r11394  
    5353   REAL(wp)          , PUBLIC ::   rn_avt_rnf      !: runoffs, value of the additional vertical mixing coef. [m2/s] 
    5454   REAL(wp)          , PUBLIC ::   rn_rfact        !: multiplicative factor for runoff 
     55   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   rn_avt_rnf0 
    5556 
    5657   LOGICAL           , PUBLIC ::   l_rnfcpl = .false.       ! runoffs recieved from oasis 
    … …  
    447448         !                                      !    - mixed upstream-centered (ln_traadv_cen2=T) 
    448449         ! 
     450         ALLOCATE ( rn_avt_rnf0(jpi,jpj) ) 
     451         rn_avt_rnf0(:,:) = rn_avt_rnf 
     452         ! 
    449453         IF ( ln_rnf_depth )   CALL ctl_warn( 'sbc_rnf_init: increased mixing turned on but effects may already',   & 
    450454            &                                              'be spread through depth by ln_rnf_depth'               ) 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/SBC/sbcssr.F90

    r6486 r11394  
    2525   USE timing         ! Timing 
    2626   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
     27   USE stopack 
     28   USE wrk_nemo       ! Memory Allocation 
    2729 
    2830   IMPLICIT NONE 
    … …  
    7577      REAL(wp) ::   zerp     ! local scalar for evaporation damping 
    7678      REAL(wp) ::   zqrp     ! local scalar for heat flux damping 
    77       REAL(wp) ::   zsrp     ! local scalar for unit conversion of rn_deds factor 
    7879      REAL(wp) ::   zerp_bnd ! local scalar for unit conversion of rn_epr_max factor 
     80      REAL(wp), POINTER, DIMENSION(:,:) :: rn_dqdt_s, zsrp 
    7981      INTEGER  ::   ierror   ! return error code 
    8082      !! 
    … …  
    9597            ! 
    9698            IF( nn_sstr == 1 ) THEN                                   !* Temperature restoring term 
     99 
     100               CALL wrk_alloc( jpi, jpj, rn_dqdt_s) 
     101               rn_dqdt_s=rn_dqdt 
     102 
     103               IF( nn_spp_dqdt > 0 ) CALL spp_gen(kt, rn_dqdt_s,nn_spp_dqdt,rn_dqdt_sd,jk_spp_dqdt ) 
    97104               DO jj = 1, jpj 
    98105                  DO ji = 1, jpi 
    99                      zqrp = rn_dqdt * ( sst_m(ji,jj) - sf_sst(1)%fnow(ji,jj,1) ) 
     106                     zqrp = rn_dqdt_s(ji,jj) * ( sst_m(ji,jj) - sf_sst(1)%fnow(ji,jj,1) ) 
    100107                     qns(ji,jj) = qns(ji,jj) + zqrp 
    101108                     qrp(ji,jj) = zqrp 
    … …  
    103110               END DO 
    104111               CALL iom_put( "qrp", qrp )                             ! heat flux damping 
     112               CALL wrk_dealloc( jpi, jpj, rn_dqdt_s ) 
    105113            ENDIF 
    106114            ! 
    107115            IF( nn_sssr == 1 ) THEN                                   !* Salinity damping term (salt flux only (sfx)) 
    108                zsrp = rn_deds / rday                                  ! from [mm/day] to [kg/m2/s] 
     116               CALL wrk_alloc( jpi, jpj, zsrp) 
     117               zsrp = rn_deds 
     118               IF( nn_spp_dedt > 0 ) CALL spp_gen(kt, zsrp, nn_spp_dedt, rn_dedt_sd, jk_spp_deds ) 
    109119!CDIR COLLAPSE 
    110120               DO jj = 1, jpj 
    111121                  DO ji = 1, jpi 
    112                      zerp = zsrp * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths 
     122                     zerp = (zsrp(ji,jj)/rday) * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths 
    113123                        &        * ( sss_m(ji,jj) - sf_sss(1)%fnow(ji,jj,1) )  
    114124                     sfx(ji,jj) = sfx(ji,jj) + zerp                 ! salt flux 
    … …  
    117127               END DO 
    118128               CALL iom_put( "erp", erp )                             ! freshwater flux damping 
     129               CALL wrk_dealloc( jpi,jpj, zsrp ) 
    119130               ! 
    120131            ELSEIF( nn_sssr == 2 ) THEN                               !* Salinity damping term (volume flux (emp) and associated heat flux (qns) 
    121                zsrp = rn_deds / rday                                  ! from [mm/day] to [kg/m2/s] 
     132               CALL wrk_alloc( jpi, jpj, zsrp) 
     133               zsrp = rn_deds 
     134               IF( nn_spp_dedt > 0 ) CALL spp_gen(kt, zsrp, nn_spp_dedt, rn_dedt_sd, jk_spp_deds ) 
    122135               zerp_bnd = rn_sssr_bnd / rday                          !       -              -     
    123136!CDIR COLLAPSE 
    124137               DO jj = 1, jpj 
    125138                  DO ji = 1, jpi                             
    126                      zerp = zsrp * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths 
     139                     zerp = (zsrp(ji,jj)/rday) * ( 1. - 2.*rnfmsk(ji,jj) )   &      ! No damping in vicinity of river mouths 
    127140                        &        * ( sss_m(ji,jj) - sf_sss(1)%fnow(ji,jj,1) )   & 
    128141                        &        / MAX(  sss_m(ji,jj), 1.e-20   ) 
    … …  
    134147               END DO 
    135148               CALL iom_put( "erp", erp )                             ! freshwater flux damping 
     149               CALL wrk_dealloc( jpi,jpj,zsrp ) 
    136150            ENDIF 
    137151            ! 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/TRA/trabbc.F90

    r6486 r11394  
    2828   USE wrk_nemo        ! Memory Allocation 
    2929   USE timing          ! Timing 
     30   USE stopack 
    3031 
    3132   IMPLICIT NONE 
    … …  
    4142 
    4243   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   qgh_trd0   ! geothermal heating trend 
     44   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   qgh_trd1   ! geothermal heating trend 
    4345   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_qgh              ! structure of input qgh (file informations, fields read) 
    4446  
    … …  
    8991      ! 
    9092      !                             !  Add the geothermal heat flux trend on temperature 
     93 
     94      IF( nn_spp_geot .GT. 0) THEN 
     95          qgh_trd1 = qgh_trd0 
     96          CALL spp_gen(kt, qgh_trd0, nn_spp_geot, rn_geot_sd, jk_spp_geot) 
     97      ENDIF 
    9198      DO jj = 2, jpjm1 
    9299         DO ji = 2, jpim1 
    93100            ik = mbkt(ji,jj) 
    94             zqgh_trd = qgh_trd0(ji,jj) / fse3t(ji,jj,ik) 
     101            zqgh_trd = qgh_trd1(ji,jj) / fse3t(ji,jj,ik) 
    95102            tsa(ji,jj,ik,jp_tem) = tsa(ji,jj,ik,jp_tem) + zqgh_trd 
    96103         END DO 
    … …  
    163170         ! 
    164171         ALLOCATE( qgh_trd0(jpi,jpj) )    ! allocation 
     172         ALLOCATE( qgh_trd1(jpi,jpj) )    ! allocation 
    165173         ! 
    166174         SELECT CASE ( nn_geoflx )        ! geothermal heat flux / (rauO * Cp) 
    … …  
    192200            ! 
    193201         END SELECT 
     202         qgh_trd1 = qgh_trd0 
    194203         ! 
    195204      ELSE 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/TRA/trabbl.F90

    r7771 r11394  
    3939   USE timing         ! Timing 
    4040   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
     41   USE stopack 
    4142 
    4243   IMPLICIT NONE 
    … …  
    6768   INTEGER , ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC ::   mgrhu    , mgrhv       ! = +/-1, sign of grad(H) in u-(v-)direction (PUBLIC for TAM) 
    6869   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)         ::   ahu_bbl_0, ahv_bbl_0   ! diffusive bbl flux coefficients at u and v-points 
     70   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)         ::   ahu_bbl_1, ahv_bbl_1   ! diffusive bbl flux coefficients at u and v-points 
    6971   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC ::   e3u_bbl_0, e3v_bbl_0   ! thichness of the bbl (e3) at u and v-points (PUBLIC for TAM) 
    7072 
    … …  
    8688         &      vtr_bbl  (jpi,jpj) , ahv_bbl  (jpi,jpj) , mbkv_d  (jpi,jpj) , mgrhv(jpi,jpj) ,     & 
    8789         &      ahu_bbl_0(jpi,jpj) , ahv_bbl_0(jpi,jpj) ,                                          & 
     90         &      ahu_bbl_1(jpi,jpj) , ahv_bbl_1(jpi,jpj) ,                                          & 
    8891         &      e3u_bbl_0(jpi,jpj) , e3v_bbl_0(jpi,jpj) ,                                      STAT=tra_bbl_alloc ) 
    8992         ! 
    … …  
    195198      ALLOCATE(zptb(1:jpi, 1:jpj)) 
    196199      ! 
     200      ahu_bbl_1 = ahu_bbl 
     201      IF( nn_spp_ahubbl .GT. 0 ) THEN 
     202          CALL spp_gen(1 , ahu_bbl_1, nn_spp_ahubbl, rn_ahubbl_sd, jk_spp_ahubbl ) 
     203      ENDIF 
     204      ahv_bbl_1 = ahv_bbl 
     205      IF( nn_spp_ahvbbl .GT. 0 ) THEN 
     206          CALL spp_gen(1 , ahv_bbl_1, nn_spp_ahvbbl, rn_ahvbbl_sd, jk_spp_ahvbbl ) 
     207      ENDIF 
     208      ! 
    197209      DO jn = 1, kjpt                                     ! tracer loop 
    198210         !                                                ! =========== 
    … …  
    209221               zbtr = r1_e12t(ji,jj)  / fse3t(ji,jj,ik) 
    210222               pta(ji,jj,ik,jn) = pta(ji,jj,ik,jn)                                                         & 
    211                   &               + (   ahu_bbl(ji  ,jj  ) * ( zptb(ji+1,jj  ) - zptb(ji  ,jj  ) )   & 
    212                   &                   - ahu_bbl(ji-1,jj  ) * ( zptb(ji  ,jj  ) - zptb(ji-1,jj  ) )   & 
    213                   &                   + ahv_bbl(ji  ,jj  ) * ( zptb(ji  ,jj+1) - zptb(ji  ,jj  ) )   & 
    214                   &                   - ahv_bbl(ji  ,jj-1) * ( zptb(ji  ,jj  ) - zptb(ji  ,jj-1) )   ) * zbtr 
     223                  &               + (   ahu_bbl_1(ji  ,jj  ) * ( zptb(ji+1,jj  ) - zptb(ji  ,jj  ) )   & 
     224                  &                   - ahu_bbl_1(ji-1,jj  ) * ( zptb(ji  ,jj  ) - zptb(ji-1,jj  ) )   & 
     225                  &                   + ahv_bbl_1(ji  ,jj  ) * ( zptb(ji  ,jj+1) - zptb(ji  ,jj  ) )   & 
     226                  &                   - ahv_bbl_1(ji  ,jj-1) * ( zptb(ji  ,jj  ) - zptb(ji  ,jj-1) )   ) * zbtr 
    215227            END DO 
    216228         END DO 
    … …  
    549561      zmbk(:,:) = REAL( mbkv_d(:,:), wp )   ;   CALL lbc_lnk(zmbk,'V',1.)   ;   mbkv_d(:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
    550562 
    551                                         !* sign of grad(H) at u- and v-points 
    552       mgrhu(jpi,:) = 0   ;   mgrhu(:,jpj) = 0   ;   mgrhv(jpi,:) = 0   ;   mgrhv(:,jpj) = 0 
     563                                        !* sign of grad(H) at u- and v-points; zero if grad(H) = 0 
     564      mgrhu(:,:) = 0   ;   mgrhv(:,:) = 0 
    553565      DO jj = 1, jpjm1 
    554566         DO ji = 1, jpim1 
    555             mgrhu(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
    556             mgrhv(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     567            IF( gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) /= 0._wp ) THEN 
     568               mgrhu(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     569            ENDIF 
     570            ! 
     571            IF( gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) /= 0._wp ) THEN 
     572               mgrhv(ji,jj) = INT(  SIGN( 1.e0, gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) )  ) 
     573            ENDIF 
    557574         END DO 
    558575      END DO 
    … …  
    569586      ahu_bbl_0(:,:) = rn_ahtbbl * e2u(:,:) * e3u_bbl_0(:,:) / e1u(:,:)  * umask(:,:,1) 
    570587      ahv_bbl_0(:,:) = rn_ahtbbl * e1v(:,:) * e3v_bbl_0(:,:) / e2v(:,:)  * vmask(:,:,1) 
    571  
    572588 
    573589      IF( cp_cfg == "orca" ) THEN   !* ORCA configuration : regional enhancement of ah_bbl 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/TRA/traldf.F90

    r6498 r11394  
    3232   USE wrk_nemo        ! Memory allocation 
    3333   USE timing          ! Timing 
     34   USE stopack 
    3435 
    3536   IMPLICIT NONE 
    … …  
    4344   REAL, SAVE, ALLOCATABLE, DIMENSION(:,:,:) ::   t0_ldf, s0_ldf   !: lateral diffusion trends of T & S for a cst profile 
    4445   !                                                               !  (key_traldf_ano only) 
     46#if defined key_traldf_c3d 
     47   REAL, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ahtu0,ahtv0,ahtw0,ahtt0 
     48#endif 
     49#if defined key_traldf_c2d 
     50   REAL, SAVE, ALLOCATABLE, DIMENSION(:,:  ) :: ahtu0,ahtv0,ahtw0,ahtt0 
     51#endif 
    4552 
    4653   !! * Substitutions 
    … …  
    7582         ztrds(:,:,:) = tsa(:,:,:,jp_sal) 
    7683      ENDIF 
     84 
     85#if defined key_traldf_c3d 
     86         IF( kt .eq. nit000 .and. (nn_spp_ahtu+nn_spp_ahtv+nn_spp_ahtw+nn_spp_ahtt) .gt. 0 ) THEN 
     87             ALLOCATE ( ahtu0(jpi,jpj,jpk), ahtv0(jpi,jpj,jpk) ) 
     88             ALLOCATE ( ahtt0(jpi,jpj,jpk), ahtw0(jpi,jpj,jpk) ) 
     89             ahtu0 = ahtu 
     90             ahtv0 = ahtv 
     91             ahtw0 = ahtw 
     92             ahtt0 = ahtt 
     93         ENDIF 
     94#endif 
     95#if defined key_traldf_c2d 
     96         IF( kt .eq. nit000 .and. (nn_spp_ahtu+nn_spp_ahtv+nn_spp_ahtw+nn_spp_ahtt) .gt. 0 ) THEN 
     97             ALLOCATE ( ahtu0(jpi,jpj), ahtv0(jpi,jpj) ) 
     98             ALLOCATE ( ahtt0(jpi,jpj), ahtw0(jpi,jpj) ) 
     99             ahtu0 = ahtu 
     100             ahtv0 = ahtv 
     101             ahtw0 = ahtw 
     102             ahtt0 = ahtt 
     103         ENDIF 
     104#endif 
     105#if defined key_traldf_c3d || defined key_traldf_c2d 
     106         IF( nn_spp_ahtu .GT. 0) THEN 
     107             ahtu = ahtu0 
     108             CALL spp_aht(kt,ahtu,nn_spp_ahtu,rn_ahtu_sd,jk_spp_ahtu) 
     109         ENDIF 
     110         IF( nn_spp_ahtv .GT. 0) THEN 
     111             ahtv = ahtv0 
     112             CALL spp_aht(kt,ahtv,nn_spp_ahtv,rn_ahtv_sd,jk_spp_ahtv) 
     113         ENDIF 
     114         IF( nn_spp_ahtw .GT. 0) THEN 
     115             ahtw = ahtw0 
     116             CALL spp_aht(kt,ahtw,nn_spp_ahtw,rn_ahtw_sd,jk_spp_ahtw) 
     117         ENDIF 
     118         IF( nn_spp_ahtt .GT. 0) THEN 
     119             ahtt = ahtt0 
     120             CALL spp_aht(kt,ahtt,nn_spp_ahtt,rn_ahtt_sd,jk_spp_ahtt) 
     121         ENDIF 
     122#endif 
    77123 
    78124      SELECT CASE ( nldf )                       ! compute lateral mixing trend and add it to the general trend 
    … …  
    215261      IF( ierr == 1 )   CALL ctl_stop( ' iso-level in z-coordinate - partial step, not allowed' ) 
    216262      IF( ierr == 2 )   CALL ctl_stop( ' isoneutral bilaplacian operator does not exist' ) 
     263      IF( ln_traldf_hor .AND. ln_traldf_grif )    & 
     264            &   CALL ctl_stop( ' horizontal operator and Griffies triads not available; sitch to isoneutral operator' ) 
    217265      IF( ln_traldf_grif .AND. ln_isfcav         )   & 
    218266           CALL ctl_stop( ' ice shelf and traldf_grif not tested') 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/TRA/traqsr.F90

    r6498 r11394  
    3333   USE wrk_nemo       ! Memory Allocation 
    3434   USE timing         ! Timing 
     35   USE stopack 
    3536 
    3637   IMPLICIT NONE 
    … …  
    5253  
    5354   ! Module variables 
    54    REAL(wp) ::   xsi0r                           !: inverse of rn_si0 
     55   REAL(wp), ALLOCATABLE ::   xsi0r(:,:)         !: inverse of rn_si0 
    5556   REAL(wp) ::   xsi1r                           !: inverse of rn_si1 
    5657   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_chl   ! structure of input Chl (file informations, fields read) 
    … …  
    182183         !                                        ! ============================================== ! 
    183184         ! 
    184          !                                                ! ------------------------- ! 
     185         !  
     186         IF( nn_spp_qsi0 > 0 ) THEN 
     187             xsi0r = rn_si0 
     188             CALL spp_gen(kt, xsi0r, nn_spp_qsi0, rn_qsi0_sd, jk_spp_qsi0 ) 
     189             xsi0r = 1.e0 / xsi0r 
     190         ENDIF 
     191         !                                               ! ------------------------- ! 
    185192         IF( ln_qsr_rgb) THEN                             !  R-G-B  light penetration ! 
    186193            !                                             ! ------------------------- ! 
    … …  
    251258!CDIR NOVERRCHK    
    252259                     DO ji = 1, jpi 
    253                         zc0 = ze0(ji,jj,jk-1) * EXP( - fse3t(ji,jj,jk-1) * xsi0r    ) 
     260                        zc0 = ze0(ji,jj,jk-1) * EXP( - fse3t(ji,jj,jk-1) * xsi0r(ji,jj) ) 
    254261                        zc1 = ze1(ji,jj,jk-1) * EXP( - fse3t(ji,jj,jk-1) * zekb(ji,jj) ) 
    255262                        zc2 = ze2(ji,jj,jk-1) * EXP( - fse3t(ji,jj,jk-1) * zekg(ji,jj) ) 
    … …  
    263270                  END DO 
    264271               END DO 
     272               ! clem: store attenuation coefficient of the first ocean level 
     273               IF ( ln_qsr_ice ) THEN 
     274                  DO jj = 1, jpj 
     275                     DO ji = 1, jpi 
     276                        zzc0 = rn_abs * EXP( - fse3t(ji,jj,1) * xsi0r(ji,jj) ) 
     277                        zzc1 = zcoef  * EXP( - fse3t(ji,jj,1) * zekb(ji,jj) ) 
     278                        zzc2 = zcoef  * EXP( - fse3t(ji,jj,1) * zekg(ji,jj) ) 
     279                        zzc3 = zcoef  * EXP( - fse3t(ji,jj,1) * zekr(ji,jj) ) 
     280                        fraqsr_1lev(ji,jj) = 1.0 - ( zzc0 + zzc1 + zzc2  + zzc3  ) * tmask(ji,jj,2)  
     281                     END DO 
     282                  END DO 
     283               ENDIF 
    265284               ! 
    266285               DO jk = 1, nksr                                        ! compute and add qsr trend to ta 
    … …  
    310329            !                                             ! ------------------------- ! 
    311330            ! 
    312             IF( lk_vvl ) THEN                                  !* variable volume 
     331            IF( lk_vvl .OR. nn_spp_qsi0 > 0 ) THEN        !* variable volume 
     332 
    313333               zz0   =        rn_abs   * r1_rau0_rcp 
    314334               zz1   = ( 1. - rn_abs ) * r1_rau0_rcp 
    … …  
    316336                  DO jj = 1, jpj 
    317337                     DO ji = 1, jpi 
    318                         zc0 = zz0 * EXP( -fsdepw(ji,jj,jk  )*xsi0r ) + zz1 * EXP( -fsdepw(ji,jj,jk  )*xsi1r ) 
    319                         zc1 = zz0 * EXP( -fsdepw(ji,jj,jk+1)*xsi0r ) + zz1 * EXP( -fsdepw(ji,jj,jk+1)*xsi1r ) 
     338                        zc0 = zz0 * EXP( -fsdepw(ji,jj,jk  )*xsi0r(ji,jj) ) + zz1 * EXP( -fsdepw(ji,jj,jk  )*xsi1r ) 
     339                        zc1 = zz0 * EXP( -fsdepw(ji,jj,jk+1)*xsi0r(ji,jj) ) + zz1 * EXP( -fsdepw(ji,jj,jk+1)*xsi1r ) 
    320340                        qsr_hc(ji,jj,jk) = qsr(ji,jj) * ( zc0*tmask(ji,jj,jk) - zc1*tmask(ji,jj,jk+1) )  
    321341                     END DO 
    … …  
    326346                  DO jj = 1, jpj 
    327347                     DO ji = 1, jpi 
    328                         zc0 = zz0 * EXP( -fsdepw(ji,jj,1)*xsi0r ) + zz1 * EXP( -fsdepw(ji,jj,1)*xsi1r ) 
    329                         zc1 = zz0 * EXP( -fsdepw(ji,jj,2)*xsi0r ) + zz1 * EXP( -fsdepw(ji,jj,2)*xsi1r ) 
     348                        zc0 = zz0 * EXP( -fsdepw(ji,jj,1)*xsi0r(ji,jj) ) + zz1 * EXP( -fsdepw(ji,jj,1)*xsi1r ) 
     349                        zc1 = zz0 * EXP( -fsdepw(ji,jj,2)*xsi0r(ji,jj) ) + zz1 * EXP( -fsdepw(ji,jj,2)*xsi1r ) 
    330350                        fraqsr_1lev(ji,jj) = ( zc0*tmask(ji,jj,1) - zc1*tmask(ji,jj,2) ) / r1_rau0_rcp 
    331351                     END DO 
    … …  
    490510         !                       ! ===================================== ! 
    491511         ! 
     512         ALLOCATE( xsi0r(jpi,jpj) ) 
    492513         xsi0r = 1.e0 / rn_si0 
    493514         xsi1r = 1.e0 / rn_si1 
    … …  
    544565!CDIR NOVERRCHK    
    545566                        DO ji = 1, jpi 
    546                            zc0 = ze0(ji,jj,jk-1) * EXP( - e3t_0(ji,jj,jk-1) * xsi0r    ) 
     567                           zc0 = ze0(ji,jj,jk-1) * EXP( - e3t_0(ji,jj,jk-1) * xsi0r(ji,jj) ) 
    547568                           zc1 = ze1(ji,jj,jk-1) * EXP( - e3t_0(ji,jj,jk-1) * zekb(ji,jj) ) 
    548569                           zc2 = ze2(ji,jj,jk-1) * EXP( - e3t_0(ji,jj,jk-1) * zekg(ji,jj) ) 
    … …  
    585606                  DO jj = 1, jpj                              ! top 400 meters 
    586607                     DO ji = 1, jpi 
    587                         zc0 = zz0 * EXP( -fsdepw(ji,jj,jk  )*xsi0r ) + zz1 * EXP( -fsdepw(ji,jj,jk  )*xsi1r ) 
    588                         zc1 = zz0 * EXP( -fsdepw(ji,jj,jk+1)*xsi0r ) + zz1 * EXP( -fsdepw(ji,jj,jk+1)*xsi1r ) 
     608                        zc0 = zz0 * EXP( -fsdepw(ji,jj,jk  )*xsi0r(ji,jj) ) + zz1 * EXP( -fsdepw(ji,jj,jk  )*xsi1r ) 
     609                        zc1 = zz0 * EXP( -fsdepw(ji,jj,jk+1)*xsi0r(ji,jj) ) + zz1 * EXP( -fsdepw(ji,jj,jk+1)*xsi1r ) 
    589610                        etot3(ji,jj,jk) = (  zc0 * tmask(ji,jj,jk) - zc1 * tmask(ji,jj,jk+1)  ) * tmask(ji,jj,1)  
    590611                     END DO 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/TRD/trddyn.F90

    r6486 r11394  
    2929   USE lib_mpp        ! MPP library 
    3030   USE wrk_nemo       ! Memory allocation 
     31   USE stopack      
    3132 
    3233   IMPLICIT NONE 
    … …  
    6869      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    6970      IF( ln_dyn_trd )   CALL trd_dyn_iom( putrd, pvtrd, ktrd, kt ) 
     71      IF( ln_dyn_trd .AND. ln_sppt_dyn ) CALL dyn_sppt_collect( putrd, pvtrd, ktrd, kt ) 
    7072          
    7173      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/TRD/trdtra.F90

    r8104 r11394  
    2929   USE iom            ! I/O manager library 
    3030   USE lib_mpp        ! MPP library 
     31   USE stopack 
    3132   USE wrk_nemo       ! Memory allocation 
    3233 
    … …  
    249250      !                   ! 3D output of tracers trends using IOM interface 
    250251      IF( ln_tra_trd )   CALL trd_tra_iom ( ptrdx, ptrdy, ktrd, kt ) 
    251  
    252       !                   ! Integral Constraints Properties for tracers trends                                       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     252      IF( ln_tra_trd .AND. ln_sppt_tra )  CALL tra_sppt_collect( ptrdx, ptrdy, ktrd, kt ) 
     253 
     254      !  Integral Constraints Properties for tracers trends                                       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    253255      IF( ln_glo_trd )   CALL trd_glo( ptrdx, ptrdy, ktrd, 'TRA', kt ) 
    254256 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfbfr.F90

    r6486 r11394  
    2626   USE wrk_nemo        ! Memory Allocation 
    2727   USE phycst, ONLY: vkarmn 
     28   USE stopack 
    2829 
    2930   IMPLICIT NONE 
    … …  
    5253   REAL(wp), PUBLIC ::   rn_tfrz0     ! bottom roughness for loglayer bfr coeff (PUBLIC for TAM) 
    5354   LOGICAL , PUBLIC ::   ln_bfrimp    ! logical switch for implicit bottom friction 
    54    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC ::  bfrcoef2d, tfrcoef2d   ! 2D bottom/top drag coefficient (PUBLIC for TAM) 
     55   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC ::  bfrcoef2d, tfrcoef2d,bfrcoef2d0   ! 2D bottom/top drag coefficient (PUBLIC for TAM) 
    5556 
    5657   !! * Substitutions 
    … …  
    6869      !!                ***  FUNCTION zdf_bfr_alloc  *** 
    6970      !!---------------------------------------------------------------------- 
    70       ALLOCATE( bfrcoef2d(jpi,jpj), tfrcoef2d(jpi,jpj), STAT=zdf_bfr_alloc ) 
     71      ALLOCATE( bfrcoef2d(jpi,jpj), tfrcoef2d(jpi,jpj), bfrcoef2d0(jpi,jpj),STAT=zdf_bfr_alloc ) 
    7172      ! 
    7273      IF( lk_mpp             )   CALL mpp_sum ( zdf_bfr_alloc ) 
    … …  
    105106         WRITE(numout,*) 'zdf_bfr : Set bottom friction coefficient (non-linear case)' 
    106107         WRITE(numout,*) '~~~~~~~~' 
     108      ENDIF 
     109      !  
     110      IF( nn_spp_bfr > 0 ) THEN 
     111         bfrcoef2d = bfrcoef2d0 
     112         CALL spp_gen(kt, bfrcoef2d, nn_spp_bfr, rn_bfr_sd, jk_spp_bfr) 
    107113      ENDIF 
    108114      ! 
    … …  
    486492      ENDIF 
    487493      ! 
     494      bfrcoef2d0(:,:) = bfrcoef2d(:,:) 
    488495      IF( nn_timing == 1 )  CALL timing_stop('zdf_bfr_init') 
    489496      ! 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfevd.F90

    r7061 r11394  
    2525   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2626   USE timing          ! Timing 
     27   USE stopack          
    2728 
    2829   IMPLICIT NONE 
    … …  
    3031 
    3132   PUBLIC   zdf_evd    ! called by step.F90 
     33   REAL, SAVE, ALLOCATABLE, DIMENSION(:,:) :: rn_avevd0 
    3234 
    3335   !! * Substitutions 
    … …  
    6971         IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
    7072         IF(lwp) WRITE(numout,*) 
     73         ALLOCATE ( rn_avevd0(jpi,jpj) ) 
     74         rn_avevd0(:,:) = rn_avevd 
    7175      ENDIF 
    7276 
    7377      zavt_evd(:,:,:) = avt(:,:,:)           ! set avt prior to evd application 
     78 
     79      IF(nn_spp_aevd.GT.0) THEN 
     80         rn_avevd0(:,:) = rn_avevd 
     81         CALL spp_gen(kt, rn_avevd0, nn_spp_aevd, rn_aevd_sd, jk_spp_aevd) 
     82      ENDIF 
    7483 
    7584      SELECT CASE ( nn_evdm ) 
    … …  
    8897                  IF(  MIN( rn2(ji,jj,jk), rn2b(ji,jj,jk) ) <= -1.e-12 ) THEN 
    8998#endif 
    90                      avt (ji  ,jj  ,jk) = rn_avevd * tmask(ji  ,jj  ,jk) 
    91                      avm (ji  ,jj  ,jk) = rn_avevd * tmask(ji  ,jj  ,jk) 
    92                      avmu(ji  ,jj  ,jk) = rn_avevd * umask(ji  ,jj  ,jk) 
    93                      avmu(ji-1,jj  ,jk) = rn_avevd * umask(ji-1,jj  ,jk) 
    94                      avmv(ji  ,jj  ,jk) = rn_avevd * vmask(ji  ,jj  ,jk) 
    95                      avmv(ji  ,jj-1,jk) = rn_avevd * vmask(ji  ,jj-1,jk) 
     99                     avt (ji  ,jj  ,jk) = rn_avevd0(ji,jj) * tmask(ji  ,jj  ,jk) 
     100                     avm (ji  ,jj  ,jk) = rn_avevd0(ji,jj) * tmask(ji  ,jj  ,jk) 
     101                     avmu(ji  ,jj  ,jk) = rn_avevd0(ji,jj) * umask(ji  ,jj  ,jk) 
     102                     avmu(ji-1,jj  ,jk) = rn_avevd0(ji,jj) * umask(ji-1,jj  ,jk) 
     103                     avmv(ji  ,jj  ,jk) = rn_avevd0(ji,jj) * vmask(ji  ,jj  ,jk) 
     104                     avmv(ji  ,jj-1,jk) = rn_avevd0(ji,jj) * vmask(ji  ,jj-1,jk) 
    96105                  ENDIF 
    97106               END DO 
    … …  
    115124                  IF(  MIN( rn2(ji,jj,jk), rn2b(ji,jj,jk) ) <= -1.e-12 )   & 
    116125#endif 
    117                      avt(ji,jj,jk) = rn_avevd * tmask(ji,jj,jk) 
     126                     avt(ji,jj,jk) = rn_avevd0(ji,jj) * tmask(ji,jj,jk) 
    118127               END DO 
    119128            END DO 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90

    r6487 r11394  
    3232   USE timing         ! Timing 
    3333   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
     34   USE stopack 
    3435 
    3536   IMPLICIT NONE 
    … …  
    806807            END DO 
    807808         END DO 
     809         IF(nn_spp_avt > 0 ) CALL spp_gen(kt,avt(:,:,jk),nn_spp_avt,rn_avt_sd,jk) 
     810         IF(nn_spp_avm > 0 ) CALL spp_gen(kt,avm(:,:,jk),nn_spp_avm,rn_avm_sd,jk) 
    808811      END DO 
    809812      ! 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90

    r6498 r11394  
    5757   USE agrif_opa_update 
    5858#endif 
    59  
    60  
     59   USE stopack 
    6160 
    6261   IMPLICIT NONE 
    … …  
    9392   REAL(wp) ::   rhftau_scl = 1.0_wp       ! scale factor applied to HF part of taum  (nn_etau=3) 
    9493 
     94   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   en             !: now turbulent kinetic energy   [m2/s2] 
     95   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   rn_lc0, rn_ediff0, rn_ediss0, rn_ebb0, rn_efr0 
     96 
    9597   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   htau           ! depth of tke penetration (nn_htau) 
    9698   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   e_niw          !: TKE budget- near-inertial waves term 
    9799   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   efr            ! surface boundary condition for nn_etau = 4 
    98100   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   dissl          ! now mixing lenght of dissipation 
     101   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avt_k , avm_k  ! not enhanced Kz 
     102   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmu_k, avmv_k ! not enhanced Kz 
    99103#if defined key_c1d 
    100104   !                                                                        !!** 1D cfg only  **   ('key_c1d') 
    … …  
    184188         avmu(:,:,:) = avmu_k(:,:,:)  
    185189         avmv(:,:,:) = avmv_k(:,:,:)  
    186       ENDIF  
     190      ENDIF 
     191      ! 
     192      IF( nn_spp_tkelc > 0 ) THEN 
     193          rn_lc0 = rn_lc 
     194          CALL spp_gen(kt,rn_lc0,nn_spp_tkelc,rn_tkelc_sd,    jk_spp_tkelc ) 
     195      ENDIF 
     196      IF( nn_spp_tkedf > 0 ) THEN 
     197          rn_ediff0 = rn_ediff 
     198          CALL spp_gen(kt,rn_ediff0,nn_spp_tkedf,rn_tkedf_sd, jk_spp_tkedf ) 
     199      ENDIF 
     200      IF( nn_spp_tkeds > 0 ) THEN 
     201          rn_ediss0 = rn_ediss 
     202          CALL spp_gen(kt,rn_ediss0,nn_spp_tkeds,rn_tkeds_sd, jk_spp_tkeds ) 
     203      ENDIF 
     204      IF( nn_spp_tkebb > 0 ) THEN 
     205          rn_ebb0 = rn_ebb 
     206          CALL spp_gen(kt,rn_ebb0,nn_spp_tkebb,rn_tkebb_sd,   jk_spp_tkebb ) 
     207      ENDIF 
     208      IF( nn_spp_tkefr > 0 ) THEN 
     209          rn_efr0 = rn_efr 
     210          CALL spp_gen(kt,rn_efr0,nn_spp_tkefr,rn_tkefr_sd,   jk_spp_tkefr ) 
     211      ENDIF 
    187212      ! 
    188213      CALL tke_tke      ! now tke (en) 
    … …  
    199224      IF( .NOT.Agrif_Root() )   CALL Agrif_Update_Tke( kt )      ! children only 
    200225#endif       
    201      !  
     226      IF (  kt .eq. nitend ) THEN 
     227        DEALLOCATE ( rn_lc0, rn_ediff0, rn_ediss0, rn_ebb0, rn_efr0 ) 
     228      ENDIF 
     229      ! 
     230 
    202231   END SUBROUTINE zdf_tke 
    203232 
    … …  
    225254      REAL(wp) ::   zrhoa  = 1.22                   ! Air density kg/m3 
    226255      REAL(wp) ::   zcdrag = 1.5e-3                 ! drag coefficient 
    227       REAL(wp) ::   zbbrau, zesh2                   ! temporary scalars 
    228       REAL(wp) ::   zfact1, zfact2, zfact3          !    -         - 
     256      REAL(wp) ::   zesh2                           ! temporary scalars 
     257      REAL(wp) ::   zfact1                          !    -         - 
    229258      REAL(wp) ::   ztx2  , zty2  , zcof            !    -         - 
    230259      REAL(wp) ::   ztau  , zdif                    !    -         - 
    … …  
    233262!!bfr      REAL(wp) ::   zebot                           !    -         - 
    234263      INTEGER , POINTER, DIMENSION(:,:  ) :: imlc 
    235       REAL(wp), POINTER, DIMENSION(:,:  ) :: zhlc 
     264      REAL(wp), POINTER, DIMENSION(:,:  ) :: zhlc, zbbrau,zfact2,zfact3 
    236265      REAL(wp), POINTER, DIMENSION(:,:,:) :: zpelc, zdiag, zd_up, zd_lw 
    237266      !!-------------------------------------------------------------------- 
    … …  
    240269      ! 
    241270      CALL wrk_alloc( jpi,jpj, imlc )    ! integer 
    242       CALL wrk_alloc( jpi,jpj, zhlc )  
     271      CALL wrk_alloc( jpi,jpj, zhlc ) 
     272      CALL wrk_alloc( jpi,jpj, zbbrau ) 
     273      CALL wrk_alloc( jpi,jpj, zfact2 ) 
     274      CALL wrk_alloc( jpi,jpj, zfact3 ) 
    243275      CALL wrk_alloc( jpi,jpj,jpk, zpelc, zdiag, zd_up, zd_lw )  
    244276      ! 
    245       zbbrau = rn_ebb / rau0       ! Local constant initialisation 
     277      zbbrau = rn_ebb0 / rau0       ! Local constant initialisation 
    246278      zfact1 = -.5_wp * rdt  
    247       zfact2 = 1.5_wp * rdt * rn_ediss 
    248       zfact3 = 0.5_wp       * rn_ediss 
     279      zfact2 = 1.5_wp * rdt * rn_ediss0 
     280      zfact3 = 0.5_wp       * rn_ediss0 
    249281      ! 
    250282      ! 
    … …  
    261293      DO jj = 2, jpjm1            ! en(1)   = rn_ebb taum / rau0  (min value rn_emin0) 
    262294         DO ji = fs_2, fs_jpim1   ! vector opt. 
    263             en(ji,jj,1) = MAX( rn_emin0, zbbrau * taum(ji,jj) ) * tmask(ji,jj,1) 
     295            en(ji,jj,1) = MAX( rn_emin0, zbbrau(ji,jj) * taum(ji,jj) ) * tmask(ji,jj,1) 
    264296         END DO 
    265297      END DO 
    … …  
    326358                  !                                           ! vertical velocity due to LC 
    327359                  zind = 0.5 - SIGN( 0.5, fsdepw(ji,jj,jk) - zhlc(ji,jj) ) 
    328                   zwlc = zind * rn_lc * zus * SIN( rpi * fsdepw(ji,jj,jk) / zhlc(ji,jj) ) 
     360                  zwlc = zind * rn_lc0(ji,jj) * zus * SIN( rpi * fsdepw(ji,jj,jk) / zhlc(ji,jj) ) 
    329361                  !                                           ! TKE Langmuir circulation source term 
    330                   en(ji,jj,jk) = en(ji,jj,jk) + rdt * ( 1._wp - fr_i(ji,jj) ) * ( zwlc * zwlc * zwlc ) /   & 
     362                  en(ji,jj,jk) = en(ji,jj,jk) + rdt * MAX(0.,1._wp - 2.*fr_i(ji,jj) ) * ( zwlc * zwlc * zwlc ) /   & 
    331363                     &   zhlc(ji,jj) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
     364 
    332365               END DO 
    333366            END DO 
    … …  
    380413               zd_up(ji,jj,jk) = zzd_up            ! Matrix (zdiag, zd_up, zd_lw) 
    381414               zd_lw(ji,jj,jk) = zzd_lw 
    382                zdiag(ji,jj,jk) = 1._wp - zzd_lw - zzd_up + zfact2 * dissl(ji,jj,jk) * tmask(ji,jj,jk) 
     415               zdiag(ji,jj,jk) = 1._wp - zzd_lw - zzd_up + zfact2(ji,jj) * dissl(ji,jj,jk) * tmask(ji,jj,jk) 
    383416               ! 
    384417               !                                   ! right hand side in en 
    385418               en(ji,jj,jk) = en(ji,jj,jk) + rdt * (  zesh2  -   avt(ji,jj,jk) * rn2(ji,jj,jk)    & 
    386                   &                                 + zfact3 * dissl(ji,jj,jk) * en (ji,jj,jk)  ) & 
     419                  &                                 + zfact3(ji,jj) * dissl(ji,jj,jk) * en (ji,jj,jk)  ) & 
    387420                  &                                 * wmask(ji,jj,jk) 
    388421            END DO 
    … …  
    455488            DO jj = 2, jpjm1 
    456489               DO ji = fs_2, fs_jpim1   ! vector opt. 
    457                   en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) )   & 
    458                      &                                 * ( 1._wp - fr_i(ji,jj) )  * wmask(ji,jj,jk) * tmask(ji,jj,1) 
     490                  en(ji,jj,jk) = en(ji,jj,jk) + rn_efr0(ji,jj) * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) )   & 
     491                     &                                 * MAX(0.,1._wp - 2.*fr_i(ji,jj) )  * wmask(ji,jj,jk) * tmask(ji,jj,1) 
    459492               END DO 
    460493            END DO 
    … …  
    464497            DO ji = fs_2, fs_jpim1   ! vector opt. 
    465498               jk = nmln(ji,jj) 
    466                en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) )   & 
    467                   &                                 * ( 1._wp - fr_i(ji,jj) )  * wmask(ji,jj,jk) * tmask(ji,jj,1) 
     499               en(ji,jj,jk) = en(ji,jj,jk) + rn_efr0(ji,jj) * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) )   & 
     500                  &                                 * MAX(0.,1._wp - 2.*fr_i(ji,jj) )  * wmask(ji,jj,jk) * tmask(ji,jj,1) 
    468501            END DO 
    469502         END DO 
    … …  
    480513                  zdif = taum(ji,jj) - ztau                            ! mean of modulus - modulus of the mean  
    481514                  zdif = rhftau_scl * MAX( 0._wp, zdif + rhftau_add )  ! apply some modifications... 
    482                   en(ji,jj,jk) = en(ji,jj,jk) + zbbrau * zdif * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) )   & 
    483                      &                        * ( 1._wp - fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
     515                  en(ji,jj,jk) = en(ji,jj,jk) + zbbrau(ji,jj) * zdif * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) )   & 
     516                     &                        * MAX(0.,1._wp - 2.*fr_i(ji,jj) ) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
    484517               END DO 
    485518            END DO 
    … …  
    516549      CALL wrk_dealloc( jpi,jpj, imlc )    ! integer 
    517550      CALL wrk_dealloc( jpi,jpj, zhlc )  
     551      CALL wrk_dealloc( jpi,jpj, zbbrau )  
     552      CALL wrk_dealloc( jpi,jpj, zfact2 )  
     553      CALL wrk_dealloc( jpi,jpj, zfact3 )  
    518554      CALL wrk_dealloc( jpi,jpj,jpk, zpelc, zdiag, zd_up, zd_lw )  
    519555      ! 
    … …  
    698734            DO ji = fs_2, fs_jpim1   ! vector opt. 
    699735               zsqen = SQRT( en(ji,jj,jk) ) 
    700                zav   = rn_ediff * zmxlm(ji,jj,jk) * zsqen 
     736               zav   = rn_ediff0(ji,jj) * zmxlm(ji,jj,jk) * zsqen 
    701737               avm  (ji,jj,jk) = MAX( zav,                  avmb(jk) ) * wmask(ji,jj,jk) 
    702738               avt  (ji,jj,jk) = MAX( zav, avtb_2d(ji,jj) * avtb(jk) ) * wmask(ji,jj,jk) 
    … …  
    704740            END DO 
    705741         END DO 
     742         IF(nn_spp_avt > 0 ) CALL spp_gen(1 ,avt(:,:,jk),nn_spp_avt,rn_avt_sd, jk_spp_avt, jk) 
     743         IF(nn_spp_avm > 0 ) CALL spp_gen(1 ,avm(:,:,jk),nn_spp_avm,rn_avm_sd, jk_spp_avm, jk) 
    706744      END DO 
    707745      CALL lbc_lnk( avm, 'W', 1. )      ! Lateral boundary conditions (sign unchanged) 
    … …  
    840878         rn_mxl0 = rmxl_min 
    841879      ENDIF 
     880 
     881      ALLOCATE( rn_lc0   (jpi,jpj) ) ; rn_lc0    = rn_lc 
     882      ALLOCATE( rn_ediff0(jpi,jpj) ) ; rn_ediff0 = rn_ediff 
     883      ALLOCATE( rn_ediss0(jpi,jpj) ) ; rn_ediss0 = rn_ediss 
     884      ALLOCATE( rn_ebb0  (jpi,jpj) ) ; rn_ebb0   = rn_ebb 
     885      ALLOCATE( rn_efr0  (jpi,jpj) ) ; rn_efr0   = rn_efr 
    842886       
    843       IF( nn_etau == 2  ) THEN 
    844           ierr = zdf_mxl_alloc() 
    845           nmln(:,:) = nlb10           ! Initialization of nmln 
    846       ENDIF 
    847887 
    848888      IF( nn_etau /= 0 .and. nn_htau == 2 ) THEN 
    … …  
    950990              CALL tke_avn                               ! recompute avt, avm, avmu, avmv and dissl (approximation) 
    951991              ! 
    952               avt_k (:,:,:) = avt (:,:,:) 
    953               avm_k (:,:,:) = avm (:,:,:) 
    954               avmu_k(:,:,:) = avmu(:,:,:) 
    955               avmv_k(:,:,:) = avmv(:,:,:) 
    956               ! 
    957992              DO jit = nit000 + 1, nit000 + 10   ;   CALL zdf_tke( jit )   ;   END DO 
    958993           ENDIF 
    959994        ELSE                                   !* Start from rest 
    960995           en(:,:,:) = rn_emin * tmask(:,:,:) 
    961            DO jk = 1, jpk                           ! set the Kz to the background value 
    962               avt (:,:,jk) = avtb(jk) * wmask (:,:,jk) 
    963               avm (:,:,jk) = avmb(jk) * wmask (:,:,jk) 
    964               avmu(:,:,jk) = avmb(jk) * wumask(:,:,jk) 
    965               avmv(:,:,jk) = avmb(jk) * wvmask(:,:,jk) 
    966            END DO 
    967996        ENDIF 
    968997        ! 
     998        avt_k (:,:,:) = avt (:,:,:) 
     999        avm_k (:,:,:) = avm (:,:,:) 
     1000        avmu_k(:,:,:) = avmu(:,:,:) 
     1001        avmv_k(:,:,:) = avmv(:,:,:) 
     1002         
    9691003     ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN   ! Create restart file 
    9701004        !                                   ! ------------------- 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90

    r10301 r11394  
    8585   USE lbcnfd, ONLY: isendto, nsndto, nfsloop, nfeloop ! Setup of north fold exchanges  
    8686   USE sbc_oce, ONLY: lk_oasis 
     87   USE stopack 
    8788   USE stopar 
    8889   USE stopts 
    … …  
    483484                            CALL dia_hsb_init   ! heat content, salt content and volume budgets 
    484485                            CALL     trd_init   ! Mixed-layer/Vorticity/Integral constraints trends 
     486                            CALL stopack_init   ! STOPACK scheme 
    485487      IF( lk_diaobs     ) THEN                  ! Observation & model comparison 
    486488                            CALL dia_obs_init            ! Initialize observational data 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/step.F90

    r9288 r11394  
    118118         CALL lbc_lnk( tsb(:,:,:,tind), 'T', 1. ) 
    119119      END DO 
    120  
     120       
     121      IF( ln_stopack )   CALL stopack_pert( kstp ) 
    121122                         CALL sbc    ( kstp )         ! Sea Boundary Condition (including sea-ice) 
    122123                                                      ! clem: moved here for bdy ice purpose 
     124 
    123125      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    124126      ! Update stochastic parameters and random T/S fluctuations 
    125127      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    126        IF( ln_sto_eos ) CALL sto_par( kstp )          ! Stochastic parameters 
    127        IF( ln_sto_eos ) CALL sto_pts( tsn  )          ! Random T/S fluctuations 
     128 
     129      IF( ln_sto_eos )   CALL sto_par( kstp )          ! Stochastic parameters 
     130      IF( ln_sto_eos )   CALL sto_pts( tsn  )          ! Random T/S fluctuations 
     131      IF( ln_skeb  )     CALL skeb_comp( kstp ) 
    128132 
    129133      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    … …  
    149153      ENDIF 
    150154      IF( ln_rnf_mouth ) THEN                         ! increase diffusivity at rivers mouths 
    151          DO jk = 2, nkrnf   ;   avt(:,:,jk) = avt(:,:,jk) + 2.e0 * rn_avt_rnf * rnfmsk(:,:) * tmask(:,:,jk)   ;   END DO 
     155         IF ( nn_spp_arnf .GT. 0 ) THEN 
     156              rn_avt_rnf0 = rn_avt_rnf 
     157              CALL spp_gen( kstp, rn_avt_rnf0,nn_spp_arnf,rn_arnf_sd,jk_spp_arnf ) 
     158         ENDIF 
     159         DO jk = 2, nkrnf   ;   avt(:,:,jk) = avt(:,:,jk) + 2.e0 * rn_avt_rnf0(:,:) * rnfmsk(:,:) * tmask(:,:,jk)   ;   END DO 
    152160      ENDIF 
    153161      IF( ln_zdfevd  )   CALL zdf_evd( kstp )         ! enhanced vertical eddy diffusivity 
    … …  
    163171      IF( lrst_oce .AND. lk_zdftke )   CALL tke_rst( kstp, 'WRITE' ) 
    164172      IF( lrst_oce .AND. lk_zdfgls )   CALL gls_rst( kstp, 'WRITE' ) 
     173      IF( lrst_oce .AND. ln_stopack)   CALL stopack_rst( kstp, 'WRITE' ) 
    165174      ! 
    166175      !  LATERAL  PHYSICS 
    167176      ! 
    168177      IF( lk_ldfslp ) THEN                            ! slope of lateral mixing 
    169          CALL eos( tsb, rhd, gdept_0(:,:,:) )               ! before in situ density 
     178         IF(ln_sto_eos ) CALL sto_pts( tsn )          ! Random T/S fluctuations 
     179                         CALL eos( tsb, rhd, gdept_0(:,:,:) )               ! before in situ density 
    170180         IF( ln_zps .AND. .NOT. ln_isfcav)                               & 
    171181            &            CALL zps_hde    ( kstp, jpts, tsb, gtsu, gtsv,  &  ! Partial steps: before horizontal gradient 
    … …  
    202212          ! Note that the computation of vertical velocity above, hence "after" sea level 
    203213          ! is necessary to compute momentum advection for the rhs of barotropic loop: 
    204             CALL eos    ( tsn, rhd, rhop, fsdept_n(:,:,:) ) ! now in situ density for hpg computation 
     214            IF(ln_sto_eos ) CALL sto_pts( tsn )                             ! Random T/S fluctuations 
     215                            CALL eos    ( tsn, rhd, rhop, fsdept_n(:,:,:) ) ! now in situ density for hpg computation 
    205216            IF( ln_zps .AND. .NOT. ln_isfcav)                               & 
    206217               &            CALL zps_hde    ( kstp, jpts, tsn, gtsu, gtsv,  &    ! Partial steps: before horizontal gradient 
    … …  
    214225                                  va(:,:,:) = 0.e0 
    215226          IF(  lk_asminc .AND. ln_asmiau .AND. & 
    216              & ln_dyninc       )  CALL dyn_asm_inc  ( kstp )   ! apply dynamics assimilation increment 
    217           IF( ln_neptsimp )       CALL dyn_nept_cor ( kstp )   ! subtract Neptune velocities (simplified) 
    218           IF( lk_bdy           )  CALL bdy_dyn3d_dmp( kstp )   ! bdy damping trends 
    219                                   CALL dyn_adv      ( kstp )   ! advection (vector or flux form) 
    220                                   CALL dyn_vor      ( kstp )   ! vorticity term including Coriolis 
    221                                   CALL dyn_ldf      ( kstp )   ! lateral mixing 
    222           IF( ln_neptsimp )       CALL dyn_nept_cor ( kstp )   ! add Neptune velocities (simplified) 
     227             & ln_dyninc       )  CALL dyn_asm_inc   ( kstp )   ! apply dynamics assimilation increment 
     228          IF( ln_sppt_dyn )       CALL dyn_sppt_apply( kstp, 0 ) 
     229          IF( ln_neptsimp )       CALL dyn_nept_cor  ( kstp )   ! subtract Neptune velocities (simplified) 
     230          IF( lk_bdy           )  CALL bdy_dyn3d_dmp ( kstp )   ! bdy damping trends 
     231                                  CALL dyn_adv       ( kstp )   ! advection (vector or flux form) 
     232                                  CALL dyn_vor       ( kstp )   ! vorticity term including Coriolis 
     233                                  CALL dyn_ldf       ( kstp )   ! lateral mixing 
     234          IF( ln_neptsimp )       CALL dyn_nept_cor  ( kstp )   ! add Neptune velocities (simplified) 
    223235#if defined key_agrif 
    224236          IF(.NOT. Agrif_Root())  CALL Agrif_Sponge_dyn        ! momentum sponge 
    … …  
    262274                             tsa(:,:,:,:) = 0.e0            ! set tracer trends to zero 
    263275 
    264       IF(  lk_asminc .AND. ln_asmiau .AND. & 
    265          & ln_trainc     )   CALL tra_asm_inc( kstp )       ! apply tracer assimilation increment 
     276      IF( ln_sppt_tra )      CALL tra_sppt_apply ( kstp, 0 ) 
     277      IF( lk_asminc .AND. ln_asmiau .AND. & 
     278        & ln_trainc      )   CALL tra_asm_inc( kstp )       ! apply tracer assimilation increment 
    266279                             CALL tra_sbc    ( kstp )       ! surface boundary condition 
    267280      IF( ln_traqsr      )   CALL tra_qsr    ( kstp )       ! penetrative solar radiation qsr 
    … …  
    280293      IF(.NOT. Agrif_Root()) CALL Agrif_Sponge_tra          ! tracers sponge 
    281294#endif 
     295      IF( ln_sppt_tra )      CALL tra_sppt_apply ( kstp, 1 ) 
    282296                             CALL tra_zdf    ( kstp )       ! vertical mixing and after tracer fields 
    283297 
    … …  
    285299         IF( ln_zdfnpc   )   CALL tra_npc( kstp )                ! update after fields by non-penetrative convection 
    286300                             CALL tra_nxt( kstp )                ! tracer fields at next time step 
     301         IF( ln_sppt_tra )   CALL tra_sppt_apply ( kstp, 2 ) 
    287302                             CALL eos    ( tsa, rhd, rhop, fsdept_n(:,:,:) )  ! Time-filtered in situ density for hpg computation 
    288303            IF( ln_zps .AND. .NOT. ln_isfcav)                                & 
    … …  
    308323                             CALL tra_nxt( kstp )                ! tracer fields at next time step 
    309324         IF( ln_bias )       CALL dyn_bias( kstp ) 
     325         IF( ln_sppt_tra )   THEN 
     326                                CALL tra_sppt_apply ( kstp, 2 ) 
     327                                CALL eos    ( tsn, rhd, rhop, fsdept_n(:,:,:) )  ! now in situ density for hpg computation 
     328                                IF( ln_zps .AND. .NOT. ln_isfcav)                & 
     329                                & CALL zps_hde    ( kstp, jpts, tsn, gtsu, gtsv,  &    ! Partial steps: before horizontal gradient 
     330                                &                   rhd, gru , grv    )  ! of t, s, rd at the last ocean level 
     331                                IF( ln_zps .AND.       ln_isfcav)                                   &  
     332                                & CALL zps_hde_isf( kstp, jpts, tsn, gtsu, gtsv,  &    ! Partial steps for top cell (ISF) 
     333                                &                   rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv ,   & 
     334                                &                   gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) ! of t, s, rd at the last ocean level 
     335         ENDIF 
    310336      ENDIF 
    311337 
    … …  
    324350 
    325351                               CALL dyn_bfr( kstp )         ! bottom friction 
     352            IF( ln_sppt_dyn )  CALL dyn_sppt_apply ( kstp, 1 ) 
    326353                               CALL dyn_zdf( kstp )         ! vertical diffusion 
    327354      ELSE 
    … …  
    331358        IF(  lk_asminc .AND. ln_asmiau .AND. & 
    332359           & ln_dyninc      )  CALL dyn_asm_inc( kstp )     ! apply dynamics assimilation increment 
     360        IF( ln_sppt_dyn )      CALL dyn_sppt_apply ( kstp, 0 ) 
    333361        IF( ln_bkgwri )        CALL asm_bkg_wri( kstp )     ! output background fields 
    334362        IF( ln_neptsimp )      CALL dyn_nept_cor( kstp )    ! subtract Neptune velocities (simplified) 
    … …  
    343371                               CALL dyn_hpg( kstp )         ! horizontal gradient of Hydrostatic pressure 
    344372                               CALL dyn_bfr( kstp )         ! bottom friction 
     373            IF( ln_sppt_dyn )  CALL dyn_sppt_apply ( kstp, 1 ) 
    345374                               CALL dyn_zdf( kstp )         ! vertical diffusion 
    346375                               CALL dyn_spg( kstp, indic )  ! surface pressure gradient 
    347376      ENDIF 
    348377                               CALL dyn_nxt( kstp )         ! lateral velocity at next time step 
     378            IF( ln_sppt_dyn )  CALL dyn_sppt_apply ( kstp, 2 ) 
     379            IF( ln_skeb )      CALL skeb_apply ( kstp ) 
    349380 
    350381                               CALL ssh_swp( kstp )         ! swap of sea surface height 
    … …  
    381412         IF( lwm.AND.numoni /= -1 ) CALL FLUSH    ( numoni )     ! flush output namelist ice 
    382413      ENDIF 
     414      IF( ln_sto_eos       )   CALL sto_rst_write( kstp )   ! write restart file for stochastic parameters 
    383415 
    384416 

  • branches/UKMO/dev_r5518_GO6_package_FOAMv14_STOPACK/NEMOGCM/NEMO/OPA_SRC/step_oce.F90

    r8400 r11394  
    9696   USE diahsb           ! heat, salt and volume budgets    (dia_hsb routine) 
    9797   USE diaharm 
     98   USE diaprod          ! ocean model: product diagnostics 
    9899   USE flo_oce          ! floats variables 
    99100   USE floats           ! floats computation               (flo_stp routine) 
    … …  
    113114   USE timing           ! Timing 
    114115 
     116   USE stopack          ! Stochastic physics 
     117 
    115118#if defined key_agrif 
    116119   USE agrif_opa_sponge ! Momemtum and tracers sponges 
Note: See TracChangeset for help on using the changeset viewer.