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

source: trunk/NEMO/OPA_SRC/TRA/trazdf_imp.F90 @ 216

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

CT : UPDATE151 : New trends organization
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 9.0 KB
Line 
1MODULE trazdf_imp
2   !!==============================================================================
3   !!                    ***  MODULE  trazdf_imp  ***
4   !! Ocean active tracers:  vertical component of the tracer mixing trend
5   !!==============================================================================
6
7   !!----------------------------------------------------------------------
8   !!   tra_zdf_imp  : update the tracer trend with the vertical diffusion
9   !!                  using an implicit time-stepping.
10   !!----------------------------------------------------------------------
11   !! * Modules used
12   USE oce             ! ocean dynamics and active tracers
13   USE dom_oce         ! ocean space and time domain
14   USE zdf_oce         ! ocean vertical physics
15   USE ldftra_oce      ! ocean active tracers: lateral physics
16   USE zdfddm          ! ocean vertical physics: double diffusion
17   USE trdmod          ! ocean active tracers trends
18   USE trdmod_oce      ! ocean variables trends
19   USE in_out_manager  ! I/O manager
20
21
22   IMPLICIT NONE
23   PRIVATE
24
25   !! * Routine accessibility
26   PUBLIC tra_zdf_imp          ! routine called by step.F90
27
28   !! * Module variable
29   REAL(wp), DIMENSION(jpk) ::   &
30      r2dt                     ! vertical profile of 2 x tracer time-step
31
32   !! * Substitutions
33#  include "domzgr_substitute.h90"
34#  include "zdfddm_substitute.h90"
35   !!----------------------------------------------------------------------
36   !!   OPA 9.0 , LODYC-IPSL (2003)
37   !!----------------------------------------------------------------------
38
39CONTAINS
40
41   SUBROUTINE tra_zdf_imp( kt )
42      !!----------------------------------------------------------------------
43      !!                  ***  ROUTINE tra_zdf_imp  ***
44      !!
45      !! ** Purpose :   Compute the trend due to the vertical tracer mixing
46      !!      using an implicit time stepping and add it to the general trend
47      !!      of the tracer equations.
48      !!
49      !! ** Method  :   The vertical diffusion of tracers (t & s) is given by:
50      !!          difft = dz( avt dz(t) ) = 1/e3t dk+1( avt/e3w dk(ta) )
51      !!      It is thus evaluated using a backward time scheme
52      !!      Surface and bottom boundary conditions: no diffusive flux on
53      !!      both tracers (bottom, applied through the masked field avt).
54      !!      Add this trend to the general trend ta,sa :
55      !!          ta = ta + dz( avt dz(t) )
56      !!         (sa = sa + dz( avs dz(t) ) if lk_zdfddm=T)
57      !!
58      !! ** Action  : - Update (ta,sa) with the before vertical diffusion trend
59      !!              - save the trends in (ttrd,strd) ('key_trdtra')
60      !!
61      !! History :
62      !!   6.0  !  90-10  (B. Blanke)  Original code
63      !!   7.0  !  91-11  (G. Madec)
64      !!        !  92-06  (M. Imbard)  correction on tracer trend loops
65      !!        !  96-01  (G. Madec)  statement function for e3
66      !!        !  97-05  (G. Madec)  vertical component of isopycnal
67      !!        !  97-07  (G. Madec)  geopotential diffusion in s-coord
68      !!        !  00-08  (G. Madec)  double diffusive mixing
69      !!   8.5  !  02-08  (G. Madec)  F90: Free form and module
70      !!   9.0  !  04-08  (C. Talandier) New trends organization
71      !!---------------------------------------------------------------------
72      !! * Modules used     
73      USE oce, ONLY :    ztdta => ua,       & ! use ua as 3D workspace   
74                         ztdsa => va          ! use va as 3D workspace   
75
76      !! * Arguments
77      INTEGER, INTENT( in ) ::   kt           ! ocean time-step index
78
79      !! * Local declarations
80      INTEGER ::   ji, jj, jk                 ! dummy loop indices
81      INTEGER ::   ikst, ikenm2, ikstp1
82      REAL(wp)  ::   zta, zsa                 ! temporary scalars
83      REAL(wp), DIMENSION(jpi,jpk) ::   &
84         zwd, zws, zwi,          &  ! ???
85         zwx, zwy, zwz, zwt         ! ???
86      !!---------------------------------------------------------------------
87
88
89      ! 0. Local constant initialization
90      ! --------------------------------
91
92      ! time step = 2 rdttra ex
93      IF( neuler == 0 .AND. kt == nit000 ) THEN
94         r2dt(:) =  rdttra(:)              ! restarting with Euler time stepping
95      ELSEIF( kt <= nit000 + 1) THEN
96         r2dt(:) = 2. * rdttra(:)          ! leapfrog
97      ENDIF
98
99      ! Save ta and sa trends
100      IF( l_trdtra )   THEN
101         ztdta(:,:,:) = ta(:,:,:) 
102         ztdsa(:,:,:) = sa(:,:,:) 
103      ENDIF
104
105      !                                                ! ===============
106      DO jj = 2, jpjm1                                 !  Vertical slab
107         !                                             ! ===============
108         ! 0. Matrix construction
109         ! ----------------------
110
111         ! Diagonal, inferior, superior (including the bottom boundary condition via avt masked)
112         DO jk = 1, jpkm1
113            DO ji = 2, jpim1
114               zwi(ji,jk) = - r2dt(jk) * avt(ji,jj,jk  )   &
115                                       / ( fse3t(ji,jj,jk) * fse3w(ji,jj,jk  ) )
116               zws(ji,jk) = - r2dt(jk) * avt(ji,jj,jk+1)   &
117                                       / ( fse3t(ji,jj,jk) * fse3w(ji,jj,jk+1) )
118               zwd(ji,jk) = 1. - zwi(ji,jk) - zws(ji,jk)
119            END DO
120         END DO
121
122         ! Surface boudary conditions
123         DO ji = 2, jpim1
124            zwi(ji,1) = 0.e0
125            zwd(ji,1) = 1. - zws(ji,1)
126         END DO
127
128         ! 1. Vertical diffusion on temperature
129         ! -------------------------===========
130
131         ! Second member construction
132         DO jk = 1, jpkm1
133            DO ji = 2, jpim1
134               zwy(ji,jk) = tb(ji,jj,jk) + r2dt(jk) * ta(ji,jj,jk)
135            END DO
136         END DO
137
138         ! Matrix inversion from the first level
139         ikst = 1
140
141#   include "zdf.matrixsolver.h90"
142
143         ! Save the masked temperature after in ta
144         ! (c a u t i o n:  temperature not its trend, Leap-frog scheme done
145         !                  it will not be done in tranxt)
146         DO jk = 1, jpkm1
147            DO ji = 2, jpim1
148               ta(ji,jj,jk) = zwx(ji,jk) * tmask(ji,jj,jk)
149            END DO
150         END DO
151
152
153         ! 2. Vertical diffusion on salinity
154         ! -------------------------========
155
156#if defined key_zdfddm
157         ! Rebuild the Matrix as avt /= avs
158
159         ! Diagonal, inferior, superior
160         ! (including the bottom boundary condition via avs masked)
161         DO jk = 1, jpkm1
162            DO ji = 2, jpim1
163               zwi(ji,jk) = - r2dt(jk) * fsavs(ji,jj,jk  )   &
164                  /( fse3t(ji,jj,jk) * fse3w(ji,jj,jk  ) )
165               zws(ji,jk) = - r2dt(jk) * fsavs(ji,jj,jk+1)   &
166                  /( fse3t(ji,jj,jk) * fse3w(ji,jj,jk+1) )
167               zwd(ji,jk) = 1. - zwi(ji,jk) - zws(ji,jk)
168            END DO
169         END DO
170
171         ! Surface boudary conditions
172         DO ji = 2, jpim1
173            zwi(ji,1) = 0.e0
174            zwd(ji,1) = 1. - zws(ji,1)
175         END DO
176#endif
177         ! Second member construction
178         DO jk = 1, jpkm1
179            DO ji = 2, jpim1
180               zwy(ji,jk) = sb(ji,jj,jk) + r2dt(jk) * sa(ji,jj,jk)
181            END DO
182         END DO
183
184         ! Matrix inversion from the first level
185         ikst = 1
186
187#   include "zdf.matrixsolver.h90"
188
189
190         ! Save the masked salinity after in sa
191         ! (c a u t i o n:  salinity not its trend, Leap-frog scheme done
192         !                  it will not be done in tranxt)
193         DO jk = 1, jpkm1
194            DO ji = 2, jpim1
195               sa(ji,jj,jk) = zwx(ji,jk)  * tmask(ji,jj,jk)
196            END DO
197         END DO
198
199         !                                             ! ===============
200      END DO                                           !   End of slab
201      !                                                ! ===============
202
203      ! save the trends for diagnostic
204      ! Compute and save the vertical diffusive temperature & salinity trends
205      IF( l_trdtra )   THEN
206         ! compute the vertical diffusive trends in substracting the previous
207         ! trends ztdta()/ztdsa() to the new one computed (dT/dt or dS/dt)
208         ! with the new temperature/salinity ta/sa
209         DO jk = 1, jpkm1
210            ztdta(:,:,jk) = ( ( ta(:,:,jk) - tb(:,:,jk) ) / r2dt(jk) )   & ! new trend
211                &           - ztdta(:,:,jk)                                ! old trend
212            ztdsa(:,:,jk) = ( ( sa(:,:,jk) - sb(:,:,jk) ) / r2dt(jk) )   & ! new trend
213                &           - ztdsa(:,:,jk)                                ! old trend
214         END DO
215
216         CALL trd_mod(ztdta, ztdsa, jpttdzdf, 'TRA', kt)
217      ENDIF
218
219      IF(l_ctl) THEN         ! print mean trends (used for debugging)
220         zta = SUM( ta(2:nictl,2:njctl,1:jpkm1) * tmask(2:nictl,2:njctl,1:jpkm1) )
221         zsa = SUM( sa(2:nictl,2:njctl,1:jpkm1) * tmask(2:nictl,2:njctl,1:jpkm1) )
222         WRITE(numout,*) ' zdf  - Ta: ', zta-t_ctl, ' Sa: ', zsa-s_ctl
223         t_ctl = zta   ;   s_ctl = zsa
224      ENDIF
225
226   END SUBROUTINE tra_zdf_imp
227
228   !!==============================================================================
229END MODULE trazdf_imp
Note: See TracBrowser for help on using the repository browser.