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.
traadv_ubs.F90 in branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/OPA_SRC/TRA – NEMO

source: branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_ubs.F90 @ 5870

Last change on this file since 5870 was 5870, checked in by acc, 8 years ago

Branch 2015/dev_r5803_NOC_WAD. Merge in trunk changes from 5803 to 5869 in preparation for merge. Also tidied and reorganised some wetting and drying code. Renamed wadlmt.F90 to wetdry.F90. Wetting drying code changes restricted to domzgr.F90, domvvl.F90 nemogcm.F90 sshwzv.F90, dynspg_ts.F90, wetdry.F90 and dynhpg.F90. Code passes full SETTE tests with ln_wd=.false.. Still awaiting test case for checking with ln_wd=.false.

  • Property svn:keywords set to Id
File size: 20.4 KB
Line 
1MODULE traadv_ubs
2   !!==============================================================================
3   !!                       ***  MODULE  traadv_ubs  ***
4   !! Ocean active tracers:  horizontal & vertical advective trend
5   !!==============================================================================
6   !! History :  1.0  !  2006-08  (L. Debreu, R. Benshila)  Original code
7   !!            3.3  !  2010-05  (C. Ethe, G. Madec)  merge TRC-TRA + switch from velocity to transport
8   !!----------------------------------------------------------------------
9
10   !!----------------------------------------------------------------------
11   !!   tra_adv_ubs : update the tracer trend with the horizontal
12   !!                 advection trends using a third order biaised scheme 
13   !!----------------------------------------------------------------------
14   USE oce            ! ocean dynamics and active tracers
15   USE dom_oce        ! ocean space and time domain
16   USE trc_oce        ! share passive tracers/Ocean variables
17   USE trd_oce        ! trends: ocean variables
18   USE traadv_fct      ! acces to routine interp_4th_cpt
19   USE trdtra         ! trends manager: tracers
20   USE dynspg_oce     ! choice/control of key cpp for surface pressure gradient
21   USE diaptr         ! poleward transport diagnostics
22   !
23   USE lib_mpp        ! I/O library
24   USE lbclnk         ! ocean lateral boundary condition (or mpp link)
25   USE in_out_manager ! I/O manager
26   USE wrk_nemo       ! Memory Allocation
27   USE timing         ! Timing
28   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
29
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC   tra_adv_ubs   ! routine called by traadv module
34
35   LOGICAL :: l_trd  ! flag to compute trends or not
36
37   !! * Substitutions
38#  include "domzgr_substitute.h90"
39#  include "vectopt_loop_substitute.h90"
40   !!----------------------------------------------------------------------
41   !! NEMO/OPA 3.7 , NEMO Consortium (2015)
42   !! $Id$
43   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
44   !!----------------------------------------------------------------------
45CONTAINS
46
47   SUBROUTINE tra_adv_ubs( kt, kit000, cdtype, p2dt, pun, pvn, pwn,          &
48      &                                                ptb, ptn, pta, kjpt, kn_ubs_v )
49      !!----------------------------------------------------------------------
50      !!                  ***  ROUTINE tra_adv_ubs  ***
51      !!                 
52      !! ** Purpose :   Compute the now trend due to the advection of tracers
53      !!      and add it to the general trend of passive tracer equations.
54      !!
55      !! ** Method  :   The 3rd order Upstream Biased Scheme (UBS) is based on an
56      !!      upstream-biased parabolic interpolation (Shchepetkin and McWilliams 2005)
57      !!      It is only used in the horizontal direction.
58      !!      For example the i-component of the advective fluxes are given by :
59      !!                !  e2u e3u un ( mi(Tn) - zltu(i  ) )   if un(i) >= 0
60      !!          ztu = !  or
61      !!                !  e2u e3u un ( mi(Tn) - zltu(i+1) )   if un(i) < 0
62      !!      where zltu is the second derivative of the before temperature field:
63      !!          zltu = 1/e3t di[ e2u e3u / e1u di[Tb] ]
64      !!        This results in a dissipatively dominant (i.e. hyper-diffusive)
65      !!      truncation error. The overall performance of the advection scheme
66      !!      is similar to that reported in (Farrow and Stevens, 1995).
67      !!        For stability reasons, the first term of the fluxes which corresponds
68      !!      to a second order centered scheme is evaluated using the now velocity
69      !!      (centered in time) while the second term which is the diffusive part
70      !!      of the scheme, is evaluated using the before velocity (forward in time).
71      !!      Note that UBS is not positive. Do not use it on passive tracers.
72      !!                On the vertical, the advection is evaluated using a FCT scheme,
73      !!      as the UBS have been found to be too diffusive.
74!!gm  !!                kn_ubs_v argument (not coded for the moment)
75      !!      controles whether the FCT is based on a 2nd order centrered scheme (kn_ubs_v=2)
76      !!      or on a 4th order compact scheme (kn_ubs_v=4).
77      !!
78      !! ** Action : - update (pta) with the now advective tracer trends
79      !!
80      !! Reference : Shchepetkin, A. F., J. C. McWilliams, 2005, Ocean Modelling, 9, 347-404.
81      !!             Farrow, D.E., Stevens, D.P., 1995, J. Phys. Ocean. 25, 1731Ð1741.
82      !!----------------------------------------------------------------------
83      INTEGER                              , INTENT(in   ) ::   kt              ! ocean time-step index
84      INTEGER                              , INTENT(in   ) ::   kit000          ! first time step index
85      CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator)
86      INTEGER                              , INTENT(in   ) ::   kjpt            ! number of tracers
87      INTEGER                              , INTENT(in   ) ::   kn_ubs_v        ! number of tracers
88      REAL(wp), DIMENSION(        jpk     ), INTENT(in   ) ::   p2dt            ! vertical profile of tracer time-step
89      REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pun, pvn, pwn   ! 3 ocean transport components
90      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb, ptn        ! before and now tracer fields
91      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta             ! tracer trend
92      !
93      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
94      REAL(wp) ::   ztra, zbtr, zcoef, z2dtt                       ! local scalars
95      REAL(wp) ::   zfp_ui, zfm_ui, zcenut, ztak, zfp_wk, zfm_wk   !   -      -
96      REAL(wp) ::   zfp_vj, zfm_vj, zcenvt, zeeu, zeev, z_hdivn    !   -      -
97      REAL(wp), POINTER, DIMENSION(:,:,:) :: ztu, ztv, zltu, zltv, zti, ztw
98      !!----------------------------------------------------------------------
99      !
100      IF( nn_timing == 1 )  CALL timing_start('tra_adv_ubs')
101      !
102      CALL wrk_alloc( jpi,jpj,jpk,   ztu, ztv, zltu, zltv, zti, ztw )
103      !
104      IF( kt == kit000 )  THEN
105         IF(lwp) WRITE(numout,*)
106         IF(lwp) WRITE(numout,*) 'tra_adv_ubs :  horizontal UBS advection scheme on ', cdtype
107         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~'
108      ENDIF
109      !
110      l_trd = .FALSE.
111      IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.
112      !
113      zltu(:,:,jpk) = 0._wp   ;   zltv(:,:,jpk) = 0._wp     ! Bottom value : set to zero one for all
114      ztw (:,:,jpk) = 0._wp   ;   zti (:,:,jpk) = 0._wp
115      IF( lk_vvl )   ztw(:,:, 1 ) = 0._wp                   ! surface value: set to zero only in vvl case
116      !
117      !                                                          ! ===========
118      DO jn = 1, kjpt                                            ! tracer loop
119         !                                                       ! ===========
120         !                                             
121         DO jk = 1, jpkm1        !==  horizontal laplacian of before tracer ==!
122            DO jj = 1, jpjm1              ! First derivative (masked gradient)
123               DO ji = 1, fs_jpim1   ! vector opt.
124                  zeeu = e2_e1u(ji,jj) * fse3u(ji,jj,jk) * umask(ji,jj,jk)
125                  zeev = e1_e2v(ji,jj) * fse3v(ji,jj,jk) * vmask(ji,jj,jk)
126                  ztu(ji,jj,jk) = zeeu * ( ptb(ji+1,jj  ,jk,jn) - ptb(ji,jj,jk,jn) )
127                  ztv(ji,jj,jk) = zeev * ( ptb(ji  ,jj+1,jk,jn) - ptb(ji,jj,jk,jn) )
128               END DO
129            END DO
130            DO jj = 2, jpjm1              ! Second derivative (divergence)
131               DO ji = fs_2, fs_jpim1   ! vector opt.
132                  zcoef = 1._wp / ( 6._wp * fse3t(ji,jj,jk) )
133                  zltu(ji,jj,jk) = (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk)  ) * zcoef
134                  zltv(ji,jj,jk) = (  ztv(ji,jj,jk) - ztv(ji,jj-1,jk)  ) * zcoef
135               END DO
136            END DO
137            !                                   
138         END DO         
139         CALL lbc_lnk( zltu, 'T', 1. )   ;    CALL lbc_lnk( zltv, 'T', 1. )   ! Lateral boundary cond. (unchanged sgn)
140         !   
141         DO jk = 1, jpkm1        !==  Horizontal advective fluxes  ==!     (UBS)
142            DO jj = 1, jpjm1
143               DO ji = 1, fs_jpim1   ! vector opt.
144                  zfp_ui = pun(ji,jj,jk) + ABS( pun(ji,jj,jk) )      ! upstream transport (x2)
145                  zfm_ui = pun(ji,jj,jk) - ABS( pun(ji,jj,jk) )
146                  zfp_vj = pvn(ji,jj,jk) + ABS( pvn(ji,jj,jk) )
147                  zfm_vj = pvn(ji,jj,jk) - ABS( pvn(ji,jj,jk) )
148                  !                                                  ! 2nd order centered advective fluxes (x2)
149                  zcenut = pun(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji+1,jj  ,jk,jn) )
150                  zcenvt = pvn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji  ,jj+1,jk,jn) )
151                  !                                                  ! UBS advective fluxes
152                  ztu(ji,jj,jk) = 0.5 * ( zcenut - zfp_ui * zltu(ji,jj,jk) - zfm_ui * zltu(ji+1,jj,jk) )
153                  ztv(ji,jj,jk) = 0.5 * ( zcenvt - zfp_vj * zltv(ji,jj,jk) - zfm_vj * zltv(ji,jj+1,jk) )
154               END DO
155            END DO
156         END DO         
157         !
158         zltu(:,:,:) = pta(:,:,:,jn)      ! store the initial trends before its update
159         !
160         DO jk = 1, jpkm1        !==  add the horizontal advective trend  ==!
161            DO jj = 2, jpjm1
162               DO ji = fs_2, fs_jpim1   ! vector opt.
163                  pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn)                        &
164                     &             - (  ztu(ji,jj,jk) - ztu(ji-1,jj  ,jk)    &
165                     &                + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk)  ) / ( e1e2t(ji,jj) * fse3t(ji,jj,jk) )
166               END DO
167            END DO
168            !                                             
169         END DO
170         !
171         zltu(:,:,:) = pta(:,:,:,jn) - zltu(:,:,:)    ! Horizontal advective trend used in vertical 2nd order FCT case
172         !                                            ! and/or in trend diagnostic (l_trd=T)
173         !               
174         IF( l_trd ) THEN                  ! trend diagnostics
175             CALL trd_tra( kt, cdtype, jn, jptra_xad, ztu, pun, ptn(:,:,:,jn) )
176             CALL trd_tra( kt, cdtype, jn, jptra_yad, ztv, pvn, ptn(:,:,:,jn) )
177         END IF
178         !                                 ! "Poleward" heat and salt transports (contribution of upstream fluxes)
179         IF( cdtype == 'TRA' .AND. ln_diaptr ) THEN 
180            IF( jn == jp_tem )  htr_adv(:) = ptr_sj( ztv(:,:,:) )
181            IF( jn == jp_sal )  str_adv(:) = ptr_sj( ztv(:,:,:) )
182         ENDIF
183         !
184         !                       !== vertical advective trend  ==!
185         !
186         SELECT CASE( kn_ubs_v )       ! select the vertical advection scheme
187         !
188         CASE(  2  )                   ! 2nd order FCT
189            !         
190            IF( l_trd )   zltv(:,:,:) = pta(:,:,:,jn)          ! store pta if trend diag.
191            !
192            !                          !*  upstream advection with initial mass fluxes & intermediate update  ==!
193            DO jk = 2, jpkm1                 ! Interior value (w-masked)
194               DO jj = 1, jpj
195                  DO ji = 1, jpi
196                     zfp_wk = pwn(ji,jj,jk) + ABS( pwn(ji,jj,jk) )
197                     zfm_wk = pwn(ji,jj,jk) - ABS( pwn(ji,jj,jk) )
198                     ztw(ji,jj,jk) = 0.5_wp * (  zfp_wk * ptb(ji,jj,jk,jn) + zfm_wk * ptb(ji,jj,jk-1,jn)  ) * wmask(ji,jj,jk)
199                  END DO
200               END DO
201            END DO
202            IF(.NOT.lk_vvl ) THEN            ! top ocean value (only in linear free surface as ztw has been w-masked)
203               IF( ln_isfcav ) THEN                ! top of the ice-shelf cavities and at the ocean surface
204                  DO jj = 1, jpj
205                     DO ji = 1, jpi
206                        ztw(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptb(ji,jj,mikt(ji,jj),jn)   ! linear free surface
207                     END DO
208                  END DO   
209               ELSE                                ! no cavities: only at the ocean surface
210                  ztw(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn)
211               ENDIF
212            ENDIF
213            !
214            DO jk = 1, jpkm1           !* trend and after field with monotonic scheme
215               z2dtt = p2dt(jk)
216               DO jj = 2, jpjm1
217                  DO ji = fs_2, fs_jpim1   ! vector opt.
218                     ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) / ( e1e2t(ji,jj) * fse3t_n(ji,jj,jk) )
219                     pta(ji,jj,jk,jn) =   pta(ji,jj,jk,jn) +  ztak 
220                     zti(ji,jj,jk)    = ( ptb(ji,jj,jk,jn) + z2dtt * ( ztak + zltu(ji,jj,jk) ) ) * tmask(ji,jj,jk)
221                  END DO
222               END DO
223            END DO
224            CALL lbc_lnk( zti, 'T', 1. )      ! Lateral boundary conditions on zti, zsi   (unchanged sign)
225            !
226            !                          !*  anti-diffusive flux : high order minus low order
227            DO jk = 2, jpkm1        ! Interior value  (w-masked)
228               DO jj = 1, jpj
229                  DO ji = 1, jpi
230                     ztw(ji,jj,jk) = (   0.5_wp * pwn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) )   &
231                        &              - ztw(ji,jj,jk)   ) * wmask(ji,jj,jk)
232                  END DO
233               END DO
234            END DO
235            !                                            ! top ocean value: high order == upstream  ==>>  zwz=0
236            IF(.NOT.lk_vvl )   ztw(:,:, 1 ) = 0._wp      ! only ocean surface as interior zwz values have been w-masked
237            !
238            CALL nonosc_z( ptb(:,:,:,jn), ztw, zti, p2dt )      !  monotonicity algorithm
239            !
240         CASE(  4  )                               ! 4th order COMPACT
241            CALL interp_4th_cpt( ptn(:,:,:,jn) , ztw )         ! 4th order compact interpolation of T at w-point
242            DO jk = 2, jpkm1
243               DO jj = 2, jpjm1
244                  DO ji = fs_2, fs_jpim1
245                     ztw(ji,jj,jk) = pwn(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk)
246                  END DO
247               END DO
248            END DO
249            IF(.NOT.lk_vvl )   ztw(:,:, 1 ) = pwn(:,:,1) * ptn(:,:,1,jn)     !!gm ISF & 4th COMPACT doesn't work
250            !
251         END SELECT
252         !
253         DO jk = 1, jpkm1        !  final trend with corrected fluxes
254            DO jj = 2, jpjm1 
255               DO ji = fs_2, fs_jpim1   ! vector opt.   
256                  pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) / ( e1e2t(ji,jj) * fse3t_n(ji,jj,jk) )
257               END DO
258            END DO
259         END DO
260         !
261         IF( l_trd )  THEN       ! vertical advective trend diagnostics
262            DO jk = 1, jpkm1                       ! (compute -w.dk[ptn]= -dk[w.ptn] + ptn.dk[w])
263               DO jj = 2, jpjm1
264                  DO ji = fs_2, fs_jpim1   ! vector opt.
265                     zltv(ji,jj,jk) = pta(ji,jj,jk,jn) - zltv(ji,jj,jk)                          &
266                        &           + ptn(ji,jj,jk,jn) * (  pwn(ji,jj,jk) - pwn(ji,jj,jk+1)  )   &
267                        &                              / ( e1e2t(ji,jj) * fse3t_n(ji,jj,jk) )
268                  END DO
269               END DO
270            END DO
271            CALL trd_tra( kt, cdtype, jn, jptra_zad, zltv )
272         ENDIF
273         !
274      END DO
275      !
276      CALL wrk_dealloc( jpi,jpj,jpk,   ztu, ztv, zltu, zltv, zti, ztw )
277      !
278      IF( nn_timing == 1 )  CALL timing_stop('tra_adv_ubs')
279      !
280   END SUBROUTINE tra_adv_ubs
281
282
283   SUBROUTINE nonosc_z( pbef, pcc, paft, p2dt )
284      !!---------------------------------------------------------------------
285      !!                    ***  ROUTINE nonosc_z  ***
286      !!     
287      !! **  Purpose :   compute monotonic tracer fluxes from the upstream
288      !!       scheme and the before field by a nonoscillatory algorithm
289      !!
290      !! **  Method  :   ... ???
291      !!       warning : pbef and paft must be masked, but the boundaries
292      !!       conditions on the fluxes are not necessary zalezak (1979)
293      !!       drange (1995) multi-dimensional forward-in-time and upstream-
294      !!       in-space based differencing for fluid
295      !!----------------------------------------------------------------------
296      REAL(wp), INTENT(in   ), DIMENSION(jpk)          ::   p2dt   ! vertical profile of tracer time-step
297      REAL(wp),                DIMENSION (jpi,jpj,jpk) ::   pbef   ! before field
298      REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) ::   paft   ! after field
299      REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) ::   pcc    ! monotonic flux in the k direction
300      !
301      INTEGER  ::   ji, jj, jk   ! dummy loop indices
302      INTEGER  ::   ikm1         ! local integer
303      REAL(wp) ::   zpos, zneg, zbt, za, zb, zc, zbig, zrtrn, z2dtt   ! local scalars
304      REAL(wp), POINTER, DIMENSION(:,:,:) :: zbetup, zbetdo
305      !!----------------------------------------------------------------------
306      !
307      IF( nn_timing == 1 )  CALL timing_start('nonosc_z')
308      !
309      CALL wrk_alloc( jpi,jpj,jpk,   zbetup, zbetdo )
310      !
311      zbig  = 1.e+40_wp
312      zrtrn = 1.e-15_wp
313      zbetup(:,:,:) = 0._wp   ;   zbetdo(:,:,:) = 0._wp
314      !
315      ! Search local extrema
316      ! --------------------
317      !                    ! large negative value (-zbig) inside land
318      pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) )
319      paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) )
320      !
321      DO jk = 1, jpkm1     ! search maximum in neighbourhood
322         ikm1 = MAX(jk-1,1)
323         DO jj = 2, jpjm1
324            DO ji = fs_2, fs_jpim1   ! vector opt.
325               zbetup(ji,jj,jk) = MAX(  pbef(ji  ,jj  ,jk  ), paft(ji  ,jj  ,jk  ),   &
326                  &                     pbef(ji  ,jj  ,ikm1), pbef(ji  ,jj  ,jk+1),   &
327                  &                     paft(ji  ,jj  ,ikm1), paft(ji  ,jj  ,jk+1)  )
328            END DO
329         END DO
330      END DO
331      !                    ! large positive value (+zbig) inside land
332      pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) )
333      paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) )
334      !
335      DO jk = 1, jpkm1     ! search minimum in neighbourhood
336         ikm1 = MAX(jk-1,1)
337         DO jj = 2, jpjm1
338            DO ji = fs_2, fs_jpim1   ! vector opt.
339               zbetdo(ji,jj,jk) = MIN(  pbef(ji  ,jj  ,jk  ), paft(ji  ,jj  ,jk  ),   &
340                  &                     pbef(ji  ,jj  ,ikm1), pbef(ji  ,jj  ,jk+1),   &
341                  &                     paft(ji  ,jj  ,ikm1), paft(ji  ,jj  ,jk+1)  )
342            END DO
343         END DO
344      END DO
345      !                    ! restore masked values to zero
346      pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:)
347      paft(:,:,:) = paft(:,:,:) * tmask(:,:,:)
348      !
349      ! Positive and negative part of fluxes and beta terms
350      ! ---------------------------------------------------
351      DO jk = 1, jpkm1
352         z2dtt = p2dt(jk)
353         DO jj = 2, jpjm1
354            DO ji = fs_2, fs_jpim1   ! vector opt.
355               ! positive & negative part of the flux
356               zpos = MAX( 0., pcc(ji  ,jj  ,jk+1) ) - MIN( 0., pcc(ji  ,jj  ,jk  ) )
357               zneg = MAX( 0., pcc(ji  ,jj  ,jk  ) ) - MIN( 0., pcc(ji  ,jj  ,jk+1) )
358               ! up & down beta terms
359               zbt = e1e2t(ji,jj) * fse3t(ji,jj,jk) / z2dtt
360               zbetup(ji,jj,jk) = ( zbetup(ji,jj,jk) - paft(ji,jj,jk) ) / (zpos+zrtrn) * zbt
361               zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zbetdo(ji,jj,jk) ) / (zneg+zrtrn) * zbt
362            END DO
363         END DO
364      END DO
365      !
366      ! monotonic flux in the k direction, i.e. pcc
367      ! -------------------------------------------
368      DO jk = 2, jpkm1
369         DO jj = 2, jpjm1
370            DO ji = fs_2, fs_jpim1   ! vector opt.
371               za = MIN( 1., zbetdo(ji,jj,jk), zbetup(ji,jj,jk-1) )
372               zb = MIN( 1., zbetup(ji,jj,jk), zbetdo(ji,jj,jk-1) )
373               zc = 0.5 * ( 1.e0 + SIGN( 1.e0, pcc(ji,jj,jk) ) )
374               pcc(ji,jj,jk) = pcc(ji,jj,jk) * ( zc * za + ( 1.e0 - zc) * zb )
375            END DO
376         END DO
377      END DO
378      !
379      CALL wrk_dealloc( jpi,jpj,jpk,   zbetup, zbetdo )
380      !
381      IF( nn_timing == 1 )  CALL timing_stop('nonosc_z')
382      !
383   END SUBROUTINE nonosc_z
384
385   !!======================================================================
386END MODULE traadv_ubs
Note: See TracBrowser for help on using the repository browser.