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_tam.F90 in branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPATAM_SRC/DYN – NEMO

source: branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPATAM_SRC/DYN/dynldf_lap_tam.F90 @ 4582

Last change on this file since 4582 was 3611, checked in by pabouttier, 12 years ago

Add TAM code and ORCA2_TAM configuration - see Ticket #1007
  • Property svn:executable set to *
File size: 10.0 KB
Line 
1MODULE dynldf_lap_tam
2#ifdef key_tam
3   !!======================================================================
4   !!                       ***  MODULE  dynldf_lap_tam  ***
5   !! Ocean dynamics:  lateral viscosity trend
6   !!                  Tangent and Adjoint Module
7   !!======================================================================
8   !!----------------------------------------------------------------------
9   !!   dyn_ldf_lap_tan  : update the momentum trend with the lateral diffusion
10   !!                      using an iso-level harmonic operator (tangent)
11   !!   dyn_ldf_lap_adj  : update the momentum trend with the lateral diffusion
12   !!                      using an iso-level harmonic operator (adjoint)
13   !!----------------------------------------------------------------------
14   !! * Modules used
15   USE par_oce
16   USE oce_tam
17   USE ldfdyn_oce
18   USE dom_oce
19   USE in_out_manager
20   USE timing          ! Timing
21
22   IMPLICIT NONE
23   PRIVATE
24
25   !! * Routine accessibility
26   PUBLIC dyn_ldf_lap_tan  ! called by dynldf_tam.F90
27   PUBLIC dyn_ldf_lap_adj  ! called by dynldf_tam.F90
28
29   !! * Substitutions
30#  include "domzgr_substitute.h90"
31#  include "ldfdyn_substitute.h90"
32#  include "vectopt_loop_substitute.h90"
33   !!----------------------------------------------------------------------
34
35CONTAINS
36
37   SUBROUTINE dyn_ldf_lap_tan( kt )
38      !!----------------------------------------------------------------------
39      !!                     ***  ROUTINE dyn_ldf_lap_tan  ***
40      !!
41      !! ** Purpose of the direct routine:
42      !!      Compute the before horizontal tracer (t & s) diffusive
43      !!      trend and add it to the general trend of tracer equation.
44      !!
45      !! ** Method of the direct routine:
46      !!      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 lk_zco=T, e3f=e3u=e3v, the vertical scale factor are simplified
53      !!      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      !!
62      !! History of the direct routine:
63      !!        !  90-09 (G. Madec) Original code
64      !!        !  91-11 (G. Madec)
65      !!        !  96-01 (G. Madec) statement function for e3 and ahm
66      !!   8.5  !  02-06 (G. Madec)  F90: Free form and module
67      !!   9.0  !  04-08 (C. Talandier) New trends organization
68      !! History of the tangent routine
69      !!   9.0  !  08-08 (A. Vidard) tangent of 9.0
70      !!   3.4  !  12-07 (P.-A. bouttier) Phasing with 3.4
71      !!----------------------------------------------------------------------
72      !! * Arguments
73      INTEGER, INTENT( in ) ::   kt       ! ocean time-step index
74      !! * Local declarations
75      INTEGER  ::   ji, jj, jk            ! dummy loop indices
76      REAL(wp) ::   &
77         zuatl, zvatl, ze2utl, ze1vtl             ! temporary scalars
78      !!----------------------------------------------------------------------
79      !
80      IF( nn_timing == 1 )  CALL timing_start('dyn_ldf_lap_tan')
81      !
82      IF( kt == nit000 ) THEN
83         IF(lwp) WRITE(numout,*)
84         IF(lwp) WRITE(numout,*) 'dyn_ldf_tan: iso-level harmonic (laplacien) operator'
85         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ '
86      ENDIF
87
88      !                                                ! ===============
89      DO jk = 1, jpkm1                                 ! Horizontal slab
90         !                                             ! ===============
91         DO jj = 2, jpjm1
92            DO ji = fs_2, fs_jpim1   ! vector opt.
93               ! horizontal diffusive trends
94               ze2utl = rotb_tl (ji,jj,jk) * fsahmf(ji,jj,jk) * fse3f(ji,jj,jk)
95               ze1vtl = hdivb_tl(ji,jj,jk) * fsahmt(ji,jj,jk)
96               zuatl = - ( ze2utl - rotb_tl(ji  ,jj-1,jk) * fsahmf(ji  ,jj-1,jk) * fse3f(ji  ,jj-1,jk) ) &
97                  &     / ( e2u(ji,jj) * fse3u(ji,jj,jk) )  &
98                  &    + ( hdivb_tl(ji+1,jj  ,jk) * fsahmt(ji+1,jj  ,jk) - ze1vtl ) / e1u(ji,jj)
99
100               zvatl = + ( ze2utl - rotb_tl(ji-1,jj  ,jk) * fsahmf(ji-1,jj  ,jk) * fse3f(ji-1,jj  ,jk) ) &
101                  &     / ( e1v(ji,jj) * fse3v(ji,jj,jk) )  &
102                  &    + ( hdivb_tl(ji  ,jj+1,jk) * fsahmt(ji  ,jj+1,jk) - ze1vtl ) / e2v(ji,jj)
103               ! add it to the general momentum trends
104               ua_tl(ji,jj,jk) = ua_tl(ji,jj,jk) + zuatl
105               va_tl(ji,jj,jk) = va_tl(ji,jj,jk) + zvatl
106            END DO
107         END DO
108         !                                             ! ===============
109      END DO                                           !   End of slab
110      !                                                ! ===============
111      IF( nn_timing == 1 )  CALL timing_stop('dyn_ldf_lap_tan')
112      !
113   END SUBROUTINE dyn_ldf_lap_tan
114
115   SUBROUTINE dyn_ldf_lap_adj( kt )
116      !!----------------------------------------------------------------------
117      !!                     ***  ROUTINE dyn_ldf_lap_adj  ***
118      !!
119      !! ** Purpose of the direct routine:
120      !!      Compute the before horizontal tracer (t & s) diffusive
121      !!      trend and add it to the general trend of tracer equation.
122      !!
123      !! ** Method of the direct routine:
124      !!      The before horizontal momentum diffusion trend is an
125      !!      harmonic operator (laplacian type) which separates the divergent
126      !!      and rotational parts of the flow.
127      !!      Its horizontal components are computed as follow:
128      !!         difu = 1/e1u di[ahmt hdivb] - 1/(e2u*e3u) dj-1[e3f ahmf rotb]
129      !!         difv = 1/e2v dj[ahmt hdivb] + 1/(e1v*e3v) di-1[e3f ahmf rotb]
130      !!      If lk_zco=T, e3f=e3u=e3v, the vertical scale factor are simplified
131      !!      in the rotational part of the diffusion.
132      !!      Add this before trend to the general trend (ua,va):
133      !!            (ua,va) = (ua,va) + (diffu,diffv)
134      !!      'key_trddyn' activated: the two components of the horizontal
135      !!                                 diffusion trend are saved.
136      !!
137      !! ** Action : - Update (ua,va) with the before iso-level harmonic
138      !!               mixing trend.
139      !!
140      !! History of the direct routine:
141      !!        !  90-09 (G. Madec) Original code
142      !!        !  91-11 (G. Madec)
143      !!        !  96-01 (G. Madec) statement function for e3 and ahm
144      !!   8.5  !  02-06 (G. Madec)  F90: Free form and module
145      !!   9.0  !  04-08 (C. Talandier) New trends organization
146      !! History of the adjoint routine
147      !!   9.0  !  08-08 (A. Vidard) adjoint of 9.0
148      !!    -   !  09-01 (A. Weaver) misc. bug fixes and reorganization
149      !!   3.4  !  12-07 (P.-A. bouttier) Phasing with 3.4
150      !!----------------------------------------------------------------------
151      !! * Arguments
152      INTEGER, INTENT( in ) ::   kt       ! ocean time-step index
153      !! * Local declarations
154      INTEGER  ::   ji, jj, jk            ! dummy loop indices
155      REAL(wp) ::   &
156           zuaad , zvaad , ze2uad, ze1vad, &        ! temporary scalars
157         & zuaad1, zvaad1, zuaad2, zvaad2           ! temporary scalars
158      !!----------------------------------------------------------------------
159      !
160      IF( nn_timing == 1 )  CALL timing_start('dyn_ldf_lap_adj')
161      !
162      IF( kt == nitend ) THEN
163         IF(lwp) WRITE(numout,*)
164         IF(lwp) WRITE(numout,*) 'dyn_ldf_adj: iso-level harmonic (laplacien) operator'
165         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ '
166      ENDIF
167      !                                                ! ===============
168      DO jk = jpkm1, 1, -1                                ! Horizontal slab
169         !                                             ! ===============
170         DO jj = jpjm1, 2, -1
171            DO ji = fs_jpim1, fs_2, -1   ! vector opt.
172               ! add it to the general momentum trends
173               zuaad = ua_ad(ji,jj,jk)
174               zvaad = va_ad(ji,jj,jk)
175               ! horizontal diffusive trends
176               zvaad1 = zvaad /   e2v(ji,jj)
177               zvaad2 = zvaad / ( e1v(ji,jj) * fse3v(ji,jj,jk) )
178               zuaad1 = zuaad /   e1u(ji,jj)
179               zuaad2 = zuaad / ( e2u(ji,jj) * fse3u(ji,jj,jk) )
180               ze1vad = - zvaad1 - zuaad1
181               ze2uad =   zvaad2 - zuaad2
182
183               rotb_ad (ji-1,jj  ,jk) = rotb_ad (ji-1,jj  ,jk) &
184                  &                   - zvaad2 * fsahmf(ji-1,jj  ,jk) * fse3f(ji-1,jj  ,jk)
185               rotb_ad (ji  ,jj-1,jk) = rotb_ad (ji  ,jj-1,jk) &
186                  &                   + zuaad2 * fsahmf(ji  ,jj-1,jk) * fse3f(ji  ,jj-1,jk)
187               rotb_ad (ji  ,jj  ,jk) = rotb_ad (ji  ,jj  ,jk) &
188                  &                   + ze2uad * fsahmf(ji  ,jj  ,jk) * fse3f(ji  ,jj  ,jk)
189
190               hdivb_ad(ji  ,jj+1,jk) = hdivb_ad(ji  ,jj+1,jk) &
191                  &                   + zvaad1 * fsahmt(ji  ,jj+1,jk)
192               hdivb_ad(ji  ,jj  ,jk) = hdivb_ad(ji  ,jj  ,jk) &
193                  &                   + ze1vad * fsahmt(ji  ,jj  ,jk)
194               hdivb_ad(ji+1,jj  ,jk) = hdivb_ad(ji+1,jj  ,jk) &
195                  &                   + zuaad1 * fsahmt(ji+1,jj  ,jk)
196            END DO
197         END DO
198         !                                             ! ===============
199      END DO                                           !   End of slab
200      !                                                ! ===============
201      IF( nn_timing == 1 )  CALL timing_stop('dyn_ldf_lap_adj')
202      !
203   END SUBROUTINE dyn_ldf_lap_adj
204
205   !!======================================================================
206#endif
207END MODULE dynldf_lap_tam
Note: See TracBrowser for help on using the repository browser.