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.
dynkeg.F90 in trunk/NEMO/OPA_SRC/DYN – NEMO

source: trunk/NEMO/OPA_SRC/DYN/dynkeg.F90 @ 3

Last change on this file since 3 was 3, checked in by opalod, 20 years ago

Initial revision

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 5.5 KB
Line 
1MODULE dynkeg
2   !!======================================================================
3   !!                       ***  MODULE  dynkeg  ***
4   !! Ocean dynamics:  kinetic energy gradient trend
5   !!======================================================================
6   
7   !!----------------------------------------------------------------------
8   !!   dyn_keg      : update the momentum trend with the horizontal tke
9   !!----------------------------------------------------------------------
10   !! * Modules used
11   USE oce             ! ocean dynamics and tracers
12   USE dom_oce         ! ocean space and time domain
13   USE in_out_manager  ! I/O manager
14   USE trddyn_oce      ! dynamics trends diagnostics
15
16   IMPLICIT NONE
17   PRIVATE
18
19   !! * Accessibility
20   PUBLIC dyn_keg                ! routine called by step.F90
21   
22   !! * Substitutions
23#  include "vectopt_loop_substitute.h90"
24   !!---------------------------------------------------------------------------------
25   !!   OPA 9.0 , LODYC-IPSL  (2003)
26   !!---------------------------------------------------------------------------------
27
28CONTAINS
29
30   SUBROUTINE dyn_keg( kt )
31      !!----------------------------------------------------------------------
32      !!                  ***  ROUTINE dyn_keg  ***
33      !!
34      !! ** Purpose :   Compute the now momentum trend due to the horizontal
35      !!      gradient of the horizontal kinetic energy and add it to the
36      !!      general momentum trend.
37      !!
38      !! ** Method  :   Compute the now horizontal kinetic energy:
39      !!         zhke = 1/2 [ mi-1( un^2 ) + mj-1( vn^2 ) ]
40      !!      Take its horizontal gradient and add it to the general momentum
41      !!      trend (ua,va).
42      !!         ua = ua - 1/e1u di[ zhke ]
43      !!         va = va - 1/e2v dj[ zhke ]
44      !!
45      !! ** Action : - Update the (ua, va) with the hor. ke gradient trend
46      !!             - Save the trends in (utrd,vtrd) ('key_trddyn')
47      !!
48      !! History :
49      !!   1.0  !  87-09  (P. Andrich, m.-a. Foujols)  Original code
50      !!   7.0  !  97-05  (G. Madec)  Split dynber into dynkeg and dynhpg
51      !!   9.0  !  02-07  (G. Madec)  F90: Free form and module
52      !!----------------------------------------------------------------------
53      !! * modules used     
54      USE oce, ONLY :   zhke => ta      ! use ta as 3D workspace   
55
56      !! * Arguments
57      INTEGER, INTENT( in ) ::   kt     ! ocean time-step index
58
59      !! * Local declarations
60      INTEGER  ::   ji, jj, jk          ! dummy loop indices
61      REAL(wp) ::   zua, zva, zu, zv    ! temporary scalars
62#if defined key_trddyn_new
63      REAL(wp) ::   zuu, zvv            ! temporary scalars
64      REAL(wp), DIMENSION(jpi,jpj) ::   &
65         zvke, zuke                     ! temporary workspace
66#endif
67      !!----------------------------------------------------------------------
68
69      IF( kt == nit000 ) THEN
70         IF(lwp) WRITE(numout,*)
71         IF(lwp) WRITE(numout,*) 'dyn_keg : kinetic energy gradient trend'
72         IF(lwp) WRITE(numout,*) '~~~~~~~'
73      ENDIF
74     
75      !                                                ! ===============
76      DO jk = 1, jpkm1                                 ! Horizontal slab
77         !                                             ! ===============
78         ! Horizontal kinetic energy at T-point
79         DO jj = 2, jpj
80            DO ji = fs_2, jpi   ! vector opt.
81               zv = 0.25 * (  vn(ji  ,jj-1,jk) * vn(ji  ,jj-1,jk)   &
82                            + vn(ji  ,jj  ,jk) * vn(ji  ,jj  ,jk)  )
83               zu = 0.25 * (  un(ji-1,jj  ,jk) * un(ji-1,jj  ,jk)   &
84                            + un(ji  ,jj  ,jk) * un(ji  ,jj  ,jk)  )
85               zhke(ji,jj,jk) = zv + zu
86#if defined key_trddyn_new
87               zvke(ji,jj) = zv
88               zuke(ji,jj) = zu
89#endif
90            END DO 
91         END DO 
92         
93         ! Horizontal gradient of Horizontal kinetic energy
94         DO jj = 2, jpjm1
95            DO ji = fs_2, fs_jpim1   ! vector opt.
96               ! gradient of kinetic energy
97               zua = -( zhke(ji+1,jj  ,jk) - zhke(ji,jj,jk) ) / e1u(ji,jj)
98               zva = -( zhke(ji  ,jj+1,jk) - zhke(ji,jj,jk) ) / e2v(ji,jj)
99               ! add to the general momentum trends
100               ua(ji,jj,jk) = ua(ji,jj,jk) + zua
101               va(ji,jj,jk) = va(ji,jj,jk) + zva
102#if defined key_trddyn
103               ! add to the general momentum trends
104               utrd(ji,jj,jk,2) = zua
105               vtrd(ji,jj,jk,2) = zva
106#endif
107#if defined key_trddyn_new
108               zuu = -( zuke(ji+1,jj  ) - zuke(ji,jj) ) / e1u(ji,jj)
109               zvv = -( zvke(ji  ,jj+1) - zvke(ji,jj) ) / e2v(ji,jj)
110               utrd(ji,jj,jk,2) = zua - zuu
111               vtrd(ji,jj,jk,3) = zva - zvv
112               utrd(ji,jj,jk,3) = zuu
113               vtrd(ji,jj,jk,2) = zvv
114#endif
115            END DO
116         END DO
117         !                                             ! ===============
118      END DO                                           !   End of slab
119      !                                                ! ===============
120
121      IF( l_ctl .AND. lwp ) THEN         ! print sum trends (used for debugging)
122         zua = SUM( ua(2:jpim1,2:jpjm1,1:jpkm1) * umask(2:jpim1,2:jpjm1,1:jpkm1) )
123         zva = SUM( va(2:jpim1,2:jpjm1,1:jpkm1) * vmask(2:jpim1,2:jpjm1,1:jpkm1) )
124         WRITE(numout,*) ' keg  - Ua: ', zua-u_ctl, ' Va: ', zva-v_ctl
125         u_ctl = zua   ;   v_ctl = zva
126      ENDIF
127
128   END SUBROUTINE dyn_keg
129
130   !!======================================================================
131END MODULE dynkeg
Note: See TracBrowser for help on using the repository browser.