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

source: trunk/NEMO/OPA_SRC/TRA/tranpc.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: 11.3 KB
Line 
1MODULE tranpc
2   !!==============================================================================
3   !!                       ***  MODULE  tranpc  ***
4   !! Ocean active tracers:  non penetrative convection scheme
5   !!==============================================================================
6
7   !!----------------------------------------------------------------------
8   !!   tra_npc      : apply the non penetrative convection scheme
9   !!   tra_npc_init : initialization and control of the scheme
10   !!----------------------------------------------------------------------
11   !! * Modules used
12   USE oce             ! ocean dynamics and active tracers
13   USE dom_oce         ! ocean space and time domain
14   USE trdmod          ! ocean active tracer trends
15   USE trdmod_oce      ! ocean variables trends
16   USE eosbn2          ! equation of state (eos routine)
17   USE lbclnk          ! lateral boundary conditions (or mpp link)
18   USE in_out_manager  ! I/O manager
19
20   IMPLICIT NONE
21   PRIVATE
22
23   !! * Routine accessibility
24   PUBLIC tra_npc      ! routine called by step.F90
25
26   !! * Module variable
27   INTEGER ::       &
28      nnpc1 =   1,  &  ! nnpc1   non penetrative convective scheme frequency
29      nnpc2 =  15      ! nnpc2   non penetrative convective scheme print frequency
30
31   !! * Substitutions
32#  include "domzgr_substitute.h90"
33   !!----------------------------------------------------------------------
34   !!   OPA 9.0 , LODYC-IPSL (2003)
35   !!----------------------------------------------------------------------
36
37CONTAINS
38
39   SUBROUTINE tra_npc( kt )
40      !!----------------------------------------------------------------------
41      !!                  ***  ROUTINE tranpc  ***
42      !!
43      !! ** Purpose :   Non penetrative convective adjustment scheme. solve
44      !!      the static instability of the water column (now, after the swap)
45      !!      while conserving heat and salt contents.
46      !!
47      !! ** Method  :   The algorithm used converges in a maximium of jpk
48      !!      iterations. instabilities are treated when the vertical density
49      !!      gradient is less than 1.e-5.
50      !!
51      !!      'key_trdtra' defined: the trend associated with this
52      !!                               algorithm is saved.
53      !!
54      !!      macro-tasked on vertical slab (jj-loop)
55      !!
56      !! ** Action  : - (tn,sn) after the application od the npc scheme
57      !!              - save the associated trends (ttrd,strd) ('key_trdtra')
58      !!
59      !! References :
60      !!      Madec, et al., 1991, JPO, 21, 9, 1349-1371.
61      !!
62      !! History :
63      !!   1.0  !  90-09  (G. Madec)  Original code
64      !!        !  91-11  (G. Madec)
65      !!        !  92-06  (M. Imbard)  periodic conditions on t and s
66      !!        !  93-03  (M. Guyon)  symetrical conditions
67      !!        !  96-01  (G. Madec)  statement function for e3
68      !!                                  suppression of common work arrays
69      !!   8.5  !  02-06  (G. Madec)  free form F90
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 ::   &
82         inpcc ,                        & ! number of statically instable water column
83         inpci ,                        & ! number of iteration for npc scheme
84         jiter, jkdown, jkp,            & ! ???
85         ikbot, ik, ikup, ikdown          ! ???
86      REAL(wp) ::   &                     ! temporary arrays
87         ze3tot, zta, zsa, zraua, ze3dwn
88      REAL(wp), DIMENSION(jpi,jpk) ::   &
89         zwx, zwy, zwz                    ! temporary arrays
90      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
91         zrhop                            ! temporary arrays
92      !!----------------------------------------------------------------------
93
94      IF( kt == nit000  )   CALL tra_npc_init
95
96
97      IF( MOD( kt, nnpc1 ) == 0 ) THEN
98
99         inpcc = 0
100         inpci = 0
101
102         ! 0. Potential density
103         ! --------------------
104
105         CALL eos( tn, sn, rhd, zrhop )
106
107         ! Save tn and sn trends
108         IF( l_trdtra )   THEN
109            ztdta(:,:,:) = tn(:,:,:) 
110            ztdsa(:,:,:) = sn(:,:,:) 
111         ENDIF
112
113         !                                                ! ===============
114         DO jj = 1, jpj                                   !  Vertical slab
115            !                                             ! ===============
116
117            ! 1. Static instability pointer
118            ! -----------------------------
119
120            DO jk = 1, jpkm1
121               DO ji = 1, jpi
122                  zwx(ji,jk) = ( zrhop(ji,jj,jk) - zrhop(ji,jj,jk+1) ) * tmask(ji,jj,jk+1)
123               END DO
124            END DO
125
126            ! 1.1 do not consider the boundary points
127
128            ! even if east-west cyclic b. c. do not considere ji=1 or jpi
129            DO jk = 1, jpkm1
130               zwx( 1 ,jk) = 0.e0
131               zwx(jpi,jk) = 0.e0
132            END DO
133            ! even if south-symmetric b. c. used, do not considere jj=1
134            IF( jj == 1 ) zwx(:,:) = 0.e0
135
136            DO jk = 1, jpkm1
137               DO ji = 1, jpi
138                  zwx(ji,jk) = 1.
139                  IF( zwx(ji,jk) < 1.e-5 ) zwx(ji,jk)=0.
140               END DO
141            END DO
142
143            zwy(:,1) = 0.
144            DO ji = 1, jpi
145               DO jk = 1, jpkm1
146                  zwy(ji,1) = zwy(ji,1) + zwx(ji,jk)
147               END DO
148            END DO
149
150            zwz(1,1) = 0.
151            DO ji = 1, jpi
152               zwz(1,1) = zwz(1,1) + zwy(ji,1)
153            END DO
154
155            inpcc = inpcc + NINT( zwz(1,1) )
156
157
158            ! 2. Vertical mixing for each instable portion of the density profil
159            ! ------------------------------------------------------------------
160
161            IF (zwz(1,1) /= 0.) THEN
162
163               ! -->> the density profil is statically instable :
164
165               DO ji = 1, jpi
166                  IF( zwy(ji,1) /= 0. ) THEN
167
168                     ! ikbot: ocean bottom level
169
170                     ikbot = mbathy(ji,jj)
171
172                     ! vertical iteration
173
174                     DO jiter = 1, jpk
175
176                        ! search of ikup : the first static instability from the sea surface
177
178                        ik = 0
179220                     CONTINUE
180                        ik = ik + 1
181                        IF( ik >= ikbot-1 ) GO TO 200
182                        zwx(ji,ik) = zrhop(ji,jj,ik) - zrhop(ji,jj,ik+1)
183                        IF( zwx(ji,ik) <= 0. ) GO TO 220
184                        ikup = ik
185                        ! the density profil is instable below ikup
186   
187                        ! ikdown : bottom of the instable portion of the density profil
188
189                        ! search of ikdown and vertical mixing from ikup to ikdown
190
191                        ze3tot= fse3t(ji,jj,ikup)
192                        zta   = tn   (ji,jj,ikup)
193                        zsa   = sn   (ji,jj,ikup)
194                        zraua = zrhop(ji,jj,ikup)
195
196                        DO jkdown = ikup+1, ikbot-1
197                           IF( zraua <= zrhop(ji,jj,jkdown) ) THEN
198                              ikdown = jkdown
199                              GO TO 240
200                           ENDIF
201                           ze3dwn =  fse3t(ji,jj,jkdown)
202                           ze3tot =  ze3tot + ze3dwn
203                           zta   = ( zta*(ze3tot-ze3dwn) + tn(ji,jj,jkdown)*ze3dwn )/ze3tot
204                           zsa   = ( zsa*(ze3tot-ze3dwn) + sn(ji,jj,jkdown)*ze3dwn )/ze3tot
205                           zraua = ( zraua*(ze3tot-ze3dwn) + zrhop(ji,jj,jkdown)*ze3dwn )/ze3tot
206                           inpci = inpci+1
207                        END DO
208                        ikdown = ikbot-1
209240                     CONTINUE
210
211                        DO jkp = ikup, ikdown-1
212                           tn(ji,jj,jkp) = zta
213                           sn(ji,jj,jkp) = zsa
214                           zrhop(ji,jj,jkp) = zraua
215                        END DO
216                        IF (ikdown == ikbot-1 .AND. zraua >= zrhop(ji,jj,ikdown) ) THEN
217                           tn(ji,jj,ikdown) = zta
218                           sn(ji,jj,ikdown) = zsa
219                           zrhop(ji,jj,ikdown) = zraua
220                        ENDIF
221
222                     END DO
223                  ENDIF
224200               CONTINUE
225               END DO
226
227               ! <<-- no more static instability on slab jj
228
229            ENDIF
230            !                                             ! ===============
231         END DO                                           !   End of slab
232         !                                                ! ===============
233
234
235         ! save the trends for diagnostic
236         ! Non penetrative mixing trends
237         IF( l_trdtra )   THEN
238            ztdta(:,:,:) = tn(:,:,:) - ztdta(:,:,:)
239            ztdsa(:,:,:) = sn(:,:,:) - ztdsa(:,:,:)
240
241            CALL trd_mod(ztdta, ztdsa, jpttdnpc, 'TRA', kt)
242         ENDIF
243     
244         ! Lateral boundary conditions on ( tn, sn )   ( Unchanged sign)
245         ! ------------------------------============
246         CALL lbc_lnk( tn, 'T', 1. )
247         CALL lbc_lnk( sn, 'T', 1. )
248     
249
250         !  2. non penetrative convective scheme statistics
251         !  -----------------------------------------------
252
253         IF( nnpc2 /= 0 .AND. MOD( kt, nnpc2 ) == 0 ) THEN
254            IF(lwp) WRITE(numout,*)' kt=',kt, ' number of statically instable',   &
255               ' water column : ',inpcc, ' number of iteration : ',inpci
256         ENDIF
257
258      ENDIF
259     
260   END SUBROUTINE tra_npc
261
262
263   SUBROUTINE tra_npc_init
264      !!----------------------------------------------------------------------
265      !!                  ***  ROUTINE tra_npc_init  ***
266      !!                   
267      !! ** Purpose :   initializations of the non-penetrative adjustment scheme
268      !!
269      !! History :
270      !!   8.5  !  02-12  (G. Madec)  F90 : free form
271      !!----------------------------------------------------------------------
272      !! * Namelist
273      NAMELIST/namnpc/ nnpc1, nnpc2
274      !!----------------------------------------------------------------------
275
276      ! Namelist namzdf : vertical diffusion
277      REWIND( numnam )
278      READ  ( numnam, namnpc )
279
280      ! Parameter print
281      ! ---------------
282      IF(lwp) THEN
283         WRITE(numout,*)
284         WRITE(numout,*) 'tra_npc_init : Non Penetrative Convection (npc) scheme'
285         WRITE(numout,*) '~~~~~~~~~~~~'
286         WRITE(numout,*) '          Namelist namnpc : set npc scheme parameters'
287         WRITE(numout,*)
288         WRITE(numout,*) '             npc scheme frequency           nnpc1  = ', nnpc1
289         WRITE(numout,*) '             npc scheme print frequency     nnpc2  = ', nnpc2
290         WRITE(numout,*)
291      ENDIF
292
293
294      ! Parameter controls
295      ! ------------------
296      IF ( nnpc1 == 0 ) THEN
297          IF(lwp) WRITE(numout,cform_war)
298          IF(lwp) WRITE(numout,*) '             nnpc1 = ', nnpc1, ' is forced to 1'
299          nnpc1 = 1
300          nwarn = nwarn + 1
301      ENDIF
302     
303   END SUBROUTINE tra_npc_init
304
305   !!======================================================================
306END MODULE tranpc
Note: See TracBrowser for help on using the repository browser.