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 6808 for branches/NERC/dev_r5549_BDY_ZEROGRAD/NEMOGCM/NEMO/OPA_SRC/TRD/trdken.F90 – NEMO

Ignore:
Timestamp:
2016-07-19T10:38:35+02:00 (8 years ago)
Author:
jamesharle
Message:

merge with trunk@6232 for consistency with SSB code

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/NERC/dev_r5549_BDY_ZEROGRAD/NEMOGCM/NEMO/OPA_SRC/TRD/trdken.F90

    r5215 r6808  
    1313   USE oce            ! ocean dynamics and tracers variables 
    1414   USE dom_oce        ! ocean space and time domain variables 
     15   USE sbc_oce        ! surface boundary condition: ocean 
    1516   USE zdf_oce        ! ocean vertical physics variables 
    1617   USE trd_oce        ! trends: ocean variables 
    1718!!gm   USE dynhpg          ! hydrostatic pressure gradient    
    1819   USE zdfbfr         ! bottom friction 
    19    USE ldftra_oce     ! ocean active tracers lateral physics 
    20    USE sbc_oce        ! surface boundary condition: ocean 
     20   USE ldftra         ! ocean active tracers lateral physics 
    2121   USE phycst         ! physical constants 
    2222   USE trdvor         ! ocean vorticity trends  
    2323   USE trdglo         ! trends:global domain averaged 
    24    USE trdmxl         ! ocean active mixed layer tracers trends  
     24   USE trdmxl         ! ocean active mixed layer tracers trends 
     25   ! 
    2526   USE in_out_manager ! I/O manager 
    2627   USE iom            ! I/O manager library 
     
    4041 
    4142   !! * Substitutions 
    42 #  include "domzgr_substitute.h90" 
    4343#  include "vectopt_loop_substitute.h90" 
    4444   !!---------------------------------------------------------------------- 
     
    9393      CALL lbc_lnk( putrd, 'U', -1. )   ;   CALL lbc_lnk( pvtrd, 'V', -1. )      ! lateral boundary conditions 
    9494      ! 
    95       IF ( lk_vvl .AND. kt /= nkstp ) THEN   ! Variable volume: set box volume at the 1st call of kt time step 
    96          nkstp = kt 
    97          DO jk = 1, jpkm1 
    98             bu   (:,:,jk) =  e1u(:,:) * e2u(:,:) * fse3u_n(:,:,jk) 
    99             bv   (:,:,jk) =  e1v(:,:) * e2v(:,:) * fse3v_n(:,:,jk) 
    100             r1_bt(:,:,jk) = 1._wp / ( e1e2t(:,:) * fse3t_n(:,:,jk) ) * tmask(:,:,jk) 
    101          END DO 
    102       ENDIF 
     95      nkstp = kt 
     96      DO jk = 1, jpkm1 
     97         bu   (:,:,jk) =    e1e2u(:,:) * e3u_n(:,:,jk) 
     98         bv   (:,:,jk) =    e1e2v(:,:) * e3v_n(:,:,jk) 
     99         r1_bt(:,:,jk) = r1_e1e2t(:,:) / e3t_n(:,:,jk) * tmask(:,:,jk) 
     100      END DO 
    103101      ! 
    104102      zke(:,:,jpk) = 0._wp 
     
    117115      ! 
    118116      SELECT CASE( ktrd ) 
    119          CASE( jpdyn_hpg )   ;   CALL iom_put( "ketrd_hpg", zke )    ! hydrostatic pressure gradient 
    120          CASE( jpdyn_spg )   ;   CALL iom_put( "ketrd_spg", zke )    ! surface pressure gradient 
    121          CASE( jpdyn_spgexp );   CALL iom_put( "ketrd_spgexp", zke ) ! surface pressure gradient (explicit) 
    122          CASE( jpdyn_spgflt );   CALL iom_put( "ketrd_spgflt", zke ) ! surface pressure gradient (filter) 
    123          CASE( jpdyn_pvo )   ;   CALL iom_put( "ketrd_pvo", zke )    ! planetary vorticity 
    124          CASE( jpdyn_rvo )   ;   CALL iom_put( "ketrd_rvo", zke )    ! relative  vorticity     (or metric term) 
    125          CASE( jpdyn_keg )   ;   CALL iom_put( "ketrd_keg", zke )    ! Kinetic Energy gradient (or had) 
    126          CASE( jpdyn_zad )   ;   CALL iom_put( "ketrd_zad", zke )    ! vertical   advection 
    127          CASE( jpdyn_ldf )   ;   CALL iom_put( "ketrd_ldf", zke )    ! lateral diffusion 
    128          CASE( jpdyn_zdf )   ;   CALL iom_put( "ketrd_zdf", zke )    ! vertical diffusion  
    129                                  !                                   ! wind stress trends 
    130                                  CALL wrk_alloc( jpi, jpj, z2dx, z2dy, zke2d ) 
    131                            z2dx(:,:) = un(:,:,1) * ( utau_b(:,:) + utau(:,:) ) * e1u(:,:) * e2u(:,:) * umask(:,:,1) 
    132                            z2dy(:,:) = vn(:,:,1) * ( vtau_b(:,:) + vtau(:,:) ) * e1v(:,:) * e2v(:,:) * vmask(:,:,1) 
    133                            zke2d(1,:) = 0._wp   ;   zke2d(:,1) = 0._wp 
    134                            DO jj = 2, jpj 
    135                               DO ji = 2, jpi 
    136                                  zke2d(ji,jj) = 0.5_wp * (   z2dx(ji,jj) + z2dx(ji-1,jj)   & 
    137                                  &                         + z2dy(ji,jj) + z2dy(ji,jj-1)   ) * r1_bt(ji,jj,1) 
    138                               END DO 
    139                            END DO 
    140                                  CALL iom_put( "ketrd_tau", zke2d ) 
    141                                  CALL wrk_dealloc( jpi, jpj     , z2dx, z2dy, zke2d ) 
    142          CASE( jpdyn_bfr )   ;   CALL iom_put( "ketrd_bfr", zke )    ! bottom friction (explicit case)  
     117         CASE( jpdyn_hpg )   ;   CALL iom_put( "ketrd_hpg"   , zke )    ! hydrostatic pressure gradient 
     118         CASE( jpdyn_spg )   ;   CALL iom_put( "ketrd_spg"   , zke )    ! surface pressure gradient 
     119         CASE( jpdyn_pvo )   ;   CALL iom_put( "ketrd_pvo"   , zke )    ! planetary vorticity 
     120         CASE( jpdyn_rvo )   ;   CALL iom_put( "ketrd_rvo"   , zke )    ! relative  vorticity     (or metric term) 
     121         CASE( jpdyn_keg )   ;   CALL iom_put( "ketrd_keg"   , zke )    ! Kinetic Energy gradient (or had) 
     122         CASE( jpdyn_zad )   ;   CALL iom_put( "ketrd_zad"   , zke )    ! vertical   advection 
     123         CASE( jpdyn_ldf )   ;   CALL iom_put( "ketrd_ldf"   , zke )    ! lateral diffusion 
     124         CASE( jpdyn_zdf )   ;   CALL iom_put( "ketrd_zdf"   , zke )    ! vertical diffusion  
     125         !                   !                                          ! wind stress trends 
     126                                 CALL wrk_alloc( jpi, jpj, z2dx, z2dy, zke2d ) 
     127                           z2dx(:,:) = un(:,:,1) * ( utau_b(:,:) + utau(:,:) ) * e1e2u(:,:) * umask(:,:,1) 
     128                           z2dy(:,:) = vn(:,:,1) * ( vtau_b(:,:) + vtau(:,:) ) * e1e2v(:,:) * vmask(:,:,1) 
     129                           zke2d(1,:) = 0._wp   ;   zke2d(:,1) = 0._wp 
     130                           DO jj = 2, jpj 
     131                              DO ji = 2, jpi 
     132                                 zke2d(ji,jj) = r1_rau0 * 0.5_wp * (   z2dx(ji,jj) + z2dx(ji-1,jj)   & 
     133                                 &                                   + z2dy(ji,jj) + z2dy(ji,jj-1)   ) * r1_bt(ji,jj,1) 
     134                              END DO 
     135                           END DO 
     136                                 CALL iom_put( "ketrd_tau"   , zke2d )  !  
     137                                 CALL wrk_dealloc( jpi, jpj     , z2dx, z2dy, zke2d ) 
     138         CASE( jpdyn_bfr )   ;   CALL iom_put( "ketrd_bfr"   , zke )    ! bottom friction (explicit case)  
    143139!!gm TO BE DONE properly 
    144140!!gm only valid if ln_bfrimp=F otherwise the bottom stress as to be recomputed at the end of the computation.... 
     
    159155!               END DO 
    160156!            END DO 
    161 !                              CALL iom_put( "ketrd_bfr", zke2d )    ! bottom friction (explicit case) 
     157!                                    CALL iom_put( "ketrd_bfr"  , zke2d )   ! bottom friction (explicit case) 
    162158!         ENDIF 
    163159!!gm end 
    164          CASE( jpdyn_atf )   ;   CALL iom_put( "ketrd_atf", zke )    ! asselin filter trends  
     160         CASE( jpdyn_atf )   ;   CALL iom_put( "ketrd_atf"   , zke )    ! asselin filter trends  
    165161!! a faire !!!!  idee changer dynnxt pour avoir un appel a jpdyn_bfr avant le swap !!! 
    166162!! reflechir a une possible sauvegarde du "vrai" un,vn pour le calcul de atf.... 
     
    171167!                  ikbu = mbku(ji,jj)          ! deepest ocean u- & v-levels 
    172168!                  ikbv = mbkv(ji,jj) 
    173 !                  z2dx(ji,jj) = un(ji,jj,ikbu) * bfrua(ji,jj) * un(ji,jj,ikbu) / fse3u(ji,jj,ikbu) 
    174 !                  z2dy(ji,jj) = un(ji,jj,ikbu) * bfrva(ji,jj) * vn(ji,jj,ikbv) / fse3v(ji,jj,ikbv) 
     169!                  z2dx(ji,jj) = un(ji,jj,ikbu) * bfrua(ji,jj) * un(ji,jj,ikbu) / e3u_n(ji,jj,ikbu) 
     170!                  z2dy(ji,jj) = un(ji,jj,ikbu) * bfrva(ji,jj) * vn(ji,jj,ikbv) / e3v_n(ji,jj,ikbv) 
    175171!               END DO 
    176172!            END DO 
     
    184180!                              CALL iom_put( "ketrd_bfri", zke2d ) 
    185181!         ENDIF 
    186          CASE( jpdyn_ken )   ;   ! kinetic energy 
    187                            ! called in dynnxt.F90 before asselin time filter 
    188                            ! with putrd=ua and pvtrd=va 
    189                            zke(:,:,:) = 0.5_wp * zke(:,:,:) 
    190                            CALL iom_put( "KE", zke ) 
    191                            ! 
    192                            CALL ken_p2k( kt , zke ) 
    193                            CALL iom_put( "ketrd_convP2K", zke )     ! conversion -rau*g*w 
     182         CASE( jpdyn_ken )   ;                                          ! kinetic energy 
     183                                 ! called in dynnxt.F90 before asselin time filter with putrd=ua and pvtrd=va 
     184                                 zke(:,:,:) = 0.5_wp * zke(:,:,:) 
     185                                 CALL iom_put( "KE", zke ) 
     186                                 ! 
     187                                 CALL ken_p2k( kt , zke ) 
     188                                 CALL iom_put( "ketrd_convP2K", zke )     ! conversion -rau*g*w 
    194189         ! 
    195190      END SELECT 
     
    226221       
    227222      !  Surface value (also valid in partial step case) 
    228       zconv(:,:,1) = zcoef * ( 2._wp * rhd(:,:,1) ) * wn(:,:,1) * fse3w(:,:,1) 
     223      zconv(:,:,1) = zcoef * ( 2._wp * rhd(:,:,1) ) * wn(:,:,1) * e3w_n(:,:,1) 
    229224 
    230225      ! interior value (2=<jk=<jpkm1) 
    231226      DO jk = 2, jpk 
    232          zconv(:,:,jk) = zcoef * ( rhd(:,:,jk) + rhd(:,:,jk-1) ) * wn(:,:,jk) * fse3w(:,:,jk) 
     227         zconv(:,:,jk) = zcoef * ( rhd(:,:,jk) + rhd(:,:,jk-1) ) * wn(:,:,jk) * e3w_n(:,:,jk) 
    233228      END DO 
    234229 
     
    237232         DO jj = 1, jpj 
    238233            DO ji = 1, jpi 
    239                zcoef = 0.5_wp / fse3t(ji,jj,jk) 
     234               zcoef = 0.5_wp / e3t_n(ji,jj,jk) 
    240235               pconv(ji,jj,jk) = zcoef * ( zconv(ji,jj,jk) + zconv(ji,jj,jk+1) ) * tmask(ji,jj,jk) 
    241236            END DO 
     
    263258      ENDIF 
    264259      !                           ! allocate box volume arrays 
    265       IF ( trd_ken_alloc() /= 0 )   CALL ctl_stop('trd_ken_alloc: failed to allocate arrays') 
    266       ! 
    267 !!gm      IF( .NOT. (ln_hpg_zco.OR.ln_hpg_zps) )   & 
    268 !!gm         &   CALL ctl_stop('trd_ken_init : only full and partial cells are coded for conversion rate') 
    269       ! 
    270       IF ( .NOT.lk_vvl ) THEN     ! constant volume: bu, bv, 1/bt computed one for all 
    271          DO jk = 1, jpkm1 
    272             bu   (:,:,jk) =  e1u(:,:) * e2u(:,:) * fse3u_n(:,:,jk) 
    273             bv   (:,:,jk) =  e1v(:,:) * e2v(:,:) * fse3v_n(:,:,jk) 
    274             r1_bt(:,:,jk) = 1._wp / ( e1e2t(:,:) * fse3t_n(:,:,jk) ) 
    275          END DO 
    276       ENDIF 
     260      IF( trd_ken_alloc() /= 0 )   CALL ctl_stop('trd_ken_alloc: failed to allocate arrays') 
    277261      ! 
    278262   END SUBROUTINE trd_ken_init 
Note: See TracChangeset for help on using the changeset viewer.