source: branches/UKMO/dev_r5518_GO6_package_text_diagnostics/NEMOGCM/NEMO/OPA_SRC/TRA/tranxt.F90 @ 10774

Last change on this file since 10774 was 10774, checked in by andmirek, 20 months ago

GMED 450 add flush after prints

File size: 21.9 KB
Line 
1MODULE tranxt
2   !!======================================================================
3   !!                       ***  MODULE  tranxt  ***
4   !! Ocean active tracers:  time stepping on temperature and salinity
5   !!======================================================================
6   !! History :  OPA  !  1991-11  (G. Madec)  Original code
7   !!            7.0  !  1993-03  (M. Guyon)  symetrical conditions
8   !!            8.0  !  1996-02  (G. Madec & M. Imbard)  opa release 8.0
9   !!             -   !  1996-04  (A. Weaver)  Euler forward step
10   !!            8.2  !  1999-02  (G. Madec, N. Grima)  semi-implicit pressure grad.
11   !!  NEMO      1.0  !  2002-08  (G. Madec)  F90: Free form and module
12   !!             -   !  2002-11  (C. Talandier, A-M Treguier) Open boundaries
13   !!             -   !  2005-04  (C. Deltel) Add Asselin trend in the ML budget
14   !!            2.0  !  2006-02  (L. Debreu, C. Mazauric) Agrif implementation
15   !!            3.0  !  2008-06  (G. Madec)  time stepping always done in trazdf
16   !!            3.1  !  2009-02  (G. Madec, R. Benshila)  re-introduce the vvl option
17   !!            3.3  !  2010-04  (M. Leclair, G. Madec)  semi-implicit hpg with asselin filter + modified LF-RA
18   !!             -   !  2010-05  (C. Ethe, G. Madec)  merge TRC-TRA
19   !!----------------------------------------------------------------------
20
21   !!----------------------------------------------------------------------
22   !!   tra_nxt       : time stepping on tracers
23   !!   tra_nxt_fix   : time stepping on tracers : fixed    volume case
24   !!   tra_nxt_vvl   : time stepping on tracers : variable volume case
25   !!----------------------------------------------------------------------
26   USE oce             ! ocean dynamics and tracers variables
27   USE dom_oce         ! ocean space and time domain variables
28   USE sbc_oce         ! surface boundary condition: ocean
29   USE sbcrnf          ! river runoffs
30   USE sbcisf          ! ice shelf melting/freezing
31   USE zdf_oce         ! ocean vertical mixing
32   USE domvvl          ! variable volume
33   USE dynspg_oce      ! surface     pressure gradient variables
34   USE dynhpg          ! hydrostatic pressure gradient
35   USE trd_oce         ! trends: ocean variables
36   USE trdtra          ! trends manager: tracers
37   USE traqsr          ! penetrative solar radiation (needed for nksr)
38   USE phycst          ! physical constant
39   USE ldftra_oce      ! lateral physics on tracers
40   USE bdy_oce         ! BDY open boundary condition variables
41   USE bdytra          ! open boundary condition (bdy_tra routine)
42   !
43   USE in_out_manager  ! I/O manager
44   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
45   USE prtctl          ! Print control
46   USE wrk_nemo        ! Memory allocation
47   USE timing          ! Timing
48#if defined key_agrif
49   USE agrif_opa_interp
50#endif
51
52
53   IMPLICIT NONE
54   PRIVATE
55
56   PUBLIC   tra_nxt       ! routine called by step.F90
57   PUBLIC   tra_nxt_fix   ! to be used in trcnxt
58   PUBLIC   tra_nxt_vvl   ! to be used in trcnxt
59
60   REAL(wp) ::   rbcp   ! Brown & Campana parameters for semi-implicit hpg
61
62   !! * Substitutions
63#  include "domzgr_substitute.h90"
64   !!----------------------------------------------------------------------
65   !! NEMO/OPA 3.3 , NEMO-Consortium (2010)
66   !! $Id$
67   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
68   !!----------------------------------------------------------------------
69CONTAINS
70
71   SUBROUTINE tra_nxt( kt )
72      !!----------------------------------------------------------------------
73      !!                   ***  ROUTINE tranxt  ***
74      !!
75      !! ** Purpose :   Apply the boundary condition on the after temperature 
76      !!             and salinity fields, achieved the time stepping by adding
77      !!             the Asselin filter on now fields and swapping the fields.
78      !!
79      !! ** Method  :   At this stage of the computation, ta and sa are the
80      !!             after temperature and salinity as the time stepping has
81      !!             been performed in trazdf_imp or trazdf_exp module.
82      !!
83      !!              - Apply lateral boundary conditions on (ta,sa)
84      !!             at the local domain   boundaries through lbc_lnk call,
85      !!             at the one-way open boundaries (lk_bdy=T),
86      !!             at the AGRIF zoom   boundaries (lk_agrif=T)
87      !!
88      !!              - Update lateral boundary conditions on AGRIF children
89      !!             domains (lk_agrif=T)
90      !!
91      !! ** Action  : - (tb,sb) and (tn,sn) ready for the next time step
92      !!              - (ta,sa) time averaged (t,s)   (ln_dynhpg_imp = T)
93      !!----------------------------------------------------------------------
94      INTEGER, INTENT(in) ::   kt    ! ocean time-step index
95      !!
96      INTEGER  ::   jk, jn    ! dummy loop indices
97      REAL(wp) ::   zfact     ! local scalars
98      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrdt, ztrds
99      !!----------------------------------------------------------------------
100      !
101      IF( nn_timing == 1 )  CALL timing_start( 'tra_nxt')
102      !
103      IF( kt == nit000 ) THEN
104         IF(lwp) WRITE(numout,*)
105         IF(lwp) WRITE(numout,*) 'tra_nxt : achieve the time stepping by Asselin filter and array swap'
106         IF(lwp) WRITE(numout,*) '~~~~~~~'
107         IF(lwp .AND. lflush) CALL flush(numout)
108         !
109         rbcp = 0.25_wp * (1._wp + atfp) * (1._wp + atfp) * ( 1._wp - atfp)      ! Brown & Campana parameter for semi-implicit hpg
110      ENDIF
111
112      ! Update after tracer on domain lateral boundaries
113      !
114#if defined key_agrif
115      CALL Agrif_tra                     ! AGRIF zoom boundaries
116#endif
117      !
118      CALL lbc_lnk( tsa(:,:,:,jp_tem), 'T', 1._wp )      ! local domain boundaries  (T-point, unchanged sign)
119      CALL lbc_lnk( tsa(:,:,:,jp_sal), 'T', 1._wp )
120      !
121#if defined key_bdy 
122      IF( lk_bdy )   CALL bdy_tra( kt )  ! BDY open boundaries
123#endif
124 
125      ! set time step size (Euler/Leapfrog)
126      IF( neuler == 0 .AND. kt == nit000 ) THEN   ;   r2dtra(:) =     rdttra(:)      ! at nit000             (Euler)
127      ELSEIF( kt <= nit000 + 1 )           THEN   ;   r2dtra(:) = 2._wp* rdttra(:)      ! at nit000 or nit000+1 (Leapfrog)
128      ENDIF
129
130      ! trends computation initialisation
131      IF( l_trdtra )   THEN                   
132         CALL wrk_alloc( jpi, jpj, jpk, ztrdt, ztrds )
133         ztrdt(:,:,jpk) = 0._wp
134         ztrds(:,:,jpk) = 0._wp
135         IF( ln_traldf_iso ) THEN              ! diagnose the "pure" Kz diffusive trend
136            CALL trd_tra( kt, 'TRA', jp_tem, jptra_zdfp, ztrdt )
137            CALL trd_tra( kt, 'TRA', jp_sal, jptra_zdfp, ztrds )
138         ENDIF
139         ! total trend for the non-time-filtered variables.
140         ! G Nurser 23 Mar 2017. Recalculate trend as Delta(e3t*T)/e3tn; e3tn cancel from tsn terms
141         IF( lk_vvl ) THEN
142            DO jk = 1, jpkm1
143               zfact = 1.0 / rdttra(jk)
144               ztrdt(:,:,jk) = ( tsa(:,:,jk,jp_tem)*fse3t_a(:,:,jk) / fse3t_n(:,:,jk) - tsn(:,:,jk,jp_tem)) * zfact
145               ztrds(:,:,jk) = ( tsa(:,:,jk,jp_sal)*fse3t_a(:,:,jk) / fse3t_n(:,:,jk) - tsn(:,:,jk,jp_sal)) * zfact
146            END DO
147         ELSE
148            DO jk = 1, jpkm1
149               zfact = 1.0 / rdttra(jk)
150               ztrdt(:,:,jk) = ( tsa(:,:,jk,jp_tem) - tsn(:,:,jk,jp_tem) ) * zfact 
151               ztrds(:,:,jk) = ( tsa(:,:,jk,jp_sal) - tsn(:,:,jk,jp_sal) ) * zfact 
152            END DO
153         END IF
154         CALL trd_tra( kt, 'TRA', jp_tem, jptra_tot, ztrdt )
155         CALL trd_tra( kt, 'TRA', jp_sal, jptra_tot, ztrds )
156         IF( .NOT.lk_vvl )  THEN
157            ! Store now fields before applying the Asselin filter
158            ! in order to calculate Asselin filter trend later.
159            ztrdt(:,:,:) = tsn(:,:,:,jp_tem) 
160            ztrds(:,:,:) = tsn(:,:,:,jp_sal)
161         END IF
162      ENDIF
163
164      IF( neuler == 0 .AND. kt == nit000 ) THEN       ! Euler time-stepping at first time-step (only swap)
165         DO jn = 1, jpts
166            DO jk = 1, jpkm1
167               tsn(:,:,jk,jn) = tsa(:,:,jk,jn)   
168            END DO
169         END DO
170         IF (l_trdtra.AND.lk_vvl) THEN      ! Zero Asselin filter contribution must be explicitly written out since for vvl
171                                            ! Asselin filter is output by tra_nxt_vvl that is not called on this time step
172            ztrdt(:,:,:) = 0._wp
173            ztrds(:,:,:) = 0._wp
174            CALL trd_tra( kt, 'TRA', jp_tem, jptra_atf, ztrdt )
175            CALL trd_tra( kt, 'TRA', jp_sal, jptra_atf, ztrds )
176         END IF
177      ELSE                                            ! Leap-Frog + Asselin filter time stepping
178         !
179         IF( lk_vvl )  THEN   ;   CALL tra_nxt_vvl( kt, nit000, rdttra, 'TRA', tsb, tsn, tsa,   &
180           &                                                              sbc_tsc, sbc_tsc_b, jpts )  ! variable volume level (vvl)
181         ELSE                 ;   CALL tra_nxt_fix( kt, nit000,         'TRA', tsb, tsn, tsa, jpts )  ! fixed    volume level
182         ENDIF
183      ENDIF     
184      !
185     ! trends computation
186      IF( l_trdtra.AND..NOT.lk_vvl) THEN      ! trend of the Asselin filter (tb filtered - tb)/dt     
187         DO jk = 1, jpkm1
188            zfact = 1._wp / r2dtra(jk)             
189            ztrdt(:,:,jk) = ( tsb(:,:,jk,jp_tem) - ztrdt(:,:,jk) ) * zfact
190            ztrds(:,:,jk) = ( tsb(:,:,jk,jp_sal) - ztrds(:,:,jk) ) * zfact
191         END DO
192         CALL trd_tra( kt, 'TRA', jp_tem, jptra_atf, ztrdt )
193         CALL trd_tra( kt, 'TRA', jp_sal, jptra_atf, ztrds )
194      END IF
195      IF( l_trdtra) CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds )
196      !
197      !                        ! control print
198      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsn(:,:,:,jp_tem), clinfo1=' nxt  - Tn: ', mask1=tmask,   &
199         &                       tab3d_2=tsn(:,:,:,jp_sal), clinfo2=       ' Sn: ', mask2=tmask )
200      !
201      IF( nn_timing == 1 )   CALL timing_stop('tra_nxt')
202      !
203   END SUBROUTINE tra_nxt
204
205
206   SUBROUTINE tra_nxt_fix( kt, kit000, cdtype, ptb, ptn, pta, kjpt )
207      !!----------------------------------------------------------------------
208      !!                   ***  ROUTINE tra_nxt_fix  ***
209      !!
210      !! ** Purpose :   fixed volume: apply the Asselin time filter and
211      !!                swap the tracer fields.
212      !!
213      !! ** Method  : - Apply a Asselin time filter on now fields.
214      !!              - save in (ta,sa) an average over the three time levels
215      !!             which will be used to compute rdn and thus the semi-implicit
216      !!             hydrostatic pressure gradient (ln_dynhpg_imp = T)
217      !!              - swap tracer fields to prepare the next time_step.
218      !!                This can be summurized for tempearture as:
219      !!             ztm = tn + rbcp * [ta -2 tn + tb ]       ln_dynhpg_imp = T
220      !!             ztm = 0                                   otherwise
221      !!                   with rbcp=1/4 * (1-atfp^4) / (1-atfp)
222      !!             tb  = tn + atfp*[ tb - 2 tn + ta ]
223      !!             tn  = ta 
224      !!             ta  = ztm       (NB: reset to 0 after eos_bn2 call)
225      !!
226      !! ** Action  : - (tb,sb) and (tn,sn) ready for the next time step
227      !!              - (ta,sa) time averaged (t,s)   (ln_dynhpg_imp = T)
228      !!----------------------------------------------------------------------
229      INTEGER         , INTENT(in   )                               ::   kt       ! ocean time-step index
230      INTEGER         , INTENT(in   )                               ::   kit000   ! first time step index
231      CHARACTER(len=3), INTENT(in   )                               ::   cdtype   ! =TRA or TRC (tracer indicator)
232      INTEGER         , INTENT(in   )                               ::   kjpt     ! number of tracers
233      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt)  ::   ptb      ! before tracer fields
234      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt)  ::   ptn      ! now tracer fields
235      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt)  ::   pta      ! tracer trend
236      !
237      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
238      LOGICAL  ::   ll_tra_hpg       ! local logical
239      REAL(wp) ::   ztn, ztd         ! local scalars
240      !!----------------------------------------------------------------------
241
242      IF( kt == kit000 )  THEN
243         IF(lwp) WRITE(numout,*)
244         IF(lwp) WRITE(numout,*) 'tra_nxt_fix : time stepping', cdtype
245         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
246         IF(lwp .AND. lflush) CALL flush(numout)
247      ENDIF
248      !
249      IF( cdtype == 'TRA' )  THEN   ;   ll_tra_hpg = ln_dynhpg_imp    ! active  tracers case  and  semi-implicit hpg   
250      ELSE                          ;   ll_tra_hpg = .FALSE.          ! passive tracers case or NO semi-implicit hpg
251      ENDIF
252      !
253      DO jn = 1, kjpt
254         !
255         DO jk = 1, jpkm1
256            DO jj = 1, jpj
257               DO ji = 1, jpi
258                  ztn = ptn(ji,jj,jk,jn)                                   
259                  ztd = pta(ji,jj,jk,jn) - 2. * ztn + ptb(ji,jj,jk,jn)      !  time laplacian on tracers
260                  !
261                  ptb(ji,jj,jk,jn) = ztn + atfp * ztd                       ! ptb <-- filtered ptn
262                  ptn(ji,jj,jk,jn) = pta(ji,jj,jk,jn)                       ! ptn <-- pta
263                  !
264                  IF( ll_tra_hpg )   pta(ji,jj,jk,jn) = ztn + rbcp * ztd    ! pta <-- Brown & Campana average
265               END DO
266           END DO
267         END DO
268         !
269      END DO
270      !
271   END SUBROUTINE tra_nxt_fix
272
273
274   SUBROUTINE tra_nxt_vvl( kt, kit000, p2dt, cdtype, ptb, ptn, pta, psbc_tc, psbc_tc_b, kjpt )
275      !!----------------------------------------------------------------------
276      !!                   ***  ROUTINE tra_nxt_vvl  ***
277      !!
278      !! ** Purpose :   Time varying volume: apply the Asselin time filter 
279      !!                and swap the tracer fields.
280      !!
281      !! ** Method  : - Apply a thickness weighted Asselin time filter on now fields.
282      !!              - save in (ta,sa) a thickness weighted average over the three
283      !!             time levels which will be used to compute rdn and thus the semi-
284      !!             implicit hydrostatic pressure gradient (ln_dynhpg_imp = T)
285      !!              - swap tracer fields to prepare the next time_step.
286      !!                This can be summurized for tempearture as:
287      !!             ztm = ( e3t_n*tn + rbcp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] )   ln_dynhpg_imp = T
288      !!                  /( e3t_n    + rbcp*[ e3t_b    - 2 e3t_n    + e3t_a    ] )   
289      !!             ztm = 0                                                       otherwise
290      !!             tb  = ( e3t_n*tn + atfp*[ e3t_b*tb - 2 e3t_n*tn + e3t_a*ta ] )
291      !!                  /( e3t_n    + atfp*[ e3t_b    - 2 e3t_n    + e3t_a    ] )
292      !!             tn  = ta
293      !!             ta  = zt        (NB: reset to 0 after eos_bn2 call)
294      !!
295      !! ** Action  : - (tb,sb) and (tn,sn) ready for the next time step
296      !!              - (ta,sa) time averaged (t,s)   (ln_dynhpg_imp = T)
297      !!----------------------------------------------------------------------
298      INTEGER         , INTENT(in   )                               ::  kt       ! ocean time-step index
299      INTEGER         , INTENT(in   )                               ::  kit000   ! first time step index
300      REAL(wp)        , INTENT(in   ), DIMENSION(jpk)               ::  p2dt     ! time-step
301      CHARACTER(len=3), INTENT(in   )                               ::  cdtype   ! =TRA or TRC (tracer indicator)
302      INTEGER         , INTENT(in   )                               ::  kjpt     ! number of tracers
303      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt)  ::  ptb      ! before tracer fields
304      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt)  ::  ptn      ! now tracer fields
305      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk,kjpt)  ::  pta      ! tracer trend
306      REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,kjpt)      ::  psbc_tc   ! surface tracer content
307      REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,kjpt)      ::  psbc_tc_b ! before surface tracer content
308
309      !!     
310      LOGICAL  ::   ll_tra_hpg, ll_traqsr, ll_rnf, ll_isf   ! local logical
311      INTEGER  ::   ji, jj, jk, jn              ! dummy loop indices
312      REAL(wp) ::   zfact, zfact1, ztc_a , ztc_n , ztc_b , ztc_f , ztc_d    ! local scalar
313      REAL(wp) ::   zfact2, ze3t_b, ze3t_n, ze3t_a, ze3t_f, ze3t_d   !   -      -
314      REAL(wp), POINTER, DIMENSION(:,:,:,:) :: ztrd_atf
315      !!----------------------------------------------------------------------
316      !
317      IF( kt == kit000 )  THEN
318         IF(lwp) WRITE(numout,*)
319         IF(lwp) WRITE(numout,*) 'tra_nxt_vvl : time stepping', cdtype
320         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
321         IF(lwp .AND. lflush) CALL flush(numout)
322      ENDIF
323      !
324      IF( cdtype == 'TRA' )  THEN   
325         ll_tra_hpg = ln_dynhpg_imp    ! active  tracers case  and  semi-implicit hpg
326         ll_traqsr  = ln_traqsr        ! active  tracers case  and  solar penetration
327         ll_rnf     = ln_rnf           ! active  tracers case  and  river runoffs
328         IF (nn_isf .GE. 1) THEN
329            ll_isf = .TRUE.            ! active  tracers case  and  ice shelf melting/freezing
330         ELSE
331            ll_isf = .FALSE.
332         END IF
333      ELSE                         
334         ll_tra_hpg = .FALSE.          ! passive tracers case or NO semi-implicit hpg
335         ll_traqsr  = .FALSE.          ! active  tracers case and NO solar penetration
336         ll_rnf     = .FALSE.          ! passive tracers or NO river runoffs
337         ll_isf     = .FALSE.          ! passive tracers or NO ice shelf melting/freezing
338      ENDIF
339      !
340      IF( ( l_trdtra .and. cdtype == 'TRA' ) .OR. ( l_trdtrc .and. cdtype == 'TRC' ) )   THEN
341         CALL wrk_alloc( jpi, jpj, jpk, kjpt, ztrd_atf )
342         ztrd_atf(:,:,:,:) = 0.0_wp
343      ENDIF
344      DO jn = 1, kjpt     
345         DO jk = 1, jpkm1
346            zfact  = 0.5_wp / p2dt(jk)
347            zfact1 = atfp * p2dt(jk)
348            zfact2 = zfact1 / rau0
349            DO jj = 1, jpj
350               DO ji = 1, jpi
351                  ze3t_b = fse3t_b(ji,jj,jk)
352                  ze3t_n = fse3t_n(ji,jj,jk)
353                  ze3t_a = fse3t_a(ji,jj,jk)
354                  !                                         ! tracer content at Before, now and after
355                  ztc_b  = ptb(ji,jj,jk,jn) * ze3t_b
356                  ztc_n  = ptn(ji,jj,jk,jn) * ze3t_n
357                  ztc_a  = pta(ji,jj,jk,jn) * ze3t_a
358                  !
359                  ze3t_d = ze3t_a - 2. * ze3t_n + ze3t_b
360                  ztc_d  = ztc_a  - 2. * ztc_n  + ztc_b
361                  !
362                  ze3t_f = ze3t_n + atfp * ze3t_d
363                  ztc_f  = ztc_n  + atfp * ztc_d
364                  !
365                  IF( jk == mikt(ji,jj) ) THEN           ! first level
366                     ze3t_f = ze3t_f - zfact2 * ( (emp_b(ji,jj)    - emp(ji,jj)   )  &
367                            &                   - (rnf_b(ji,jj)    - rnf(ji,jj)   )  &
368                            &                   + (fwfisf_b(ji,jj) - fwfisf(ji,jj))  )
369                     ztc_f  = ztc_f  - zfact1 * ( psbc_tc(ji,jj,jn) - psbc_tc_b(ji,jj,jn) )
370                  ENDIF
371
372                  ! solar penetration (temperature only)
373                  IF( ll_traqsr .AND. jn == jp_tem .AND. jk <= nksr )                            & 
374                     &     ztc_f  = ztc_f  - zfact1 * ( qsr_hc(ji,jj,jk) - qsr_hc_b(ji,jj,jk) ) 
375
376                  ! river runoff
377                  IF( ll_rnf .AND. jk <= nk_rnf(ji,jj) )                                          &
378                     &     ztc_f  = ztc_f  - zfact1 * ( rnf_tsc(ji,jj,jn) - rnf_tsc_b(ji,jj,jn) ) & 
379                     &                              * fse3t_n(ji,jj,jk) / h_rnf(ji,jj)
380
381                  ! ice shelf
382                  IF( ll_isf ) THEN
383                     ! level fully include in the Losch_2008 ice shelf boundary layer
384                     IF ( jk >= misfkt(ji,jj) .AND. jk < misfkb(ji,jj) )                          &
385                        ztc_f  = ztc_f  - zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) )  &
386                               &                 * fse3t_n(ji,jj,jk) * r1_hisf_tbl (ji,jj)
387                     ! level partially include in Losch_2008 ice shelf boundary layer
388                     IF ( jk == misfkb(ji,jj) )                                                   &
389                        ztc_f  = ztc_f  - zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) )  &
390                               &                 * fse3t_n(ji,jj,jk) * r1_hisf_tbl (ji,jj) * ralpha(ji,jj)
391                  END IF
392
393                  ze3t_f = 1.e0 / ze3t_f
394                  ptb(ji,jj,jk,jn) = ztc_f * ze3t_f       ! ptb <-- ptn filtered
395                  ptn(ji,jj,jk,jn) = pta(ji,jj,jk,jn)     ! ptn <-- pta
396                  !
397                  IF( ll_tra_hpg ) THEN        ! semi-implicit hpg (T & S only)
398                     ze3t_d           = 1.e0   / ( ze3t_n + rbcp * ze3t_d )
399                     pta(ji,jj,jk,jn) = ze3t_d * ( ztc_n  + rbcp * ztc_d  )   ! ta <-- Brown & Campana average
400                  ENDIF
401                  IF( ( l_trdtra .and. cdtype == 'TRA' ) .OR. ( l_trdtrc .and. cdtype == 'TRC' ) ) THEN
402                     ztrd_atf(ji,jj,jk,jn) = (ztc_f - ztc_n) * zfact/ze3t_n
403                  ENDIF
404               END DO
405            END DO
406         END DO
407         !
408      END DO
409      !
410      IF( l_trdtra .and. cdtype == 'TRA' ) THEN
411         CALL trd_tra( kt, cdtype, jp_tem, jptra_atf, ztrd_atf(:,:,:,jp_tem) )
412         CALL trd_tra( kt, cdtype, jp_sal, jptra_atf, ztrd_atf(:,:,:,jp_sal) )
413         CALL wrk_dealloc( jpi, jpj, jpk, kjpt, ztrd_atf )
414      ENDIF
415      IF( l_trdtrc .and. cdtype == 'TRC' ) THEN
416         DO jn = 1, kjpt
417            CALL trd_tra( kt, cdtype, jn, jptra_atf, ztrd_atf(:,:,:,jn) )
418         END DO
419         CALL wrk_dealloc( jpi, jpj, jpk, kjpt, ztrd_atf )
420      ENDIF
421
422   END SUBROUTINE tra_nxt_vvl
423
424   !!======================================================================
425END MODULE tranxt
Note: See TracBrowser for help on using the repository browser.