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.
trcldf_iso.F90 in trunk/NEMO/TOP_SRC/TRP – NEMO

source: trunk/NEMO/TOP_SRC/TRP/trcldf_iso.F90 @ 941

Last change on this file since 941 was 941, checked in by cetlod, 16 years ago

phasing the passive tracer transport module to the new version of NEMO, see ticket 143
  • Property svn:executable set to *
File size: 11.5 KB
RevLine 
[941]1MODULE trcldf_iso
2   !!==============================================================================
3   !!                    ***  MODULE  trcldf_iso  ***
4   !! Ocean passive tracers:  horizontal component of the lateral tracer mixing trend
5   !!==============================================================================
6#if defined key_top && defined key_ldfslp 
7   !!----------------------------------------------------------------------
8   !!   'key_top'          and                                   TOP models
9   !!   'key_ldfslp'                  rotation of the lateral mixing tensor
10   !!----------------------------------------------------------------------
11   !!   trc_ldf_iso : update the tracer trend with the horizontal component
12   !!                 of iso neutral laplacian operator or horizontal
13   !!                 laplacian operator in s-coordinate
14   !!----------------------------------------------------------------------
15   !! * Modules used
16   USE oce_trc      ! ocean dynamics and tracers variables
17   USE trc          ! ocean passive tracers variables
18   USE prtctl_trc   ! Print control for debbuging
19
20   IMPLICIT NONE
21   PRIVATE
22
23   !! * Routine accessibility
24   PUBLIC trc_ldf_iso  ! routine called by step.F90
25
26   !! * Substitutions
27#  include "top_substitute.h90"
28   !!----------------------------------------------------------------------
29   !!   TOP 1.0 , LOCEAN-IPSL (2005)
30   !! $Header: /home/opalod/NEMOCVSROOT/NEMO/TOP_SRC/TRP/trcldf_iso.F90,v 1.10 2007/10/12 09:26:30 opalod Exp $
31   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
32   !!----------------------------------------------------------------------
33
34CONTAINS
35
36   SUBROUTINE trc_ldf_iso( kt )
37      !!----------------------------------------------------------------------
38      !!                  ***  ROUTINE trc_ldf_iso  ***
39      !!
40      !! ** Purpose :   Compute the before horizontal tracer  diffusive
41      !!      trend and add it to the general trend of tracer equation.
42      !!
43      !! ** Method  :   The horizontal component of the lateral diffusive trends
44      !!      is provided by a 2nd order operator rotated along neural or geopo-
45      !!      tential surfaces to which an eddy induced advection can be added
46      !!      It is computed using before fields (forward in time) and isopyc-
47      !!      nal or geopotential slopes computed in routine ldfslp.
48      !!
49      !!      horizontal fluxes associated with the rotated lateral mixing:
50      !!         zftu = (aht+ahtb0) e2u*e3u/e1u di[ tb ]
51      !!               - aht       e2u*uslp    dk[ mi(mk(tb)) ]
52      !!         zftv = (aht+ahtb0) e1v*e3v/e2v dj[ tb ]
53      !!               - aht       e2u*vslp    dk[ mj(mk(tb)) ]
54      !!      add horizontal Eddy Induced advective fluxes (lk_traldf_eiv=T):
55      !!         zftu = zftu - dk-1[ aht e2u mi(wslpi) ] mi( tb )
56      !!         zftv = zftv - dk-1[ aht e1v mj(wslpj) ] mj( tb )
57      !!      take the horizontal divergence of the fluxes:
58      !!         difft = 1/(e1t*e2t*e3t) {  di-1[ zftu ] +  dj-1[ zftv ]  }
59      !!      Add this trend to the general trend tra :
60      !!         tra = tra + difft
61      !!
62      !! ** Action  : - Update tra arrays with the before isopycnal or
63      !!                geopotential s-coord harmonic mixing trend.
64      !!              - Save the trends in trtrd ('key_trc_diatrd')
65      !!
66      !! History :
67      !!        !  94-08  (G. Madec, M. Imbard)
68      !!        !  97-05  (G. Madec)  split into traldf and trazdf
69      !!        !  98-03  (L. Bopp, MA Foujols) passive tracer generalisation
70      !!        !  00-10  (MA Foujols E Kestenare) USE passive tracer coefficient
71      !!   8.5  !  02-08  (G. Madec)  Free form, F90
72      !!   9.0  !  04-03  (C. Ethe)  Free form, F90
73      !!----------------------------------------------------------------------
74      !! * Modules used
75      USE oce_trc       , zftu => ua,  &  ! use ua as workspace
76         &                zfsu => va      ! use va as workspace
77
78      !! * Arguments
79      INTEGER, INTENT( in ) ::   kt       ! ocean time-step index
80
81      !! * Local declarations
82      INTEGER ::   ji, jj, jk,jn             ! dummy loop indices
83      REAL(wp) ::   &
84         zabe1, zabe2, zcof1, zcof2,   &  ! temporary scalars
85         zmsku, zmskv, zbtr,           &
86#if defined key_trcldf_eiv
87         zcg1, zcg2, zuwk, zvwk,       &
88         zuwk1, zvwk1,                 &
89#endif
90         ztra
91
92      REAL(wp), DIMENSION(jpi,jpj) ::   &
93         zdkt, zdk1t            ! workspace
94
95#if defined key_trcldf_eiv
96      REAL(wp), DIMENSION(jpi,jpj) ::   &
97         zftug, zftvg
98
99#if defined key_trc_diatrd
100      REAL(wp) ::   &
101         ztagu, ztagv
102#endif
103
104#endif
105
106      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
107         zftv                       ! workspace
108      CHARACTER (len=22) :: charout
109      !!----------------------------------------------------------------------
110
111      IF( kt == nittrc000 ) THEN
112         IF(lwp) WRITE(numout,*)
113         IF(lwp) WRITE(numout,*) 'trc_ldf_iso : iso neutral lateral diffusion or'
114         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~   horizontal laplacian diffusion in s-coordinate'
115#if defined key_trcldf_eiv && defined key_diaeiv
116         u_trc_eiv(:,:,:) = 0.e0
117         v_trc_eiv(:,:,:) = 0.e0
118#endif
119      ENDIF
120
121
122      DO jn = 1, jptra
123
124         !                                                ! ===============
125         DO jk = 1, jpkm1                                 ! Horizontal slab
126            !                                             ! ===============
127            ! 1. Vertical tracer gradient at level jk and jk+1
128            ! ------------------------------------------------
129            ! surface boundary condition: zdkt(jk=1)=zdkt(jk=2)
130
131            zdk1t(:,:) = ( trb(:,:,jk,jn) - trb(:,:,jk+1,jn) ) * tmask(:,:,jk+1)
132
133            IF( jk == 1 ) THEN
134               zdkt(:,:) = zdk1t(:,:)
135            ELSE
136               zdkt(:,:) = ( trb(:,:,jk-1,jn) - trb(:,:,jk,jn) ) * tmask(:,:,jk)
137            ENDIF
138
139
140            ! 2. Horizontal fluxes
141            ! --------------------
142
143            DO jj = 1 , jpjm1
144               DO ji = 1, fs_jpim1   ! vector opt.
145                  zabe1 = ( fsahtru(ji,jj,jk) + ahtrb0 ) * e2u(ji,jj) * fse3u(ji,jj,jk) / e1u(ji,jj)
146                  zabe2 = ( fsahtrv(ji,jj,jk) + ahtrb0 ) * e1v(ji,jj) * fse3v(ji,jj,jk) / e2v(ji,jj)
147
148                  zmsku = 1. / MAX(   tmask(ji+1,jj,jk  ) + tmask(ji,jj,jk+1)   &
149                     + tmask(ji+1,jj,jk+1) + tmask(ji,jj,jk  ), 1. )
150
151                  zmskv = 1. / MAX(   tmask(ji,jj+1,jk  ) + tmask(ji,jj,jk+1)   &
152                     + tmask(ji,jj+1,jk+1) + tmask(ji,jj,jk  ), 1. )
153
154                  zcof1 = -fsahtru(ji,jj,jk) * e2u(ji,jj) * uslp(ji,jj,jk) * zmsku
155                  zcof2 = -fsahtrv(ji,jj,jk) * e1v(ji,jj) * vslp(ji,jj,jk) * zmskv
156
157                  zftu(ji,jj,jk) = umask(ji,jj,jk) * (   zabe1 * (   trb(ji+1,jj,jk,jn) - trb(ji,jj,jk,jn)  )   &
158                     &                              + zcof1 * (   zdkt (ji+1,jj) + zdk1t(ji,jj)      &
159                     &                                          + zdk1t(ji+1,jj) + zdkt (ji,jj)  )  )
160
161                  zftv(ji,jj,jk) = vmask(ji,jj,jk) * (   zabe2 * (   trb(ji,jj+1,jk,jn) - trb(ji,jj,jk,jn)  )   &
162                     &                              + zcof2 * (   zdkt (ji,jj+1) + zdk1t(ji,jj)      &
163                     &                                          + zdk1t(ji,jj+1) + zdkt (ji,jj)  )  )
164
165               END DO
166            END DO
167
168#   if defined key_trcldf_eiv
169            !                              ! ---------------------------------------!
170            !                              ! Eddy induced vertical advective fluxes !
171            !                              ! ---------------------------------------!
172            DO jj = 1, jpjm1
173               DO ji = 1, fs_jpim1   ! vector opt.
174                  zuwk = ( wslpi(ji,jj,jk  ) + wslpi(ji+1,jj,jk  ) ) * fsaeitru(ji,jj,jk  ) * umask(ji,jj,jk  )
175                  zuwk1= ( wslpi(ji,jj,jk+1) + wslpi(ji+1,jj,jk+1) ) * fsaeitru(ji,jj,jk+1) * umask(ji,jj,jk+1)
176                  zvwk = ( wslpj(ji,jj,jk  ) + wslpj(ji,jj+1,jk  ) ) * fsaeitrv(ji,jj,jk  ) * vmask(ji,jj,jk  )
177                  zvwk1= ( wslpj(ji,jj,jk+1) + wslpj(ji,jj+1,jk+1) ) * fsaeitrv(ji,jj,jk+1) * vmask(ji,jj,jk+1)
178
179                  zcg1= -0.25 * e2u(ji,jj) * umask(ji,jj,jk) * ( zuwk-zuwk1 )
180                  zcg2= -0.25 * e1v(ji,jj) * vmask(ji,jj,jk) * ( zvwk-zvwk1 )
181
182                  zftug(ji,jj) = zcg1 * ( trb(ji+1,jj,jk,jn) + trb(ji,jj,jk,jn) )
183                  zftvg(ji,jj) = zcg2 * ( trb(ji,jj+1,jk,jn) + trb(ji,jj,jk,jn) )
184
185                  zftu(ji,jj,jk) = zftu(ji,jj,jk) + zftug(ji,jj)
186                  zftv(ji,jj,jk) = zftv(ji,jj,jk) + zftvg(ji,jj)
187
188#   if defined key_diaeiv
189                  u_trc_eiv(ji,jj,jk) = -2. * zcg1 / ( e2u(ji,jj) * fse3u(ji,jj,jk) )
190                  v_trc_eiv(ji,jj,jk) = -2. * zcg2 / ( e1v(ji,jj) * fse3v(ji,jj,jk) )
191#   endif
192               END DO
193            END DO
194#   endif
195
196            ! II.4 Second derivative (divergence) and add to the general trend
197            ! ----------------------------------------------------------------
198
199            DO jj = 2 , jpjm1
200               DO ji = fs_2, fs_jpim1   ! vector opt.
201                  zbtr= 1. / ( e1t(ji,jj)*e2t(ji,jj)*fse3t(ji,jj,jk) )
202                  ztra = zbtr * (  zftu(ji,jj,jk) - zftu(ji-1,jj  ,jk)   &
203                     &          + zftv(ji,jj,jk) - zftv(ji  ,jj-1,jk)  )
204                  tra (ji,jj,jk,jn) = tra (ji,jj,jk,jn) + ztra
205#if defined key_trc_diatrd
206                  IF (luttrd(jn)) trtrd (ji,jj,jk,ikeep(jn),4) = ( zftu(ji,jj,jk) - zftu(ji-1,jj,jk  ) ) * zbtr
207                  IF (luttrd(jn)) trtrd (ji,jj,jk,ikeep(jn),5) = ( zftv(ji,jj,jk) - zftv(ji,jj-1,jk  ) ) * zbtr
208#endif
209               END DO
210            END DO
211
212#if defined key_trc_diatrd
213#   if defined key_trcldf_eiv
214            DO jj = 2 , jpjm1
215               DO ji = fs_2, fs_jpim1   ! vector opt.
216                  zbtr= 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
217                  ztagu = ( zftug(ji,jj) - zftug(ji-1,jj  ) ) * zbtr
218                  ztagv = ( zftvg(ji,jj) - zftvg(ji  ,jj-1) ) * zbtr
219                  IF (luttrd(jn)) trtrd (ji,jj,jk,ikeep(jn),4) = trtrd(ji,jj,jk,ikeep(jn),4) - ztagu
220                  IF (luttrd(jn)) trtrd (ji,jj,jk,ikeep(jn),5) = trtrd(ji,jj,jk,ikeep(jn),5) - ztagv
221                  IF (luttrd(jn)) trtrd (ji,jj,jk,ikeep(jn),7) = ztagu
222                  IF (luttrd(jn)) trtrd (ji,jj,jk,ikeep(jn),8) = ztagv
223               END DO
224             END DO
225#   endif
226#endif
227
228            !                                          ! ===============
229         END DO                                        !   End of slab 
230         !                                             ! ===============
231
232      END DO
233
234      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
235         WRITE(charout, FMT="('ldf - iso')")
236         CALL prt_ctl_trc_info(charout)
237         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd')
238      ENDIF
239
240   END SUBROUTINE trc_ldf_iso
241
242#else
243   !!----------------------------------------------------------------------
244   !!   Dummy module :             No rotation of the lateral mixing tensor
245   !!----------------------------------------------------------------------
246CONTAINS
247   SUBROUTINE trc_ldf_iso( kt )               ! Empty routine
248      INTEGER, INTENT(in) :: kt
249      WRITE(*,*) 'trc_ldf_iso: You should not have seen this print! error?', kt
250   END SUBROUTINE trc_ldf_iso
251#endif
252
253   !!==============================================================================
254END MODULE trcldf_iso
Note: See TracBrowser for help on using the repository browser.