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.
dynzdf_exp.F90 in tags/nemo_dev_x3/NEMO/OPA_SRC/DYN – NEMO

source: tags/nemo_dev_x3/NEMO/OPA_SRC/DYN/dynzdf_exp.F90 @ 105

Last change on this file since 105 was 105, checked in by cvs2svn, 20 years ago

This commit was manufactured by cvs2svn to create tag 'nemo_dev_x3'.
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 7.3 KB
Line 
1MODULE dynzdf_exp
2   !!==============================================================================
3   !!                     ***  MODULE  dynzdf_exp  ***
4   !! Ocean dynamics:  vertical component(s) of the momentum mixing trend
5   !!==============================================================================
6
7   !!----------------------------------------------------------------------
8   !!   dyn_zdf_exp  : update the momentum trend with the vertical diffu-
9   !!                  sion using an explicit time-stepping scheme.
10   !!----------------------------------------------------------------------
11   !! * Modules used
12   USE oce             ! ocean dynamics and tracers
13   USE dom_oce         ! ocean space and time domain
14   USE phycst          ! physical constants
15   USE zdf_oce         ! ocean vertical physics
16   USE in_out_manager  ! I/O manager
17   USE taumod          ! surface ocean stress
18   USE trddyn_oce     ! dynamics trends diagnostics variables
19
20   IMPLICIT NONE
21   PRIVATE
22
23   !! * Routine accessibility
24   PUBLIC dyn_zdf_exp    ! called by step.F90
25
26   !! * Substitutions
27#  include "domzgr_substitute.h90"
28#  include "vectopt_loop_substitute.h90"
29   !!----------------------------------------------------------------------
30   !!   OPA 9.0 , LODYC-IPSL  (2003)
31   !!----------------------------------------------------------------------
32
33CONTAINS
34
35   SUBROUTINE dyn_zdf_exp( kt )
36      !!----------------------------------------------------------------------
37      !!                  ***  ROUTINE dyn_zdf_exp  ***
38      !!                   
39      !! ** Purpose :   Compute the trend due to the vert. momentum diffusion
40      !!
41      !! ** Method  :   Explicit forward time stepping with a time splitting
42      !!      technique. The vertical diffusion of momentum is given by:
43      !!         diffu = dz( avmu dz(u) ) = 1/e3u dk+1( avmu/e3uw dk(ub) )
44      !!      Surface boundary conditions: wind stress input
45      !!      Bottom boundary conditions : bottom stress (cf zdfbfr.F90)
46      !!      Add this trend to the general trend ua :
47      !!         ua = ua + dz( avmu dz(u) )
48      !!
49      !! ** Action : - Update (ua,va) with the vertical diffusive trend
50      !!             - Save the trends in (utrd,vtrd) ('key_diatrends')
51      !!
52      !! History :
53      !!        !  90-10  (B. Blanke)  Original code
54      !!        !  97-05  (G. Madec)  vertical component of isopycnal
55      !!   8.5  !  02-08  (G. Madec)  F90: Free form and module
56      !!---------------------------------------------------------------------
57      !! * Arguments
58      INTEGER, INTENT( in ) ::   kt       ! ocean time-step index
59
60      !! * Local declarations
61      INTEGER ::   ji, jj, jk, jl         ! dummy loop indices
62      REAL(wp) ::   &
63         zrau0r, zlavmr, z2dt, zua, zva   ! temporary scalars
64      REAL(wp), DIMENSION(jpi,jpk) ::   &
65         zwx, zwy, zwz, zww               ! temporary workspace arrays
66#if defined key_trddyn
67      INTEGER ::   &
68         ikbu, ikbum1 , ikbv, ikbvm1      ! temporary integers
69#endif
70      !!----------------------------------------------------------------------
71
72      IF( kt == nit000 ) THEN
73         IF(lwp) WRITE(numout,*)
74         IF(lwp) WRITE(numout,*) 'dyn_zdf_exp : vertical momentum diffusion explicit operator'
75         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ '
76      ENDIF
77
78      ! Local constant initialization
79      ! -----------------------------
80      zrau0r = 1. / rau0                                   ! inverse of the reference density
81      zlavmr = 1. / float( n_zdfexp )                      ! inverse of the number of sub time step
82      z2dt = 2. * rdt                                      ! Leap-frog environnement
83      IF( neuler == 0 .AND. kt == nit000 )   z2dt = rdt    ! Euler time stepping when starting from rest
84
85      !                                                ! ===============
86      DO jj = 2, jpjm1                                 !  Vertical slab
87         !                                             ! ===============
88
89         ! Surface boundary condition
90         DO ji = 2, jpim1
91            zwy(ji,1) = taux(ji,jj) * zrau0r
92            zww(ji,1) = tauy(ji,jj) * zrau0r
93         END DO 
94
95         ! Initialization of x, z and contingently trends array
96         DO jk = 1, jpk
97            DO ji = 2, jpim1
98               zwx(ji,jk) = ub(ji,jj,jk)
99               zwz(ji,jk) = vb(ji,jj,jk)
100#if defined key_trddyn
101               utrd(ji,jj,jk,7) = ua(ji,jj,jk)
102               vtrd(ji,jj,jk,7) = va(ji,jj,jk)
103#endif
104            END DO 
105         END DO 
106
107         ! Time splitting loop
108         DO jl = 1, n_zdfexp
109
110            ! First vertical derivative
111            DO jk = 2, jpk
112               DO ji = 2, jpim1
113                  zwy(ji,jk) = avmu(ji,jj,jk) * ( zwx(ji,jk-1) - zwx(ji,jk) ) / fse3uw(ji,jj,jk) 
114                  zww(ji,jk) = avmv(ji,jj,jk) * ( zwz(ji,jk-1) - zwz(ji,jk) ) / fse3vw(ji,jj,jk)
115               END DO 
116            END DO 
117
118            ! Second vertical derivative and trend estimation at kt+l*rdt/n_zdfexp
119            DO jk = 1, jpkm1
120               DO ji = 2, jpim1
121                  zua = zlavmr*( zwy(ji,jk) - zwy(ji,jk+1) ) / fse3u(ji,jj,jk)
122                  zva = zlavmr*( zww(ji,jk) - zww(ji,jk+1) ) / fse3v(ji,jj,jk)
123                  ua(ji,jj,jk) = ua(ji,jj,jk) + zua
124                  va(ji,jj,jk) = va(ji,jj,jk) + zva
125
126                  zwx(ji,jk) = zwx(ji,jk) + z2dt*zua*umask(ji,jj,jk)
127                  zwz(ji,jk) = zwz(ji,jk) + z2dt*zva*vmask(ji,jj,jk)
128               END DO 
129            END DO 
130
131         END DO 
132
133#if defined key_trddyn
134         ! diagnose the vertical diffusive momentum trends
135
136         ! save the total vertical momentum diffusive trend
137         DO jk = 1, jpkm1
138            DO ji = 2, jpim1
139               utrd(ji,jj,jk,7) = ua(ji,jj,jk) - utrd(ji,jj,jk,7)
140               vtrd(ji,jj,jk,7) = va(ji,jj,jk) - vtrd(ji,jj,jk,7)
141            END DO
142         END DO 
143 
144         ! subtract and save surface and momentum fluxes
145         DO ji = 2, jpim1
146            ! save the surface momentum fluxes
147            tautrd(ji,jj,1) = zwy(ji,1) / fse3u(ji,jj,1)
148            tautrd(ji,jj,2) = zww(ji,1) / fse3v(ji,jj,1)
149            ! save bottom friction momentum fluxes
150            ikbu   = MIN( mbathy(ji+1,jj), mbathy(ji,jj) )
151            ikbum1 = MAX( ikbu-1, 1 )
152            ikbv   = MIN( mbathy(ji,jj+1), mbathy(ji,jj) )
153            ikbvm1 = MAX( ikbv-1, 1 )
154            tautrd(ji,jj,3) = avmu(ji,jj,ikbu) * zwx(ji,ikbum1)   &
155                            / ( fse3u(ji,jj,ikbum1) * fse3uw(ji,jj,ikbu) )
156            tautrd(ji,jj,4) = avmv(ji,jj,ikbv) * zwz(ji,ikbvm1)   &
157                            / ( fse3v(ji,jj,ikbvm1) * fse3vw(ji,jj,ikbv) )
158            ! subtract surface forcing and bottom friction trend from vertical
159            ! diffusive momentum trend
160            utrd(ji,jj,1     ,7) = utrd(ji,jj,1     ,7) - tautrd(ji,jj,1)
161            utrd(ji,jj,ikbum1,7) = utrd(ji,jj,ikbum1,7) - tautrd(ji,jj,3)
162            vtrd(ji,jj,1     ,7) = vtrd(ji,jj,1     ,7) - tautrd(ji,jj,2)
163            vtrd(ji,jj,ikbvm1,7) = vtrd(ji,jj,ikbvm1,7) - tautrd(ji,jj,4)
164         END DO
165#endif
166         !                                             ! ===============
167      END DO                                           !   End of slab
168      !                                                ! ===============
169
170   END SUBROUTINE dyn_zdf_exp
171
172   !!==============================================================================
173END MODULE dynzdf_exp
Note: See TracBrowser for help on using the repository browser.