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

source: trunk/NEMO/OPA_SRC/DYN/dynldf_lap.F90 @ 247

Last change on this file since 247 was 247, checked in by opalod, 19 years ago

CL : Add CVS Header and CeCILL licence information
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 6.8 KB
Line 
1MODULE dynldf_lap
2   !!======================================================================
3   !!                       ***  MODULE  dynldf_lap  ***
4   !! Ocean dynamics:  lateral viscosity trend
5   !!======================================================================
6
7   !!----------------------------------------------------------------------
8   !!   dyn_ldf_lap  : update the momentum trend with the lateral diffusion
9   !!                  using an iso-level harmonic operator
10   !!----------------------------------------------------------------------
11   !! * Modules used
12   USE oce             ! ocean dynamics and tracers
13   USE dom_oce         ! ocean space and time domain
14   USE ldfdyn_oce      ! ocean dynamics: lateral physics
15   USE zdf_oce         ! ocean vertical physics
16   USE in_out_manager  ! I/O manager
17   USE trdmod          ! ocean dynamics trends
18   USE trdmod_oce      ! ocean variables trends
19   USE ldfslp          ! iso-neutral slopes
20
21   IMPLICIT NONE
22   PRIVATE
23
24   !! * Routine accessibility
25   PUBLIC dyn_ldf_lap  ! called by step.F90
26
27   !! * Substitutions
28#  include "domzgr_substitute.h90"
29#  include "ldfdyn_substitute.h90"
30#  include "vectopt_loop_substitute.h90"
31   !!----------------------------------------------------------------------
32   !!   OPA 9.0 , LOCEAN-IPSL (2005)
33   !! $Header$
34   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
35   !!----------------------------------------------------------------------
36
37CONTAINS
38
39   SUBROUTINE dyn_ldf_lap( kt )
40      !!----------------------------------------------------------------------
41      !!                     ***  ROUTINE dyn_ldf_lap  ***
42      !!                       
43      !! ** Purpose :   Compute the before horizontal tracer (t & s) diffusive
44      !!      trend and add it to the general trend of tracer equation.
45      !!
46      !! ** Method  :   The before horizontal momentum diffusion trend is an
47      !!      harmonic operator (laplacian type) which separates the divergent
48      !!      and rotational parts of the flow.
49      !!      Its horizontal components are computed as follow:
50      !!         difu = 1/e1u di[ahmt hdivb] - 1/(e2u*e3u) dj-1[e3f ahmf rotb]
51      !!         difv = 1/e2v dj[ahmt hdivb] + 1/(e1v*e3v) di-1[e3f ahmf rotb]
52      !!      If 'key_s_coord' key is not activated, the vertical scale factor
53      !!      is simplified in the rotational part of the diffusion.
54      !!      Add this before trend to the general trend (ua,va):
55      !!            (ua,va) = (ua,va) + (diffu,diffv)
56      !!      'key_trddyn' activated: the two components of the horizontal
57      !!                                 diffusion trend are saved.
58      !!
59      !! ** Action : - Update (ua,va) with the before iso-level harmonic
60      !!               mixing trend.
61      !!             - Save in (ztdua,ztdva) arrays the trends ('key_trddyn')
62      !!
63      !! History :
64      !!        !  90-09 (G. Madec) Original code
65      !!        !  91-11 (G. Madec)
66      !!        !  96-01 (G. Madec) statement function for e3 and ahm
67      !!   8.5  !  02-06 (G. Madec)  F90: Free form and module
68      !!   9.0  !  04-08 (C. Talandier) New trends organization
69      !!----------------------------------------------------------------------
70      !! * Modules used     
71      USE oce, ONLY :    ztdua => ta,   & ! use ta as 3D workspace   
72                         ztdva => sa      ! use sa as 3D workspace   
73
74      !! * Arguments
75      INTEGER, INTENT( in ) ::   kt       ! ocean time-step index
76
77      !! * Local declarations
78      INTEGER  ::   ji, jj, jk            ! dummy loop indices
79      REAL(wp) ::   &
80         zua, zva, ze2u, ze1v             ! temporary scalars
81      !!----------------------------------------------------------------------
82
83      IF( kt == nit000 ) THEN
84         IF(lwp) WRITE(numout,*)
85         IF(lwp) WRITE(numout,*) 'dyn_ldf : iso-level harmonic (laplacien) operator'
86         IF(lwp) WRITE(numout,*) '~~~~~~~ '
87      ENDIF
88
89      ! Save ua and va trends
90      IF( l_trddyn )   THEN
91         ztdua(:,:,:) = ua(:,:,:) 
92         ztdva(:,:,:) = va(:,:,:) 
93      ENDIF
94
95      !                                                ! ===============
96      DO jk = 1, jpkm1                                 ! Horizontal slab
97         !                                             ! ===============
98         DO jj = 2, jpjm1
99            DO ji = fs_2, fs_jpim1   ! vector opt.
100#if defined key_s_coord || defined key_partial_steps
101               ze2u = rotb (ji,jj,jk)*fsahmf(ji,jj,jk)*fse3f(ji,jj,jk)
102               ze1v = hdivb(ji,jj,jk)*fsahmt(ji,jj,jk)
103               ! horizontal diffusive trends
104               zua = - ( ze2u - rotb (ji,jj-1,jk)*fsahmf(ji,jj-1,jk)*fse3f(ji,jj-1,jk) ) / ( e2u(ji,jj) * fse3u(ji,jj,jk) )   &
105                     + ( hdivb(ji+1,jj,jk)*fsahmt(ji+1,jj,jk) - ze1v                   ) / e1u(ji,jj)
106
107               zva = + ( ze2u - rotb (ji-1,jj,jk)*fsahmf(ji-1,jj,jk)*fse3f(ji-1,jj,jk) ) / ( e1v(ji,jj) * fse3v(ji,jj,jk) )   &
108                     + ( hdivb(ji,jj+1,jk)*fsahmt(ji,jj+1,jk) - ze1v                   ) / e2v(ji,jj)
109#else
110               ! horizontal diffusive trends
111               ze2u = rotb (ji,jj,jk)*fsahmf(ji,jj,jk)
112               ze1v = hdivb(ji,jj,jk)*fsahmt(ji,jj,jk)
113               zua = - (                ze2u                  - rotb (ji,jj-1,jk)*fsahmf(ji,jj-1,jk) ) / e2u(ji,jj)   &
114                     + ( hdivb(ji+1,jj,jk)*fsahmt(ji+1,jj,jk) -                ze1v                  ) / e1u(ji,jj)
115
116               zva = + (                ze2u                  - rotb (ji-1,jj,jk)*fsahmf(ji-1,jj,jk) ) / e1v(ji,jj)   &
117                     + ( hdivb(ji,jj+1,jk)*fsahmt(ji,jj+1,jk) -                ze1v                  ) / e2v(ji,jj)
118#endif
119
120               ! add it to the general momentum trends
121               ua(ji,jj,jk) = ua(ji,jj,jk) + zua
122               va(ji,jj,jk) = va(ji,jj,jk) + zva
123            END DO
124         END DO
125         !                                             ! ===============
126      END DO                                           !   End of slab
127      !                                                ! ===============
128
129      ! save the lateral diffusion trends for diagnostic
130      ! momentum trends
131      IF( l_trddyn )   THEN
132         ztdua(:,:,:) = ua(:,:,:) - ztdua(:,:,:)
133         ztdva(:,:,:) = va(:,:,:) - ztdva(:,:,:)
134
135         CALL trd_mod(ztdua, ztdva, jpdtdldf, 'DYN', kt)
136      ENDIF
137
138      IF(l_ctl) THEN         ! print sum trends (used for debugging)
139         zua = SUM( ua(2:nictl,2:njctl,1:jpkm1) * umask(2:nictl,2:njctl,1:jpkm1) )
140         zva = SUM( va(2:nictl,2:njctl,1:jpkm1) * vmask(2:nictl,2:njctl,1:jpkm1) )
141         WRITE(numout,*) ' ldf  - Ua: ', zua-u_ctl, ' Va: ', zva-v_ctl
142         u_ctl = zua   ;   v_ctl = zva
143      ENDIF
144
145   END SUBROUTINE dyn_ldf_lap
146
147   !!======================================================================
148END MODULE dynldf_lap
Note: See TracBrowser for help on using the repository browser.