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.
traatf.F90 in NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA – NEMO

source: NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traatf.F90 @ 11480

Last change on this file since 11480 was 11480, checked in by davestorkey, 5 years ago

2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps : Merge in changes from branch of branch.
Main changes:

  1. "nxt" modules renamed as "atf" and now just do Asselin time filtering. The time level swapping is achieved by swapping indices.
  2. Some additional prognostic grid variables changed to use a time dimension.

Notes:

  1. This merged branch passes SETTE tests but does not identical results to the SETTE tests with the trunk@10721 unless minor bugs to do with Euler timestepping and the OFF timestepping are fixed in the trunk (NEMO tickets #2310 and #2311).
  2. The nn_dttrc > 1 option for TOP (TOP has a different timestep to OCE) doesn't work. But it doesn't work in the trunk or NEMO 4.0 release either.
  • Property svn:keywords set to Id
File size: 19.2 KB
Line 
1MODULE traatf
2   !!======================================================================
3   !!                       ***  MODULE  traatf  ***
4   !! Ocean active tracers:  Asselin time filtering for 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   !!            4.1  !  2019-08  (A. Coward, D. Storkey) rename tranxt.F90 -> traatf.F90. Now only does time filtering.
20   !!----------------------------------------------------------------------
21
22   !!----------------------------------------------------------------------
23   !!   tra_atf       : time filtering on tracers
24   !!   tra_atf_fix   : time filtering on tracers : fixed    volume case
25   !!   tra_atf_vvl   : time filtering on tracers : variable volume case
26   !!----------------------------------------------------------------------
27   USE oce             ! ocean dynamics and tracers variables
28   USE dom_oce         ! ocean space and time domain variables
29   USE sbc_oce         ! surface boundary condition: ocean
30   USE sbcrnf          ! river runoffs
31   USE sbcisf          ! ice shelf melting
32   USE zdf_oce         ! ocean vertical mixing
33   USE domvvl          ! variable volume
34   USE trd_oce         ! trends: ocean variables
35   USE trdtra          ! trends manager: tracers
36   USE traqsr          ! penetrative solar radiation (needed for nksr)
37   USE phycst          ! physical constant
38   USE ldftra          ! lateral physics : tracers
39   USE ldfslp          ! lateral physics : slopes
40   USE bdy_oce  , ONLY : ln_bdy
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 timing          ! Timing
47#if defined key_agrif
48   USE agrif_oce_interp
49#endif
50
51   IMPLICIT NONE
52   PRIVATE
53
54   PUBLIC   tra_atf       ! routine called by step.F90
55   PUBLIC   tra_atf_fix   ! to be used in trcnxt
56   PUBLIC   tra_atf_vvl   ! to be used in trcnxt
57
58   !! * Substitutions
59#  include "vectopt_loop_substitute.h90"
60   !!----------------------------------------------------------------------
61   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
62   !! $Id$
63   !! Software governed by the CeCILL license (see ./LICENSE)
64   !!----------------------------------------------------------------------
65CONTAINS
66
67   SUBROUTINE tra_atf( kt, Kbb, Kmm, Kaa, pts )
68      !!----------------------------------------------------------------------
69      !!                   ***  ROUTINE traatf  ***
70      !!
71      !! ** Purpose :   Apply the boundary condition on the after temperature 
72      !!             and salinity fields and add the Asselin time filter on now fields.
73      !!
74      !! ** Method  :   At this stage of the computation, ta and sa are the
75      !!             after temperature and salinity as the time stepping has
76      !!             been performed in trazdf_imp or trazdf_exp module.
77      !!
78      !!              - Apply lateral boundary conditions on (ta,sa)
79      !!             at the local domain   boundaries through lbc_lnk call,
80      !!             at the one-way open boundaries (ln_bdy=T),
81      !!             at the AGRIF zoom   boundaries (lk_agrif=T)
82      !!
83      !!              - Update lateral boundary conditions on AGRIF children
84      !!             domains (lk_agrif=T)
85      !!
86      !! ** Action  : - ts(Kmm) time filtered
87      !!----------------------------------------------------------------------
88      INTEGER                                  , INTENT(in   ) :: kt             ! ocean time-step index
89      INTEGER                                  , INTENT(in   ) :: Kbb, Kmm, Kaa  ! time level indices
90      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts            ! active tracers
91      !!
92      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
93      REAL(wp) ::   zfact            ! local scalars
94      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrdt, ztrds
95      !!----------------------------------------------------------------------
96      !
97      IF( ln_timing )   CALL timing_start( 'tra_atf')
98      !
99      IF( kt == nit000 ) THEN
100         IF(lwp) WRITE(numout,*)
101         IF(lwp) WRITE(numout,*) 'tra_atf : apply Asselin time filter to "now" fields'
102         IF(lwp) WRITE(numout,*) '~~~~~~~'
103      ENDIF
104
105      ! Update after tracer on domain lateral boundaries
106      !
107#if defined key_agrif
108      CALL Agrif_tra                     ! AGRIF zoom boundaries
109#endif
110      !                                              ! local domain boundaries  (T-point, unchanged sign)
111      CALL lbc_lnk_multi( 'traatf', pts(:,:,:,jp_tem,Kaa), 'T', 1., pts(:,:,:,jp_sal,Kaa), 'T', 1. )
112      !
113      IF( ln_bdy )   CALL bdy_tra( kt, Kbb, pts, Kaa )  ! BDY open boundaries
114 
115      ! set time step size (Euler/Leapfrog)
116      IF( neuler == 0 .AND. kt == nit000 ) THEN   ;   r2dt =        rdt   ! at nit000             (Euler)
117      ELSEIF( kt <= nit000 + 1 )           THEN   ;   r2dt = 2._wp* rdt   ! at nit000 or nit000+1 (Leapfrog)
118      ENDIF
119
120      ! trends computation initialisation
121      IF( l_trdtra )   THEN                   
122         ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) )
123         ztrdt(:,:,jpk) = 0._wp
124         ztrds(:,:,jpk) = 0._wp
125         IF( ln_traldf_iso ) THEN              ! diagnose the "pure" Kz diffusive trend
126            CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_tem, jptra_zdfp, ztrdt )
127            CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_sal, jptra_zdfp, ztrds )
128         ENDIF
129         ! total trend for the non-time-filtered variables.
130         zfact = 1.0 / rdt
131         ! G Nurser 23 Mar 2017. Recalculate trend as Delta(e3t*T)/e3tn; e3tn cancel from pts(Kmm) terms
132         DO jk = 1, jpkm1
133            ztrdt(:,:,jk) = ( pts(:,:,jk,jp_tem,Kaa)*e3t(:,:,jk,Kaa) / e3t(:,:,jk,Kmm) - pts(:,:,jk,jp_tem,Kmm)) * zfact
134            ztrds(:,:,jk) = ( pts(:,:,jk,jp_sal,Kaa)*e3t(:,:,jk,Kaa) / e3t(:,:,jk,Kmm) - pts(:,:,jk,jp_sal,Kmm)) * zfact
135         END DO
136         CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_tem, jptra_tot, ztrdt )
137         CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_sal, jptra_tot, ztrds )
138         IF( ln_linssh ) THEN       ! linear sea surface height only
139            ! Store now fields before applying the Asselin filter
140            ! in order to calculate Asselin filter trend later.
141            ztrdt(:,:,:) = pts(:,:,:,jp_tem,Kmm) 
142            ztrds(:,:,:) = pts(:,:,:,jp_sal,Kmm)
143         ENDIF
144      ENDIF
145
146      IF( neuler == 0 .AND. kt == nit000 ) THEN       ! Euler time-stepping
147         !
148         IF (l_trdtra .AND. .NOT. ln_linssh ) THEN   ! Zero Asselin filter contribution must be explicitly written out since for vvl
149            !                                        ! Asselin filter is output by tra_atf_vvl that is not called on this time step
150            ztrdt(:,:,:) = 0._wp
151            ztrds(:,:,:) = 0._wp
152            CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_tem, jptra_atf, ztrdt )
153            CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_sal, jptra_atf, ztrds )
154         END IF
155         !
156      ELSE                                            ! Leap-Frog + Asselin filter time stepping
157         !
158         IF( ln_linssh ) THEN   ;   CALL tra_atf_fix( kt, Kbb, Kmm, Kaa, nit000,      'TRA', pts, jpts )  ! linear free surface
159         ELSE                   ;   CALL tra_atf_vvl( kt, Kbb, Kmm, Kaa, nit000, rdt, 'TRA', pts, sbc_tsc, sbc_tsc_b, jpts )  ! non-linear free surface
160         ENDIF
161         !
162         CALL lbc_lnk_multi( 'traatf', pts(:,:,:,jp_tem,Kbb) , 'T', 1., pts(:,:,:,jp_sal,Kbb) , 'T', 1., &
163                  &                    pts(:,:,:,jp_tem,Kmm) , 'T', 1., pts(:,:,:,jp_sal,Kmm) , 'T', 1., &
164                  &                    pts(:,:,:,jp_tem,Kaa), 'T', 1., pts(:,:,:,jp_sal,Kaa), 'T', 1.  )
165         !
166      ENDIF     
167      !
168      IF( l_trdtra .AND. ln_linssh ) THEN      ! trend of the Asselin filter (tb filtered - tb)/dt     
169         zfact = 1._wp / r2dt             
170         DO jk = 1, jpkm1
171            ztrdt(:,:,jk) = ( pts(:,:,jk,jp_tem,Kmm) - ztrdt(:,:,jk) ) * zfact
172            ztrds(:,:,jk) = ( pts(:,:,jk,jp_sal,Kmm) - ztrds(:,:,jk) ) * zfact
173         END DO
174         CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_tem, jptra_atf, ztrdt )
175         CALL trd_tra( kt, Kmm, Kaa, 'TRA', jp_sal, jptra_atf, ztrds )
176      END IF
177      IF( l_trdtra )   DEALLOCATE( ztrdt , ztrds )
178      !
179      !                        ! control print
180      IF(ln_ctl)   CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Kmm), clinfo1=' nxt  - Tn: ', mask1=tmask,   &
181         &                       tab3d_2=pts(:,:,:,jp_sal,Kmm), clinfo2=       ' Sn: ', mask2=tmask )
182      !
183      IF( ln_timing )   CALL timing_stop('tra_atf')
184      !
185   END SUBROUTINE tra_atf
186
187
188   SUBROUTINE tra_atf_fix( kt, Kbb, Kmm, Kaa, kit000, cdtype, pt, kjpt )
189      !!----------------------------------------------------------------------
190      !!                   ***  ROUTINE tra_atf_fix  ***
191      !!
192      !! ** Purpose :   fixed volume: apply the Asselin time filter to the "now" field
193      !!
194      !! ** Method  : - Apply a Asselin time filter on now fields.
195      !!
196      !! ** Action  : - pt(Kmm) ready for the next time step
197      !!----------------------------------------------------------------------
198      INTEGER                                  , INTENT(in   ) ::  kt            ! ocean time-step index
199      INTEGER                                  , INTENT(in   ) ::  Kbb, Kmm, Kaa ! time level indices
200      INTEGER                                  , INTENT(in   ) ::  kit000        ! first time step index
201      CHARACTER(len=3)                         , INTENT(in   ) ::  cdtype        ! =TRA or TRC (tracer indicator)
202      INTEGER                                  , INTENT(in   ) ::  kjpt          ! number of tracers
203      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::  pt            ! tracer fields
204      !
205      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
206      REAL(wp) ::   ztn, ztd         ! local scalars
207      !!----------------------------------------------------------------------
208      !
209      IF( kt == kit000 )  THEN
210         IF(lwp) WRITE(numout,*)
211         IF(lwp) WRITE(numout,*) 'tra_atf_fix : time filtering', cdtype
212         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
213      ENDIF
214      !
215      DO jn = 1, kjpt
216         !
217         DO jk = 1, jpkm1
218            DO jj = 2, jpjm1
219               DO ji = fs_2, fs_jpim1
220                  ztn = pt(ji,jj,jk,jn,Kmm)                                   
221                  ztd = pt(ji,jj,jk,jn,Kaa) - 2._wp * ztn + pt(ji,jj,jk,jn,Kbb)  ! time laplacian on tracers
222                  !
223                  pt(ji,jj,jk,jn,Kmm) = ztn + atfp * ztd                      ! pt <-- filtered pt
224               END DO
225           END DO
226         END DO
227         !
228      END DO
229      !
230   END SUBROUTINE tra_atf_fix
231
232
233   SUBROUTINE tra_atf_vvl( kt, Kbb, Kmm, Kaa, kit000, p2dt, cdtype, pt, psbc_tc, psbc_tc_b, kjpt )
234      !!----------------------------------------------------------------------
235      !!                   ***  ROUTINE tra_atf_vvl  ***
236      !!
237      !! ** Purpose :   Time varying volume: apply the Asselin time filter 
238      !!
239      !! ** Method  : - Apply a thickness weighted Asselin time filter on now fields.
240      !!             pt(Kmm)  = ( e3t(Kmm)*pt(Kmm) + atfp*[ e3t(Kbb)*pt(Kbb) - 2 e3t(Kmm)*pt(Kmm) + e3t_a*pt(Kaa) ] )
241      !!                       /( e3t(Kmm)         + atfp*[ e3t(Kbb)         - 2 e3t(Kmm)         + e3t(Kaa)    ] )
242      !!
243      !! ** Action  : - pt(Kmm) ready for the next time step
244      !!----------------------------------------------------------------------
245      INTEGER                                  , INTENT(in   ) ::  kt        ! ocean time-step index
246      INTEGER                                  , INTENT(in   ) ::  Kbb, Kmm, Kaa ! time level indices
247      INTEGER                                  , INTENT(in   ) ::  kit000    ! first time step index
248      REAL(wp)                                 , INTENT(in   ) ::  p2dt      ! time-step
249      CHARACTER(len=3)                         , INTENT(in   ) ::  cdtype    ! =TRA or TRC (tracer indicator)
250      INTEGER                                  , INTENT(in   ) ::  kjpt      ! number of tracers
251      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::  pt        ! tracer fields
252      REAL(wp), DIMENSION(jpi,jpj    ,kjpt)    , INTENT(in   ) ::  psbc_tc   ! surface tracer content
253      REAL(wp), DIMENSION(jpi,jpj    ,kjpt)    , INTENT(in   ) ::  psbc_tc_b ! before surface tracer content
254      !
255      LOGICAL  ::   ll_traqsr, ll_rnf, ll_isf   ! local logical
256      INTEGER  ::   ji, jj, jk, jn              ! dummy loop indices
257      REAL(wp) ::   zfact, zfact1, ztc_a , ztc_n , ztc_b , ztc_f , ztc_d    ! local scalar
258      REAL(wp) ::   zfact2, ze3t_b, ze3t_n, ze3t_a, ze3t_f, ze3t_d   !   -      -
259      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) ::   ztrd_atf
260      !!----------------------------------------------------------------------
261      !
262      IF( kt == kit000 )  THEN
263         IF(lwp) WRITE(numout,*)
264         IF(lwp) WRITE(numout,*) 'tra_atf_vvl : time filtering', cdtype
265         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
266      ENDIF
267      !
268      IF( cdtype == 'TRA' )  THEN   
269         ll_traqsr  = ln_traqsr        ! active  tracers case  and  solar penetration
270         ll_rnf     = ln_rnf           ! active  tracers case  and  river runoffs
271         ll_isf     = ln_isf           ! active  tracers case  and  ice shelf melting
272      ELSE                          ! passive tracers case
273         ll_traqsr  = .FALSE.          ! NO solar penetration
274         ll_rnf     = .FALSE.          ! NO river runoffs ????          !!gm BUG ? 
275         ll_isf     = .FALSE.          ! NO ice shelf melting/freezing  !!gm BUG ??
276      ENDIF
277      !
278      IF( ( l_trdtra .AND. cdtype == 'TRA' ) .OR. ( l_trdtrc .AND. cdtype == 'TRC' ) )   THEN
279         ALLOCATE( ztrd_atf(jpi,jpj,jpk,kjpt) )
280         ztrd_atf(:,:,:,:) = 0.0_wp
281      ENDIF
282      zfact = 1._wp / p2dt
283      zfact1 = atfp * p2dt
284      zfact2 = zfact1 * r1_rau0
285      DO jn = 1, kjpt     
286         DO jk = 1, jpkm1
287            DO jj = 2, jpjm1
288               DO ji = fs_2, fs_jpim1
289                  ze3t_b = e3t(ji,jj,jk,Kbb)
290                  ze3t_n = e3t(ji,jj,jk,Kmm)
291                  ze3t_a = e3t(ji,jj,jk,Kaa)
292                  !                                         ! tracer content at Before, now and after
293                  ztc_b  = pt(ji,jj,jk,jn,Kbb) * ze3t_b
294                  ztc_n  = pt(ji,jj,jk,jn,Kmm) * ze3t_n
295                  ztc_a  = pt(ji,jj,jk,jn,Kaa) * ze3t_a
296                  !
297                  ze3t_d = ze3t_a - 2. * ze3t_n + ze3t_b
298                  ztc_d  = ztc_a  - 2. * ztc_n  + ztc_b
299                  !
300                  ze3t_f = ze3t_n + atfp * ze3t_d
301                  ztc_f  = ztc_n  + atfp * ztc_d
302                  !
303                  IF( jk == mikt(ji,jj) ) THEN           ! first level
304                     ze3t_f = ze3t_f - zfact2 * ( (emp_b(ji,jj)    - emp(ji,jj)   )  &
305                            &                   + (fwfisf_b(ji,jj) - fwfisf(ji,jj))  )
306                     ztc_f  = ztc_f  - zfact1 * ( psbc_tc(ji,jj,jn) - psbc_tc_b(ji,jj,jn) )
307                  ENDIF
308                  IF( ln_rnf_depth ) THEN
309                     ! Rivers are not just at the surface must go down to nk_rnf(ji,jj)
310                     IF( mikt(ji,jj) <=jk .and. jk <= nk_rnf(ji,jj)  ) THEN
311                        ze3t_f = ze3t_f - zfact2 * ( - (rnf_b(ji,jj) - rnf(ji,jj)   )  ) &
312                    &                            * ( e3t(ji,jj,jk,Kmm) / h_rnf(ji,jj) ) 
313                     ENDIF
314                  ELSE
315                     IF( jk == mikt(ji,jj) ) THEN           ! first level
316                        ze3t_f = ze3t_f - zfact2 * ( - (rnf_b(ji,jj)    - rnf(ji,jj)   ) ) 
317                     ENDIF
318                  ENDIF
319
320                  !
321                  ! solar penetration (temperature only)
322                  IF( ll_traqsr .AND. jn == jp_tem .AND. jk <= nksr )                            & 
323                     &     ztc_f  = ztc_f  - zfact1 * ( qsr_hc(ji,jj,jk) - qsr_hc_b(ji,jj,jk) ) 
324                     !
325                  ! river runoff
326                  IF( ll_rnf .AND. jk <= nk_rnf(ji,jj) )                                          &
327                     &     ztc_f  = ztc_f  - zfact1 * ( rnf_tsc(ji,jj,jn) - rnf_tsc_b(ji,jj,jn) ) & 
328                     &                              * e3t(ji,jj,jk,Kmm) / h_rnf(ji,jj)
329                     !
330                  ! ice shelf
331                  IF( ll_isf ) THEN
332                     ! level fully include in the Losch_2008 ice shelf boundary layer
333                     IF ( jk >= misfkt(ji,jj) .AND. jk < misfkb(ji,jj) )                          &
334                        ztc_f  = ztc_f  - zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) )  &
335                               &                 * e3t(ji,jj,jk,Kmm) * r1_hisf_tbl (ji,jj)
336                     ! level partially include in Losch_2008 ice shelf boundary layer
337                     IF ( jk == misfkb(ji,jj) )                                                   &
338                        ztc_f  = ztc_f  - zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) )  &
339                               &                 * e3t(ji,jj,jk,Kmm) * r1_hisf_tbl (ji,jj) * ralpha(ji,jj)
340                  END IF
341                  !
342                  ze3t_f = 1.e0 / ze3t_f
343                  pt(ji,jj,jk,jn,Kmm) = ztc_f * ze3t_f    ! time filtered "now" field
344                  !
345                  IF( ( l_trdtra .and. cdtype == 'TRA' ) .OR. ( l_trdtrc .and. cdtype == 'TRC' ) ) THEN
346                     ztrd_atf(ji,jj,jk,jn) = (ztc_f - ztc_n) * zfact/ze3t_n
347                  ENDIF
348                  !
349               END DO
350            END DO
351         END DO
352         !
353      END DO
354      !
355      IF( ( l_trdtra .AND. cdtype == 'TRA' ) .OR. ( l_trdtrc .AND. cdtype == 'TRC' ) )   THEN
356         IF( l_trdtra .AND. cdtype == 'TRA' ) THEN
357            CALL trd_tra( kt, Kmm, Kaa, cdtype, jp_tem, jptra_atf, ztrd_atf(:,:,:,jp_tem) )
358            CALL trd_tra( kt, Kmm, Kaa, cdtype, jp_sal, jptra_atf, ztrd_atf(:,:,:,jp_sal) )
359         ENDIF
360         IF( l_trdtrc .AND. cdtype == 'TRC' ) THEN
361            DO jn = 1, kjpt
362               CALL trd_tra( kt, Kmm, Kaa, cdtype, jn, jptra_atf, ztrd_atf(:,:,:,jn) )
363            END DO
364         ENDIF
365         DEALLOCATE( ztrd_atf )
366      ENDIF
367      !
368   END SUBROUTINE tra_atf_vvl
369
370   !!======================================================================
371END MODULE traatf
Note: See TracBrowser for help on using the repository browser.