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.
dynadv_cen2.F90 in branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/DYN – NEMO

source: branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/DYN/dynadv_cen2.F90 @ 11643

Last change on this file since 11643 was 11101, checked in by frrh, 5 years ago

Merge changes from Met Office GMED ticket 450 to reduce unnecessary
text output from NEMO.
This output, which is typically not switchable, is rarely of interest
in normal (non-debugging) runs and simply redunantley consumes extra
file space.
Further, the presence of this text output has been shown to
significantly degrade performance of models which are run during
Met Office HPC RAID (disk) checks.
The new code introduces switches which are configurable via the
changes made in the associated Met Office MOCI ticket 399.
File size: 8.4 KB
RevLine 
[643]1MODULE dynadv_cen2
2   !!======================================================================
3   !!                       ***  MODULE  dynadv  ***
4   !! Ocean dynamics: Update the momentum trend with the flux form advection
5   !!                 using a 2nd order centred scheme
6   !!======================================================================
[1566]7   !! History :  2.0  ! 2006-08  (G. Madec, S. Theetten)  Original code
8   !!            3.2  ! 2009-07  (R. Benshila)  Suppression of rigid-lid option
[643]9   !!----------------------------------------------------------------------
10
11   !!----------------------------------------------------------------------
12   !!   dyn_adv_cen2       : flux form momentum advection (ln_dynadv_cen2=T)
13   !!                        trends using a 2nd order centred scheme 
14   !!----------------------------------------------------------------------
15   USE oce            ! ocean dynamics and tracers
16   USE dom_oce        ! ocean space and time domain
[4990]17   USE trd_oce        ! trends: ocean variables
18   USE trddyn         ! trend manager: dynamics
19   !
[643]20   USE in_out_manager ! I/O manager
[2715]21   USE lib_mpp        ! MPP library
[1129]22   USE prtctl         ! Print control
[4990]23   USE wrk_nemo       ! Memory Allocation
24   USE timing         ! Timing
[643]25
26   IMPLICIT NONE
27   PRIVATE
28
[1566]29   PUBLIC   dyn_adv_cen2   ! routine called by step.F90
[643]30
31   !! * Substitutions
32#  include "domzgr_substitute.h90"
33#  include "vectopt_loop_substitute.h90"
34   !!----------------------------------------------------------------------
[2715]35   !! NEMO/OPA 4.0 , NEMO Consortium (2011)
[1152]36   !! $Id$
[2715]37   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
[643]38   !!----------------------------------------------------------------------
39CONTAINS
40
41   SUBROUTINE dyn_adv_cen2( kt )
42      !!----------------------------------------------------------------------
43      !!                  ***  ROUTINE dyn_adv_cen2  ***
44      !!
45      !! ** Purpose :   Compute the now momentum advection trend in flux form
[1566]46      !!              and the general trend of the momentum equation.
[643]47      !!
48      !! ** Method  :   Trend evaluated using now fields (centered in time)
49      !!
[1566]50      !! ** Action  :   (ua,va) updated with the now vorticity term trend
[643]51      !!----------------------------------------------------------------------
[1566]52      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
[2715]53      !
[1566]54      INTEGER  ::   ji, jj, jk   ! dummy loop indices
[2715]55      REAL(wp) ::   zbu, zbv     ! local scalars
[3294]56      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zfu_t, zfv_t, zfu_f, zfv_f, zfu_uw, zfv_vw, zfw
57      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zfu, zfv
[643]58      !!----------------------------------------------------------------------
[3294]59      !
60      IF( nn_timing == 1 )  CALL timing_start('dyn_adv_cen2')
61      !
62      CALL wrk_alloc( jpi, jpj, jpk, zfu_t, zfv_t, zfu_f, zfv_f, zfu_uw, zfv_vw, zfu, zfv, zfw )
63      !
[2715]64      IF( kt == nit000 .AND. lwp ) THEN
65         WRITE(numout,*)
66         WRITE(numout,*) 'dyn_adv_cen2 : 2nd order flux form momentum advection'
67         WRITE(numout,*) '~~~~~~~~~~~~'
[11101]68         IF(lflush) CALL flush(numout)
[643]69      ENDIF
[3294]70      !
[1129]71      IF( l_trddyn ) THEN           ! Save ua and va trends
72         zfu_uw(:,:,:) = ua(:,:,:)
73         zfv_vw(:,:,:) = va(:,:,:)
74      ENDIF
[643]75
[1566]76      !                                      ! ====================== !
77      !                                      !  Horizontal advection  !
78      DO jk = 1, jpkm1                       ! ====================== !
79         !                                         ! horizontal volume fluxes
[643]80         zfu(:,:,jk) = 0.25 * e2u(:,:) * fse3u(:,:,jk) * un(:,:,jk)
81         zfv(:,:,jk) = 0.25 * e1v(:,:) * fse3v(:,:,jk) * vn(:,:,jk)
[1566]82         !
83         DO jj = 1, jpjm1                          ! horizontal momentum fluxes at T- and F-point
[643]84            DO ji = 1, fs_jpim1   ! vector opt.
85               zfu_t(ji+1,jj  ,jk) = ( zfu(ji,jj,jk) + zfu(ji+1,jj  ,jk) ) * ( un(ji,jj,jk) + un(ji+1,jj  ,jk) )
86               zfv_f(ji  ,jj  ,jk) = ( zfv(ji,jj,jk) + zfv(ji+1,jj  ,jk) ) * ( un(ji,jj,jk) + un(ji  ,jj+1,jk) )
87               zfu_f(ji  ,jj  ,jk) = ( zfu(ji,jj,jk) + zfu(ji  ,jj+1,jk) ) * ( vn(ji,jj,jk) + vn(ji+1,jj  ,jk) )
88               zfv_t(ji  ,jj+1,jk) = ( zfv(ji,jj,jk) + zfv(ji  ,jj+1,jk) ) * ( vn(ji,jj,jk) + vn(ji  ,jj+1,jk) )
89            END DO
90         END DO
[1566]91         DO jj = 2, jpjm1                          ! divergence of horizontal momentum fluxes
[643]92            DO ji = fs_2, fs_jpim1   ! vector opt.
93               zbu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk)
94               zbv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk)
[1566]95               !
96               ua(ji,jj,jk) = ua(ji,jj,jk) - (  zfu_t(ji+1,jj  ,jk) - zfu_t(ji  ,jj  ,jk)    &
97                  &                           + zfv_f(ji  ,jj  ,jk) - zfv_f(ji  ,jj-1,jk)  ) / zbu
98               va(ji,jj,jk) = va(ji,jj,jk) - (  zfu_f(ji  ,jj  ,jk) - zfu_f(ji-1,jj  ,jk)    &
99                  &                           + zfv_t(ji  ,jj+1,jk) - zfv_t(ji  ,jj  ,jk)  ) / zbv
[643]100            END DO
101         END DO
[1566]102      END DO
103      !
104      IF( l_trddyn ) THEN                          ! save the horizontal advection trend for diagnostic
[1129]105         zfu_uw(:,:,:) = ua(:,:,:) - zfu_uw(:,:,:)
106         zfv_vw(:,:,:) = va(:,:,:) - zfv_vw(:,:,:)
[4990]107         CALL trd_dyn( zfu_uw, zfv_vw, jpdyn_keg, kt )
[1129]108         zfu_t(:,:,:) = ua(:,:,:)
109         zfv_t(:,:,:) = va(:,:,:)
110      ENDIF
[1566]111      !
[1129]112
[1566]113      !                                      ! ==================== !
114      !                                      !  Vertical advection  !
115      DO jk = 1, jpkm1                       ! ==================== !
116         !                                         ! Vertical volume fluxesÊ
[643]117         zfw(:,:,jk) = 0.25 * e1t(:,:) * e2t(:,:) * wn(:,:,jk)
[1566]118         !
119         IF( jk == 1 ) THEN                        ! surface/bottom advective fluxes                   
120            zfu_uw(:,:,jpk) = 0.e0                      ! Bottom  value : flux set to zero
[643]121            zfv_vw(:,:,jpk) = 0.e0
[1566]122            !                                           ! Surface value :
123            IF( lk_vvl ) THEN                                ! variable volume : flux set to zero
[643]124               zfu_uw(:,:, 1 ) = 0.e0   
125               zfv_vw(:,:, 1 ) = 0.e0
[1566]126            ELSE                                             ! constant volume : advection through the surface
[643]127               DO jj = 2, jpjm1
128                  DO ji = fs_2, fs_jpim1
129                     zfu_uw(ji,jj, 1 ) = 2.e0 * ( zfw(ji,jj,1) + zfw(ji+1,jj  ,1) ) * un(ji,jj,1)
130                     zfv_vw(ji,jj, 1 ) = 2.e0 * ( zfw(ji,jj,1) + zfw(ji  ,jj+1,1) ) * vn(ji,jj,1)
131                  END DO
132               END DO
133            ENDIF
[1566]134         ELSE                                      ! interior fluxes
[643]135            DO jj = 2, jpjm1
136               DO ji = fs_2, fs_jpim1   ! vector opt.
137                  zfu_uw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji+1,jj  ,jk) ) * ( un(ji,jj,jk) + un(ji,jj,jk-1) )
138                  zfv_vw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji  ,jj+1,jk) ) * ( vn(ji,jj,jk) + vn(ji,jj,jk-1) )
139               END DO
140            END DO
141         ENDIF
142      END DO
[1566]143      DO jk = 1, jpkm1                             ! divergence of vertical momentum flux divergence
[643]144         DO jj = 2, jpjm1 
145            DO ji = fs_2, fs_jpim1   ! vector opt.
[1566]146               ua(ji,jj,jk) =  ua(ji,jj,jk) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) )    &
[643]147                  &  / ( e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) )
[1566]148               va(ji,jj,jk) =  va(ji,jj,jk) - ( zfv_vw(ji,jj,jk) - zfv_vw(ji,jj,jk+1) )    &
[643]149                  &  / ( e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) )
150            END DO
151         END DO
152      END DO
[1566]153      !
154      IF( l_trddyn ) THEN                          ! save the vertical advection trend for diagnostic
[1129]155         zfu_t(:,:,:) = ua(:,:,:) - zfu_t(:,:,:)
156         zfv_t(:,:,:) = va(:,:,:) - zfv_t(:,:,:)
[4990]157         CALL trd_dyn( zfu_t, zfv_t, jpdyn_zad, kt )
[1129]158      ENDIF
[1566]159      !                                            ! Control print
[1129]160      IF(ln_ctl)   CALL prt_ctl( tab3d_1=ua, clinfo1=' cen2 adv - Ua: ', mask1=umask,   &
161         &                       tab3d_2=va, clinfo2=           ' Va: ', mask2=vmask, clinfo3='dyn' )
162      !
[3294]163      CALL wrk_dealloc( jpi, jpj, jpk, zfu_t, zfv_t, zfu_f, zfv_f, zfu_uw, zfv_vw, zfu, zfv, zfw )
[2715]164      !
[3294]165      IF( nn_timing == 1 )  CALL timing_stop('dyn_adv_cen2')
166      !
[643]167   END SUBROUTINE dyn_adv_cen2
168
169   !!==============================================================================
170END MODULE dynadv_cen2
Note: See TracBrowser for help on using the repository browser.