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.
dtadyn.F90 in trunk/NEMO/OFF_SRC – NEMO

source: trunk/NEMO/OFF_SRC/dtadyn.F90 @ 1735

Last change on this file since 1735 was 1735, checked in by cetlod, 14 years ago

suppress useless variables in TOP modules, see ticket:602

  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 40.6 KB
Line 
1MODULE dtadyn
2   !!======================================================================
3   !!                       ***  MODULE  dtadyn  ***
4   !! OFFLINE : interpolation of the physical fields
5   !!=====================================================================
6
7   !!----------------------------------------------------------------------
8   !!   dta_dyn_init : initialization, namelist read, and parameters control
9   !!   dta_dyn      : Interpolation of the fields
10   !!----------------------------------------------------------------------
11   !! * Modules used
12   USE oce             ! ocean dynamics and tracers variables
13   USE dom_oce         ! ocean space and time domain variables
14   USE zdf_oce         ! ocean vertical physics
15   USE in_out_manager  ! I/O manager
16   USE phycst          ! physical constants
17   USE sbc_oce
18   USE ldfslp
19   USE ldfeiv          ! eddy induced velocity coef.      (ldf_eiv routine)
20   USE ldftra_oce      ! ocean tracer   lateral physics
21   USE zdfmxl
22   USE trabbl
23   USE eosbn2
24   USE zdfddm          ! vertical  physics: double diffusion
25   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
26   USE zpshde
27   USE lib_mpp         ! distributed memory computing library
28
29   IMPLICIT NONE
30   PRIVATE
31
32   !! *  Routine accessibility
33   PUBLIC dta_dyn_init   ! called by opa.F90
34   PUBLIC dta_dyn        ! called by step.F90
35
36   LOGICAL , PUBLIC :: &
37      lperdyn = .TRUE. , & ! boolean for periodic fields or not
38      lfirdyn = .TRUE.     ! boolean for the first call or not
39
40   INTEGER , PUBLIC :: &
41      ndtadyn = 73 ,  & ! Number of dat in one year
42      ndtatot = 73 ,  & ! Number of data in the input field
43      nsptint = 1 ,   & ! type of spatial interpolation
44      nficdyn = 2       ! number of dynamical fields
45
46   CHARACTER(len=45)  ::  &
47      cfile_grid_T = 'dyna_grid_T.nc', &   !: name of the grid_T file
48      cfile_grid_U = 'dyna_grid_U.nc', &   !: name of the grid_U file
49      cfile_grid_V = 'dyna_grid_V.nc', &   !: name of the grid_V file
50      cfile_grid_W = 'dyna_grid_W.nc'      !: name of the grid_W file
51   
52   REAL(wp)      ::   &
53      rnspdta  ,       &  !: number of time step per 2 consecutives data
54      rnspdta2            !: rnspdta * 0.5
55
56   INTEGER ::     &
57      ndyn1, ndyn2 , &
58      nlecoff = 0  , & ! switch for the first read
59      numfl_t, numfl_u, &
60      numfl_v, numfl_w
61
62   REAL(wp), DIMENSION(jpi,jpj,jpk,2) ::   &
63      tdta   ,   & ! temperature at two consecutive times
64      sdta   ,   & ! salinity at two consecutive times
65      udta   ,   & ! zonal velocity at two consecutive times
66      vdta   ,   & ! meridional velocity at two consecutive times
67      wdta   ,   & ! vertical velocity at two consecutive times
68#if defined key_trc_diatrd
69      hdivdta,   & ! horizontal divergence
70#endif
71      avtdta       ! vertical diffusivity coefficient
72
73   REAL(wp), DIMENSION(jpi,jpj,2) ::       &
74      hmlddta,   & ! mixed layer depth at two consecutive times
75      wspddta,   & ! wind speed at two consecutive times
76      frlddta,   & ! sea-ice fraction at two consecutive times
77      empdta ,   & ! E-P at two consecutive times
78      qsrdta       ! short wave heat flux at two consecutive times
79
80#if defined key_ldfslp
81   REAL(wp), DIMENSION(jpi,jpj,jpk,2) ::   &
82      uslpdta ,  & ! zonal isopycnal slopes
83      vslpdta ,  & ! meridional isopycnal slopes
84      wslpidta , & ! zonal diapycnal slopes
85      wslpjdta     ! meridional diapycnal slopes
86#endif
87
88#if ! defined key_off_degrad &&  defined key_traldf_c2d
89   REAL(wp), DIMENSION(jpi,jpj,2) ::   &
90      ahtwdta      ! Lateral diffusivity
91# if defined key_trcldf_eiv 
92   REAL(wp), DIMENSION(jpi,jpj,2) ::   &
93      aeiwdta      ! G&M coefficient
94# endif
95#endif
96
97#if defined key_off_degrad
98   REAL(wp), DIMENSION(jpi,jpj,jpk,2) ::   &
99      ahtudta, ahtvdta, ahtwdta  !  Lateral diffusivity
100# if defined key_trcldf_eiv
101   REAL(wp), DIMENSION(jpi,jpj,jpk,2) ::   &
102      aeiudta, aeivdta, aeiwdta  ! G&M coefficient
103# endif
104
105#endif
106
107#if defined key_trcbbl_dif   ||   defined key_trcbbl_adv
108   REAL(wp), DIMENSION(jpi,jpj,2) ::       &
109      bblxdta ,  & ! frequency of bbl in the x direction at 2 consecutive times
110      bblydta      ! frequency of bbl in the y direction at 2 consecutive times
111#endif
112
113   !! * Substitutions
114#  include "domzgr_substitute.h90"
115#  include "vectopt_loop_substitute.h90"
116   !!----------------------------------------------------------------------
117   !!   OPA 9.0 , LOCEAN-IPSL  (2005)
118   !!   $Id$
119   !!   This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
120   !!----------------------------------------------------------------------
121
122CONTAINS
123
124   SUBROUTINE dta_dyn( kt )
125      !!----------------------------------------------------------------------
126      !!                  ***  ROUTINE dta_dyn  ***
127      !!
128      !! ** Purpose : Prepares dynamics and physics fields from an
129      !!              OPA9 simulation  for an off-line simulation
130      !!               for passive tracer
131      !!
132      !! ** Method : calculates the position of DATA to read READ DATA
133      !!             (example month changement) computes slopes IF needed
134      !!             interpolates DATA IF needed
135      !!
136      !! ** History :
137      !!   ! original  : 92-01 (M. Imbard: sub domain)
138      !!   ! addition  : 98-04 (L.Bopp MA Foujols: slopes for isopyc.)
139      !!   ! addition  : 98-05 (L. Bopp read output of coupled run)
140      !!   ! addition  : 05-03 (O. Aumont and A. El Moussaoui) F90
141      !!   ! addition  : 05-12 (C. Ethe) Adapted for DEGINT
142      !!----------------------------------------------------------------------
143      !! * Arguments
144      INTEGER, INTENT( in ) ::   kt       ! ocean time-step index
145
146      !! * Local declarations
147      INTEGER ::   iper, iperm1, iswap, izt   
148
149      REAL(wp) :: zpdtan, zpdtpe, zdemi, zt
150      REAL(wp) :: zweigh
151
152      ! 0. Initialization
153      ! -----------------
154
155      IF( lfirdyn ) THEN
156         ! first time step MUST BE nit000
157         IF( kt /= nit000 ) THEN
158            IF (lwp) THEN
159               WRITE (numout,*) ' kt MUST BE EQUAL to nit000. kt = ',kt ,' nit000 = ',nit000 
160              STOP 'dtadyn'
161            ENDIF
162          ENDIF 
163          ! Initialize the parameters of the interpolation
164          CALL dta_dyn_init
165      ENDIF
166
167      zt       = ( FLOAT (kt) + rnspdta2 ) / rnspdta
168      izt      = INT( zt )
169      zweigh   = zt - FLOAT( INT(zt) )
170
171      IF( lperdyn ) THEN
172         iperm1 = MOD( izt, ndtadyn )
173      ELSE
174         iperm1 = MOD( izt, ndtatot - 1 ) + 1
175      ENDIF
176
177      iper = iperm1 + 1
178      IF( iperm1 == 0 ) THEN
179          IF( lperdyn ) THEN
180              iperm1 = ndtadyn
181          ELSE
182              IF( lfirdyn ) THEN
183                  IF (lwp) WRITE (numout,*) & 
184                      &   ' dynamic file is not periodic with or without interpolation  &
185                      &   we take the first value for the previous period iperm1 = 0  '
186              END IF
187          END IF
188      END IF
189
190      iswap  = 0
191
192      ! 1. First call lfirdyn = true
193      ! ----------------------------
194
195      IF( lfirdyn ) THEN
196         ! store the information of the period read
197         ndyn1 = iperm1
198         ndyn2 = iper
199         
200         IF (lwp) THEN
201            WRITE (numout,*) ' dynamics data read for the period ndyn1 =',ndyn1, &
202               &             ' and for the period ndyn2 = ',ndyn2
203            WRITE (numout,*) ' time step is : ', kt
204            WRITE (numout,*) ' we have ndtadyn = ',ndtadyn,' records in the dynamic file for one year'
205         END IF
206         !
207         IF( iperm1 /= 0 ) THEN         ! data read for the iperm1 period
208            CALL dynrea( kt, iperm1 ) 
209         ELSE
210            CALL dynrea( kt, 1 )
211         ENDIF
212         
213#if defined key_ldfslp
214         ! Computes slopes
215         ! Caution : here tn, sn and avt are used as workspace
216         tn (:,:,:) = tdta  (:,:,:,2)
217         sn (:,:,:) = sdta  (:,:,:,2)
218         avt(:,:,:) = avtdta(:,:,:,2)
219         
220         CALL eos( tn, sn, rhd, rhop )   ! Time-filtered in situ density
221         CALL bn2( tn, sn, rn2 )         ! before Brunt-Vaisala frequency
222         IF( ln_zps )   &
223            &   CALL zps_hde( kt, tn , sn , rhd,  &  ! Partial steps: before Horizontal DErivative
224            &                 gtu, gsu, gru,  &  ! of t, s, rd at the bottom ocean level
225            &                 gtv, gsv, grv )
226         CALL zdf_mxl( kt )              ! mixed layer depth
227         CALL ldf_slp( kt, rhd, rn2 )
228         
229         uslpdta (:,:,:,2) = uslp (:,:,:)
230         vslpdta (:,:,:,2) = vslp (:,:,:)
231         wslpidta(:,:,:,2) = wslpi(:,:,:)
232         wslpjdta(:,:,:,2) = wslpj(:,:,:)
233#endif
234         
235         ! swap from record 2 to 1
236         CALL swap_dyn_data
237         
238         iswap = 1        !  indicates swap
239         
240         CALL dynrea( kt, iper )    ! data read for the iper period
241         
242#if defined key_ldfslp
243         ! Computes slopes
244         ! Caution : here tn, sn and avt are used as workspace
245         tn (:,:,:) = tdta  (:,:,:,2)
246         sn (:,:,:) = sdta  (:,:,:,2)
247         avt(:,:,:) = avtdta(:,:,:,2)
248         
249         CALL eos( tn, sn, rhd, rhop )   ! Time-filtered in situ density
250         CALL bn2( tn, sn, rn2 )         ! before Brunt-Vaisala frequency
251         IF( ln_zps )   &
252            &   CALL zps_hde( kt, tn , sn , rhd,  &  ! Partial steps: before Horizontal DErivative
253            &                 gtu, gsu, gru,  &  ! of t, s, rd at the bottom ocean level
254            &                 gtv, gsv, grv )
255         CALL zdf_mxl( kt )              ! mixed layer depth
256         CALL ldf_slp( kt, rhd, rn2 )
257         
258         uslpdta (:,:,:,2) = uslp (:,:,:)
259         vslpdta (:,:,:,2) = vslp (:,:,:)
260         wslpidta(:,:,:,2) = wslpi(:,:,:)
261         wslpjdta(:,:,:,2) = wslpj(:,:,:)
262#endif
263         !
264         lfirdyn=.FALSE.    ! trace the first call
265      ENDIF
266      !
267      ! And now what we have to do at every time step
268      ! check the validity of the period in memory
269      !
270      IF( iperm1 /= ndyn1 ) THEN
271
272         IF( iperm1 == 0. ) THEN
273            IF (lwp) THEN
274               WRITE (numout,*) ' dynamic file is not periodic with periodic interpolation'
275               WRITE (numout,*) ' we take the last value for the last period '
276               WRITE (numout,*) ' iperm1 = 12,   iper = 13  '
277            ENDIF
278            iperm1 = 12
279            iper   = 13
280         ENDIF
281         !
282         ! We have to prepare a new read of data : swap from record 2 to 1
283         !
284         CALL swap_dyn_data
285
286         iswap = 1        !  indicates swap
287         
288         CALL dynrea( kt, iper )    ! data read for the iper period
289
290#if defined key_ldfslp
291         ! Computes slopes
292         ! Caution : here tn, sn and avt are used as workspace
293         tn (:,:,:) = tdta  (:,:,:,2)
294         sn (:,:,:) = sdta  (:,:,:,2)
295         avt(:,:,:) = avtdta(:,:,:,2)
296         
297         CALL eos( tn, sn, rhd, rhop )   ! Time-filtered in situ density
298         CALL bn2( tn, sn, rn2 )         ! before Brunt-Vaisala frequency
299         IF( ln_zps )   &
300            &   CALL zps_hde( kt, tn , sn , rhd,  &  ! Partial steps: before Horizontal DErivative
301            &                 gtu, gsu, gru,  &  ! of t, s, rd at the bottom ocean level
302            &                 gtv, gsv, grv )
303         CALL zdf_mxl( kt )              ! mixed layer depth
304         CALL ldf_slp( kt, rhd, rn2 )
305         
306         uslpdta (:,:,:,2) = uslp (:,:,:)
307         vslpdta (:,:,:,2) = vslp (:,:,:)
308         wslpidta(:,:,:,2) = wslpi(:,:,:)
309         wslpjdta(:,:,:,2) = wslpj(:,:,:)
310#endif
311       
312         ! store the information of the period read
313         ndyn1 = ndyn2
314         ndyn2 = iper
315
316         IF (lwp) THEN
317            WRITE (numout,*) ' dynamics data read for the period ndyn1 =',ndyn1, &
318               &             ' and for the period ndyn2 = ',ndyn2
319            WRITE (numout,*) ' time step is : ', kt
320         END IF
321         !
322      END IF
323      !
324      ! Compute the data at the given time step
325      !----------------------------------------     
326
327      IF( nsptint == 0 ) THEN   
328         ! No spatial interpolation, data are probably correct
329         ! We have to initialize data if we have changed the period         
330         CALL assign_dyn_data         
331      ELSE IF( nsptint == 1 ) THEN
332         ! linear interpolation
333         CALL linear_interp_dyn_data( zweigh )
334      ELSE 
335         ! other interpolation
336         WRITE (numout,*) ' this kind of interpolation do not exist at the moment : we stop'
337         STOP 'dtadyn'         
338      END IF
339     
340      ! In any case, we need rhop
341      CALL eos( tn, sn, rhd, rhop ) 
342     
343#if ! defined key_off_degrad && defined key_traldf_c2d
344      ! In case of 2D varying coefficients, we need aeiv and aeiu
345      IF( lk_traldf_eiv )   CALL ldf_eiv( kt )      ! eddy induced velocity coefficient
346#endif
347   
348   END SUBROUTINE dta_dyn
349
350   SUBROUTINE dynrea( kt, kenr )
351      !!----------------------------------------------------------------------
352      !!                  ***  ROUTINE dynrea  ***
353      !!
354      !! ** Purpose : READ dynamics fiels from OPA9 netcdf output
355      !!
356      !! ** Method : READ the kenr records of DATA and store in
357      !!             in udta(...,2), .... 
358      !!
359      !! ** History : additions : M. Levy et M. Benjelloul jan 2001
360      !!              (netcdf FORMAT)
361      !!              05-03 (O. Aumont and A. El Moussaoui) F90
362      !!              06-07 : (C. Ethe) use of iom module
363      !!----------------------------------------------------------------------
364      !! * Modules used
365      USE iom
366
367      !! * Arguments
368      INTEGER, INTENT( in ) ::   kt, kenr       ! time index
369      !! * Local declarations
370      INTEGER ::  jkenr
371
372      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
373        zu, zv, zw, zt, zs, zavt ,   &     ! 3-D dynamical fields
374        zhdiv                              ! horizontal divergence
375
376      REAL(wp), DIMENSION(jpi,jpj) :: &
377         zemp, zqsr, zmld, zice, zwspd, &
378         ztaux, ztauy
379#if defined key_trcbbl_dif   ||   defined key_trcbbl_adv
380      REAL(wp), DIMENSION(jpi,jpj) :: zbblx, zbbly
381#endif
382
383#if ! defined key_off_degrad && defined key_traldf_c2d
384      REAL(wp), DIMENSION(jpi,jpj) :: zahtw 
385#   if defined key_trcldf_eiv
386      REAL(wp), DIMENSION(jpi,jpj) :: zaeiw 
387#  endif
388#endif
389
390#if defined key_off_degrad
391   REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
392      zahtu, zahtv, zahtw  !  Lateral diffusivity
393# if defined key_trcldf_eiv
394   REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
395      zaeiu, zaeiv, zaeiw  ! G&M coefficient
396# endif
397#endif
398
399      !---------------------------------------------------------------
400      ! 0. Initialization
401     
402      ! cas d'un fichier non periodique : on utilise deux fois le premier et
403      ! le dernier champ temporel
404
405      jkenr = kenr
406
407      IF(lwp) THEN
408         WRITE(numout,*)
409         WRITE(numout,*) 'Dynrea : reading dynamical fields, kenr = ', jkenr
410         WRITE(numout,*) ' ~~~~~~~'
411#if defined key_off_degrad
412         WRITE(numout,*) ' Degraded fields'
413#endif
414         WRITE(numout,*)
415      ENDIF
416
417
418      IF( kt == nit000 .AND. nlecoff == 0 ) THEN
419         nlecoff = 1
420         CALL  iom_open ( cfile_grid_T, numfl_t )
421         CALL  iom_open ( cfile_grid_U, numfl_u )
422         CALL  iom_open ( cfile_grid_V, numfl_v )
423         CALL  iom_open ( cfile_grid_W, numfl_w )
424      ENDIF
425
426      ! file grid-T
427      !---------------
428      CALL iom_get( numfl_t, jpdom_data, 'votemper', zt   (:,:,:), jkenr )
429      CALL iom_get( numfl_t, jpdom_data, 'vosaline', zs   (:,:,:), jkenr )
430      CALL iom_get( numfl_t, jpdom_data, 'somixhgt', zmld (:,:  ), jkenr )
431      CALL iom_get( numfl_t, jpdom_data, 'sowaflcd', zemp (:,:  ), jkenr )
432      CALL iom_get( numfl_t, jpdom_data, 'soshfldo', zqsr (:,:  ), jkenr )
433      CALL iom_get( numfl_t, jpdom_data, 'soicecov', zice (:,:  ), jkenr )
434      IF( iom_varid( numfl_t, 'sowindsp', ldstop = .FALSE. ) > 0 ) THEN
435         CALL iom_get( numfl_t, jpdom_data, 'sowindsp', zwspd(:,:), jkenr ) 
436      ELSE
437         CALL iom_get( numfl_u, jpdom_data, 'sozotaux', ztaux(:,:), jkenr )
438         CALL iom_get( numfl_v, jpdom_data, 'sometauy', ztauy(:,:), jkenr )
439         CALL tau2wnd( ztaux, ztauy, zwspd )
440      ENDIF
441      ! files grid-U / grid_V
442      CALL iom_get( numfl_u, jpdom_data, 'vozocrtx', zu   (:,:,:), jkenr )
443      CALL iom_get( numfl_v, jpdom_data, 'vomecrty', zv   (:,:,:), jkenr )
444
445#if defined key_trcbbl_dif || defined key_trcbbl_adv
446      IF( iom_varid( numfl_u, 'sobblcox', ldstop = .FALSE. ) > 0  .AND. &
447      &   iom_varid( numfl_v, 'sobblcoy', ldstop = .FALSE. ) > 0 ) THEN
448         CALL iom_get( numfl_u, jpdom_data, 'sobblcox', zbblx(:,:), jkenr )
449         CALL iom_get( numfl_v, jpdom_data, 'sobblcoy', zbbly(:,:), jkenr )
450      ELSE
451         CALL bbl_sign( zt, zs, zbblx, zbbly )   
452      ENDIF
453#endif
454
455      ! file grid-W
456!!      CALL iom_get ( numfl_w, jpdom_data, 'vovecrtz', zw   (:,:,:), jkenr )
457      ! Computation of vertical velocity using horizontal divergence
458      CALL wzv( zu, zv, zw, zhdiv )
459
460# if defined key_zdfddm
461      CALL iom_get( numfl_w, jpdom_data, 'voddmavs', zavt (:,:,:), jkenr )
462#else
463      CALL iom_get( numfl_w, jpdom_data, 'votkeavt', zavt (:,:,:), jkenr )
464#endif 
465
466#if ! defined key_off_degrad && defined key_traldf_c2d
467      CALL iom_get( numfl_w, jpdom_data, 'soleahtw', zahtw (:,: ), jkenr )
468#  if   defined key_trcldf_eiv 
469      CALL iom_get( numfl_w, jpdom_data, 'soleaeiw', zaeiw (:,: ), jkenr )
470#  endif
471#endif
472
473#if defined key_off_degrad
474      CALL iom_get( numfl_u, jpdom_data, 'vozoahtu', zahtu(:,:,:), jkenr )
475      CALL iom_get( numfl_v, jpdom_data, 'vomeahtv', zahtv(:,:,:), jkenr )
476      CALL iom_get( numfl_w, jpdom_data, 'voveahtw', zahtw(:,:,:), jkenr )
477#  if defined key_trcldf_eiv
478      CALL iom_get( numfl_u, jpdom_data, 'vozoaeiu', zaeiu(:,:,:), jkenr )
479      CALL iom_get( numfl_v, jpdom_data, 'vomeaeiv', zaeiv(:,:,:), jkenr )
480      CALL iom_get( numfl_w, jpdom_data, 'voveaeiw', zaeiw(:,:,:), jkenr )
481#  endif
482#endif
483
484      udta(:,:,:,2) = zu(:,:,:) * umask(:,:,:)
485      vdta(:,:,:,2) = zv(:,:,:) * vmask(:,:,:) 
486      wdta(:,:,:,2) = zw(:,:,:) * tmask(:,:,:)
487
488#if defined key_trc_diatrd
489      hdivdta(:,:,:,2) = zhdiv(:,:,:) * tmask(:,:,:)
490#endif
491
492      tdta(:,:,:,2)   = zt  (:,:,:) * tmask(:,:,:)
493      sdta(:,:,:,2)   = zs  (:,:,:) * tmask(:,:,:)
494      avtdta(:,:,:,2) = zavt(:,:,:) * tmask(:,:,:)
495
496#if ! defined key_off_degrad && defined key_traldf_c2d
497      ahtwdta(:,:,2)  = zahtw(:,:) * tmask(:,:,1)
498#if defined key_trcldf_eiv
499      aeiwdta(:,:,2)  = zaeiw(:,:) * tmask(:,:,1)
500#endif
501#endif
502
503#if defined key_off_degrad
504        ahtudta(:,:,:,2) = zahtu(:,:,:) * umask(:,:,:)
505        ahtvdta(:,:,:,2) = zahtv(:,:,:) * vmask(:,:,:)
506        ahtwdta(:,:,:,2) = zahtw(:,:,:) * tmask(:,:,:)
507#  if defined key_trcldf_eiv
508        aeiudta(:,:,:,2) = zaeiu(:,:,:) * umask(:,:,:)
509        aeivdta(:,:,:,2) = zaeiv(:,:,:) * vmask(:,:,:)
510        aeiwdta(:,:,:,2) = zaeiw(:,:,:) * tmask(:,:,:)
511#  endif
512#endif
513
514      ! fluxes
515      !
516      wspddta(:,:,2)  = zwspd(:,:) * tmask(:,:,1)
517      frlddta(:,:,2)  = MIN( 1., zice(:,:) ) * tmask(:,:,1)
518      empdta (:,:,2)  = zemp(:,:) * tmask(:,:,1)
519      qsrdta (:,:,2)  = zqsr(:,:) * tmask(:,:,1)
520      hmlddta(:,:,2)  = zmld(:,:) * tmask(:,:,1)
521     
522#if defined key_trcbbl_dif   ||   defined key_trcbbl_adv
523      bblxdta(:,:,2) = MAX( 0., zbblx(:,:) )
524      bblydta(:,:,2) = MAX( 0., zbbly(:,:) )
525
526      WHERE( bblxdta(:,:,2) > 2. ) bblxdta(:,:,2) = 0.
527      WHERE( bblydta(:,:,2) > 2. ) bblydta(:,:,2) = 0.
528#endif
529
530      IF( kt == nitend ) THEN
531         CALL iom_close ( numfl_t )
532         CALL iom_close ( numfl_u )
533         CALL iom_close ( numfl_v )
534         CALL iom_close ( numfl_w )
535      ENDIF
536     
537   END SUBROUTINE dynrea
538
539   SUBROUTINE dta_dyn_init
540      !!----------------------------------------------------------------------
541      !!                  ***  ROUTINE dta_dyn_init  ***
542      !!
543      !! ** Purpose :   initializations of parameters for the interpolation
544      !!
545      !! ** Method :
546      !!
547      !! History :
548      !!    ! original  : 92-01 (M. Imbard: sub domain)
549      !!    ! 98-04 (L.Bopp MA Foujols: slopes for isopyc.)
550      !!    ! 98-05 (L. Bopp read output of coupled run)
551      !!    ! 05-03 (O. Aumont and A. El Moussaoui) F90
552      !!----------------------------------------------------------------------
553      !! * Modules used
554
555      !! * Local declarations
556
557      REAL(wp) ::   znspyr   !: number of time step per year
558
559      NAMELIST/namdyn/ ndtadyn, ndtatot, nsptint, nficdyn, lperdyn,  &
560      &                cfile_grid_T, cfile_grid_U, cfile_grid_V, cfile_grid_W
561      !!----------------------------------------------------------------------
562
563      !  Define the dynamical input parameters
564      ! ======================================
565
566      ! Read Namelist namdyn : Lateral physics on tracers
567      REWIND( numnam )
568      READ  ( numnam, namdyn )
569
570      IF(lwp) THEN
571         WRITE(numout,*)
572         WRITE(numout,*) 'namdyn : offline dynamical selection'
573         WRITE(numout,*) '~~~~~~~'
574         WRITE(numout,*) '  Namelist namdyn : set parameters for the lecture of the dynamical fields'
575         WRITE(numout,*) 
576         WRITE(numout,*) ' number of elements in the FILE for a year  ndtadyn = ' , ndtadyn
577         WRITE(numout,*) ' total number of elements in the FILE       ndtatot = ' , ndtatot
578         WRITE(numout,*) ' type of interpolation                      nsptint = ' , nsptint
579         WRITE(numout,*) ' number of dynamics FILE                    nficdyn = ' , nficdyn
580         WRITE(numout,*) ' loop on the same FILE                      lperdyn = ' , lperdyn
581         WRITE(numout,*) '  '
582         WRITE(numout,*) ' name of grid_T file                   cfile_grid_T = ', TRIM(cfile_grid_T)   
583         WRITE(numout,*) ' name of grid_U file                   cfile_grid_U = ', TRIM(cfile_grid_U) 
584         WRITE(numout,*) ' name of grid_V file                   cfile_grid_V = ', TRIM(cfile_grid_V) 
585         WRITE(numout,*) ' name of grid_W file                   cfile_grid_W = ', TRIM(cfile_grid_W)     
586         WRITE(numout,*) ' '
587      ENDIF
588
589      znspyr   = nyear_len(1) * rday / rdt 
590      rnspdta  = znspyr / FLOAT( ndtadyn )
591      rnspdta2 = rnspdta * 0.5 
592
593   END SUBROUTINE dta_dyn_init
594
595   SUBROUTINE wzv( pu, pv, pw, phdiv )
596      !!----------------------------------------------------------------------
597      !!                    ***  ROUTINE wzv  ***
598      !!
599      !! ** Purpose :   Compute the now vertical velocity after the array swap
600      !!
601      !! ** Method  :
602      !! ** Method  : - Divergence:
603      !!      - compute the now divergence given by :
604      !!         * z-coordinate
605      !!         hdiv = 1/(e1t*e2t) [ di(e2u  u) + dj(e1v  v) ]
606      !!     - Using the incompressibility hypothesis, the vertical
607      !!      velocity is computed by integrating the horizontal divergence
608      !!      from the bottom to the surface.
609      !!        The boundary conditions are w=0 at the bottom (no flux) and,
610      !!      in regid-lid case, w=0 at the sea surface.
611      !!
612      !!
613      !! History :
614      !!   9.0  !  02-07  (G. Madec)  Vector optimization
615      !!----------------------------------------------------------------------
616      !! * Arguments
617      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in  )   :: pu, pv    !:  horizontal velocities
618      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( out )   :: pw        !:  verticla velocity
619      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( inout ) :: phdiv     !:  horizontal divergence
620
621      !! * Local declarations
622      INTEGER  ::  ji, jj, jk
623      REAL(wp) ::  zu, zu1, zv, zv1, zet
624
625
626      ! Computation of vertical velocity using horizontal divergence
627      phdiv(:,:,:) = 0.
628      DO jk = 1, jpkm1
629         DO jj = 2, jpjm1
630            DO ji = fs_2, fs_jpim1   ! vector opt.
631#if defined key_zco
632               zu  = pu(ji  ,jj  ,jk) * umask(ji  ,jj  ,jk) * e2u(ji  ,jj  )
633               zu1 = pu(ji-1,jj  ,jk) * umask(ji-1,jj  ,jk) * e2u(ji-1,jj  )
634               zv  = pv(ji  ,jj  ,jk) * vmask(ji  ,jj  ,jk) * e1v(ji  ,jj  )
635               zv1 = pv(ji  ,jj-1,jk) * vmask(ji  ,jj-1,jk) * e1v(ji  ,jj-1)
636               zet = 1. / ( e1t(ji,jj) * e2t(ji,jj) )
637#else
638               zu  = pu(ji  ,jj  ,jk) * umask(ji  ,jj  ,jk) * e2u(ji  ,jj  ) * fse3u(ji  ,jj  ,jk)
639               zu1 = pu(ji-1,jj  ,jk) * umask(ji-1,jj  ,jk) * e2u(ji-1,jj  ) * fse3u(ji-1,jj  ,jk)
640               zv  = pv(ji  ,jj  ,jk) * vmask(ji  ,jj  ,jk) * e1v(ji  ,jj  ) * fse3v(ji  ,jj  ,jk)
641               zv1 = pv(ji  ,jj-1,jk) * vmask(ji  ,jj-1,jk) * e1v(ji  ,jj-1) * fse3v(ji  ,jj-1,jk)
642               zet = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
643#endif
644               phdiv(ji,jj,jk) = ( zu - zu1 + zv - zv1 ) * zet 
645            END DO
646         END DO
647      ENDDO
648
649      ! Lateral boundary conditions on phdiv
650      CALL lbc_lnk( phdiv, 'T', 1. )
651
652
653      ! computation of vertical velocity from the bottom
654      pw(:,:,jpk) = 0.
655      DO jk = jpkm1, 1, -1
656         pw(:,:,jk) = pw(:,:,jk+1) - fse3t(:,:,jk) * phdiv(:,:,jk)
657      END DO
658
659   END SUBROUTINE wzv
660
661   SUBROUTINE tau2wnd( ptaux, ptauy, pwspd )
662      !!---------------------------------------------------------------------
663      !!                    ***  ROUTINE sbc_tau2wnd  ***
664      !!
665      !! ** Purpose : Estimation of wind speed as a function of wind stress
666      !!
667      !! ** Method  : |tau|=rhoa*Cd*|U|^2
668      !!---------------------------------------------------------------------
669      !! * Arguments
670      REAL(wp), DIMENSION(jpi,jpj), INTENT( in  ) ::  &
671         ptaux, ptauy                              !: wind stress in i-j direction resp.
672      REAL(wp), DIMENSION(jpi,jpj), INTENT( out ) ::  &
673         pwspd                                     !: wind speed
674      REAL(wp) ::   zrhoa  = 1.22         ! Air density kg/m3
675      REAL(wp) ::   zcdrag = 1.5e-3       ! drag coefficient
676      REAL(wp) ::   ztx, zty, ztau, zcoef ! temporary variables
677      INTEGER  ::   ji, jj                ! dummy indices
678      !!---------------------------------------------------------------------
679      zcoef = 1. / ( zrhoa * zcdrag )
680!CDIR NOVERRCHK
681      DO jj = 2, jpjm1
682!CDIR NOVERRCHK
683         DO ji = fs_2, fs_jpim1   ! vector opt.
684            ztx = ptaux(ji,jj) * umask(ji,jj,1) + ptaux(ji-1,jj  ) * umask(ji-1,jj  ,1)
685            zty = ptauy(ji,jj) * vmask(ji,jj,1) + ptauy(ji  ,jj-1) * vmask(ji  ,jj-1,1)
686            ztau = 0.5 * SQRT( ztx * ztx + zty * zty )
687            pwspd(ji,jj) = SQRT ( ztau * zcoef ) * tmask(ji,jj,1)
688         END DO
689      END DO
690      CALL lbc_lnk( pwspd(:,:), 'T', 1. )
691
692   END SUBROUTINE tau2wnd
693
694#if defined key_trcbbl_dif   ||   defined key_trcbbl_adv
695
696   SUBROUTINE bbl_sign( ptn, psn, pbblx, pbbly )
697      !!----------------------------------------------------------------------
698      !!                    ***  ROUTINE bbl_sign  ***
699      !!
700      !! ** Purpose :   Compute the sign of local gradient of density multiplied by the slope
701      !!                along the bottom slope gradient : grad( rho) * grad(h)
702      !!                Need to compute the diffusive bottom boundary layer
703      !!
704      !! ** Method  :   When the product grad( rho) * grad(h) < 0 (where grad
705      !!      is an along bottom slope gradient) an additional lateral diffu-
706      !!      sive trend along the bottom slope is added to the general tracer
707      !!      trend, otherwise nothing is done. See trcbbl.F90
708      !!
709      !!
710      !! History :
711      !!   9.0  !  02-07  (G. Madec)  Vector optimization
712      !!----------------------------------------------------------------------
713      !! * Arguments
714      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in  ) ::  &
715         ptn             ,  &                           !: temperature
716         psn                                            !: salinity
717      REAL(wp), DIMENSION(jpi,jpj), INTENT( out ) ::  &
718         pbblx , pbbly                                  !: sign of bbl in i-j direction resp.
719     
720      !! * Local declarations
721      INTEGER  ::   ji, jj                   ! dummy loop indices
722      INTEGER  ::   ik
723      REAL(wp) ::   &
724         ztx, zsx, zhx, zalbetx, zgdrhox,     &  ! temporary scalars
725         zty, zsy, zhy, zalbety, zgdrhoy 
726      REAL(wp), DIMENSION(jpi,jpj) ::    &
727        ztnb, zsnb, zdep
728      REAL(wp) ::    fsalbt, pft, pfs, pfh   ! statement function
729      !!----------------------------------------------------------------------
730      ! ratio alpha/beta
731      ! ================
732      !  fsalbt: ratio of thermal over saline expension coefficients
733      !       pft :  potential temperature in degrees celcius
734      !       pfs :  salinity anomaly (s-35) in psu
735      !       pfh :  depth in meters
736
737      fsalbt( pft, pfs, pfh ) =                                              &
738         ( ( ( -0.255019e-07 * pft + 0.298357e-05 ) * pft                    &
739                                   - 0.203814e-03 ) * pft                    &
740                                   + 0.170907e-01 ) * pft                    &
741                                   + 0.665157e-01                            &
742         +(-0.678662e-05 * pfs - 0.846960e-04 * pft + 0.378110e-02 ) * pfs   &
743         +  ( ( - 0.302285e-13 * pfh                                         &
744                - 0.251520e-11 * pfs                                         &
745                + 0.512857e-12 * pft * pft          ) * pfh                  &
746                                     - 0.164759e-06   * pfs                  &
747             +(   0.791325e-08 * pft - 0.933746e-06 ) * pft                  &
748                                     + 0.380374e-04 ) * pfh
749
750      ! 0. 2D fields of bottom temperature and salinity, and bottom slope
751      ! -----------------------------------------------------------------
752      ! mbathy= number of w-level, minimum value=1 (cf domrea.F90)
753#  if defined key_vectopt_loop
754      jj = 1
755      DO ji = 1, jpij   ! vector opt. (forced unrolling)
756#  else
757      DO jj = 1, jpj
758         DO ji = 1, jpi
759#  endif
760            ik          =  MAX( mbathy(ji,jj) - 1, 1 )    ! vertical index of the bottom ocean T-level
761            ztnb(ji,jj) = ptn(ji,jj,ik) * tmask(ji,jj,1)  ! masked T and S at ocean bottom
762            zsnb(ji,jj) = psn(ji,jj,ik) * tmask(ji,jj,1)
763            zdep(ji,jj) = fsdept(ji,jj,ik)                ! depth of the ocean bottom T-level
764#  if ! defined key_vectopt_loop
765         END DO
766#  endif
767      END DO
768
769      !!----------------------------------------------------------------------
770      ! 1. Criteria of additional bottom diffusivity: grad(rho).grad(h)<0
771      ! --------------------------------------------
772      ! Sign of the local density gradient along the i- and j-slopes
773      ! multiplied by the slope of the ocean bottom
774
775      SELECT CASE ( neos )
776
777      CASE ( 0 )                 ! Jackett and McDougall (1994) formulation
778
779#  if defined key_vectopt_loop
780      jj = 1
781      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
782#  else
783      DO jj = 1, jpjm1
784         DO ji = 1, jpim1
785#  endif
786            ! temperature, salinity anomalie and depth
787            ztx = 0.5 * ( ztnb(ji,jj) + ztnb(ji+1,jj) )
788            zsx = 0.5 * ( zsnb(ji,jj) + zsnb(ji+1,jj) ) - 35.0
789            zhx = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) )
790            !
791            zty = 0.5 * ( ztnb(ji,jj+1) + ztnb(ji,jj) )
792            zsy = 0.5 * ( zsnb(ji,jj+1) + zsnb(ji,jj) ) - 35.0
793            zhy = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) )
794            ! masked ratio alpha/beta
795            zalbetx = fsalbt( ztx, zsx, zhx ) * umask(ji,jj,1)
796            zalbety = fsalbt( zty, zsy, zhy ) * vmask(ji,jj,1)
797            ! local density gradient along i-bathymetric slope
798            zgdrhox = zalbetx * ( ztnb(ji+1,jj) - ztnb(ji,jj) )   &
799                   -            ( zsnb(ji+1,jj) - zsnb(ji,jj) )
800            ! local density gradient along j-bathymetric slope
801            zgdrhoy = zalbety * ( ztnb(ji,jj+1) - ztnb(ji,jj) )   &
802                   -            ( zsnb(ji,jj+1) - zsnb(ji,jj) )
803            ! sign of local i-gradient of density multiplied by the i-slope
804            pbblx(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhox * ( zdep(ji+1,jj) - zdep(ji,jj) ) )
805            ! sign of local j-gradient of density multiplied by the j-slope
806            pbbly(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhoy * ( zdep(ji,jj+1) - zdep(ji,jj) ) )
807#  if ! defined key_vectopt_loop
808         END DO
809#  endif
810      END DO
811
812      CASE ( 1 )               ! Linear formulation function of temperature only
813                               !
814#  if defined key_vectopt_loop
815      jj = 1
816      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
817#  else
818      DO jj = 1, jpjm1
819         DO ji = 1, jpim1
820#  endif
821            ! local 'density/temperature' gradient along i-bathymetric slope
822            zgdrhox =  ztnb(ji+1,jj) - ztnb(ji,jj)
823            ! local density gradient along j-bathymetric slope
824            zgdrhoy =  ztnb(ji,jj+1) - ztnb(ji,jj)
825            ! sign of local i-gradient of density multiplied by the i-slope
826            pbblx(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhox * ( zdep(ji+1,jj) - zdep(ji,jj) ) )
827            ! sign of local j-gradient of density multiplied by the j-slope
828            pbbly(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhoy * ( zdep(ji,jj+1) - zdep(ji,jj) ) )
829#  if ! defined key_vectopt_loop
830         END DO
831#  endif
832      END DO
833
834      CASE ( 2 )               ! Linear formulation function of temperature and salinity
835
836#  if defined key_vectopt_loop
837      jj = 1
838      DO ji = 1, jpij-jpi   ! vector opt. (forced unrolling)
839#  else
840      DO jj = 1, jpjm1
841         DO ji = 1, jpim1
842#  endif     
843            ! local density gradient along i-bathymetric slope
844            zgdrhox = - ( rbeta*( zsnb(ji+1,jj) - zsnb(ji,jj) )   &
845                      -  ralpha*( ztnb(ji+1,jj) - ztnb(ji,jj) ) )
846            ! local density gradient along j-bathymetric slope
847            zgdrhoy = - ( rbeta*( zsnb(ji,jj+1) - zsnb(ji,jj) )   &
848                      -  ralpha*( ztnb(ji,jj+1) - ztnb(ji,jj) ) )
849            ! sign of local i-gradient of density multiplied by the i-slope
850            pbblx(ji,jj) = 0.5 - SIGN( 0.5, - zgdrhox * ( zdep(ji+1,jj) - zdep(ji,jj) ) )
851            ! sign of local j-gradient of density multiplied by the j-slope
852            pbbly(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhoy * ( zdep(ji,jj+1) - zdep(ji,jj) ) )
853#  if ! defined key_vectopt_loop
854         END DO
855#  endif
856      END DO
857
858      CASE DEFAULT
859
860         WRITE(ctmp1,*) '          bad flag value for neos = ', neos
861         CALL ctl_stop(ctmp1)
862
863      END SELECT
864   
865      ! Lateral boundary conditions
866      CALL lbc_lnk( pbblx, 'U', 1. )
867      CALL lbc_lnk( pbbly, 'V', 1. )
868
869   END SUBROUTINE bbl_sign
870
871#endif
872
873   SUBROUTINE swap_dyn_data
874      !!----------------------------------------------------------------------
875      !!                    ***  ROUTINE swap_dyn_data  ***
876      !!
877      !! ** Purpose :   swap array data
878      !!
879      !! History :
880      !!   9.0  !  07-09  (C. Ethe)
881      !!----------------------------------------------------------------------
882
883
884      ! swap from record 2 to 1
885      tdta   (:,:,:,1) = tdta   (:,:,:,2)
886      sdta   (:,:,:,1) = sdta   (:,:,:,2)
887      avtdta (:,:,:,1) = avtdta (:,:,:,2)
888      udta   (:,:,:,1) = udta   (:,:,:,2)
889      vdta   (:,:,:,1) = vdta   (:,:,:,2)
890      wdta   (:,:,:,1) = wdta   (:,:,:,2)
891#if defined key_trc_diatrd
892      hdivdta(:,:,:,1) = hdivdta(:,:,:,2)
893#endif
894
895#if defined key_ldfslp
896      uslpdta (:,:,:,1) = uslpdta (:,:,:,2)
897      vslpdta (:,:,:,1) = vslpdta (:,:,:,2)
898      wslpidta(:,:,:,1) = wslpidta(:,:,:,2)
899      wslpjdta(:,:,:,1) = wslpjdta(:,:,:,2)
900#endif
901      hmlddta(:,:,1) = hmlddta(:,:,2) 
902      wspddta(:,:,1) = wspddta(:,:,2) 
903      frlddta(:,:,1) = frlddta(:,:,2) 
904      empdta (:,:,1) = empdta (:,:,2) 
905      qsrdta (:,:,1) = qsrdta (:,:,2) 
906
907#if ! defined key_off_degrad && defined key_traldf_c2d
908      ahtwdta(:,:,1) = ahtwdta(:,:,2)
909#  if defined key_trcldf_eiv
910      aeiwdta(:,:,1) = aeiwdta(:,:,2)
911#  endif
912#endif
913
914#if defined key_off_degrad
915      ahtudta(:,:,:,1) = ahtudta(:,:,:,2)
916      ahtvdta(:,:,:,1) = ahtvdta(:,:,:,2)
917      ahtwdta(:,:,:,1) = ahtwdta(:,:,:,2)
918#  if defined key_trcldf_eiv
919      aeiudta(:,:,:,1) = aeiudta(:,:,:,2)
920      aeivdta(:,:,:,1) = aeivdta(:,:,:,2)
921      aeiwdta(:,:,:,1) = aeiwdta(:,:,:,2)
922#  endif
923#endif
924
925#if defined key_trcbbl_dif || defined key_trcbbl_adv
926      bblxdta(:,:,1) = bblxdta(:,:,2)
927      bblydta(:,:,1) = bblydta(:,:,2)
928#endif
929
930   END SUBROUTINE swap_dyn_data
931
932   SUBROUTINE assign_dyn_data
933      !!----------------------------------------------------------------------
934      !!                    ***  ROUTINE assign_dyn_data  ***
935      !!
936      !! ** Purpose :   Assign dynamical data to the data that have been read
937      !!                without time interpolation
938      !!
939      !!----------------------------------------------------------------------
940     
941      tn (:,:,:) = tdta  (:,:,:,2)
942      sn (:,:,:) = sdta  (:,:,:,2)
943      avt(:,:,:) = avtdta(:,:,:,2)
944     
945      un (:,:,:) = udta  (:,:,:,2) 
946      vn (:,:,:) = vdta  (:,:,:,2)
947      wn (:,:,:) = wdta  (:,:,:,2)
948
949#if defined key_trc_diatrd
950      hdivn(:,:,:) = hdivdta(:,:,:,2)
951#endif
952
953#if defined key_zdfddm
954      avs(:,:,:)   = avtdta (:,:,:,2)
955#endif
956
957     
958#if defined key_ldfslp
959      uslp (:,:,:) = uslpdta (:,:,:,2) 
960      vslp (:,:,:) = vslpdta (:,:,:,2) 
961      wslpi(:,:,:) = wslpidta(:,:,:,2) 
962      wslpj(:,:,:) = wslpjdta(:,:,:,2) 
963#endif
964
965      hmld(:,:) = hmlddta(:,:,2) 
966      wndm(:,:) = wspddta(:,:,2) 
967      fr_i(:,:) = frlddta(:,:,2) 
968      emp (:,:) = empdta (:,:,2) 
969      emps(:,:) = emp(:,:) 
970      qsr (:,:) = qsrdta (:,:,2) 
971
972#if ! defined key_off_degrad && defined key_traldf_c2d   
973      ahtw(:,:) = ahtwdta(:,:,2)
974#  if defined key_trcldf_eiv
975      aeiw(:,:) = aeiwdta(:,:,2)
976#  endif
977#endif
978     
979#if defined key_off_degrad
980      ahtu(:,:,:) = ahtudta(:,:,:,2)
981      ahtv(:,:,:) = ahtvdta(:,:,:,2)
982      ahtw(:,:,:) = ahtwdta(:,:,:,2)
983#  if defined key_trcldf_eiv
984      aeiu(:,:,:) = aeiudta(:,:,:,2)
985      aeiv(:,:,:) = aeivdta(:,:,:,2)
986      aeiw(:,:,:) = aeiwdta(:,:,:,2)
987#  endif
988     
989#endif
990     
991#if defined key_trcbbl_dif ||  defined key_trcbbl_adv
992      bblx(:,:) = bblxdta(:,:,2)
993      bbly(:,:) = bblydta(:,:,2)
994#endif
995
996   END SUBROUTINE assign_dyn_data
997
998   SUBROUTINE linear_interp_dyn_data( pweigh )
999      !!----------------------------------------------------------------------
1000      !!                    ***  ROUTINE linear_interp_dyn_data  ***
1001      !!
1002      !! ** Purpose :   linear interpolation of data
1003      !!
1004      !!----------------------------------------------------------------------
1005      !! * Argument
1006      REAL(wp), INTENT( in ) ::   pweigh       ! weigh
1007
1008      !! * Local declarations
1009      REAL(wp) :: zweighm1
1010      !!----------------------------------------------------------------------
1011
1012      zweighm1 = 1. - pweigh
1013     
1014      tn (:,:,:) = zweighm1 * tdta  (:,:,:,1) + pweigh * tdta  (:,:,:,2)
1015      sn (:,:,:) = zweighm1 * sdta  (:,:,:,1) + pweigh * sdta  (:,:,:,2)
1016      avt(:,:,:) = zweighm1 * avtdta(:,:,:,1) + pweigh * avtdta(:,:,:,2)
1017     
1018      un (:,:,:) = zweighm1 * udta  (:,:,:,1) + pweigh * udta  (:,:,:,2) 
1019      vn (:,:,:) = zweighm1 * vdta  (:,:,:,1) + pweigh * vdta  (:,:,:,2)
1020      wn (:,:,:) = zweighm1 * wdta  (:,:,:,1) + pweigh * wdta  (:,:,:,2)
1021
1022#if defined key_trc_diatrd
1023      hdivn(:,:,:) = zweighm1 * hdivdta(:,:,:,1) + pweigh * hdivdta(:,:,:,2)
1024#endif
1025
1026#if defined key_zdfddm
1027      avs(:,:,:)   = zweighm1 * avtdta (:,:,:,1) + pweigh * avtdta (:,:,:,2)
1028#endif
1029
1030     
1031#if defined key_ldfslp
1032      uslp (:,:,:) = zweighm1 * uslpdta (:,:,:,1) + pweigh * uslpdta (:,:,:,2) 
1033      vslp (:,:,:) = zweighm1 * vslpdta (:,:,:,1) + pweigh * vslpdta (:,:,:,2) 
1034      wslpi(:,:,:) = zweighm1 * wslpidta(:,:,:,1) + pweigh * wslpidta(:,:,:,2) 
1035      wslpj(:,:,:) = zweighm1 * wslpjdta(:,:,:,1) + pweigh * wslpjdta(:,:,:,2) 
1036#endif
1037
1038      hmld(:,:) = zweighm1 * hmlddta(:,:,1) + pweigh  * hmlddta(:,:,2) 
1039      wndm(:,:) = zweighm1 * wspddta(:,:,1) + pweigh  * wspddta(:,:,2) 
1040      fr_i(:,:) = zweighm1 * frlddta(:,:,1) + pweigh  * frlddta(:,:,2) 
1041      emp (:,:) = zweighm1 * empdta (:,:,1) + pweigh  * empdta (:,:,2) 
1042      emps(:,:) = emp(:,:) 
1043      qsr (:,:) = zweighm1 * qsrdta (:,:,1) + pweigh  * qsrdta (:,:,2) 
1044
1045#if ! defined key_off_degrad && defined key_traldf_c2d   
1046      ahtw(:,:) = zweighm1 * ahtwdta(:,:,1) + pweigh * ahtwdta(:,:,2)
1047#  if defined key_trcldf_eiv
1048      aeiw(:,:) = zweighm1 * aeiwdta(:,:,1) + pweigh * aeiwdta(:,:,2)
1049#  endif
1050#endif
1051     
1052#if defined key_off_degrad
1053      ahtu(:,:,:) = zweighm1 * ahtudta(:,:,:,1) + pweigh * ahtudta(:,:,:,2)
1054      ahtv(:,:,:) = zweighm1 * ahtvdta(:,:,:,1) + pweigh * ahtvdta(:,:,:,2)
1055      ahtw(:,:,:) = zweighm1 * ahtwdta(:,:,:,1) + pweigh * ahtwdta(:,:,:,2)
1056#  if defined key_trcldf_eiv
1057      aeiu(:,:,:) = zweighm1 * aeiudta(:,:,:,1) + pweigh * aeiudta(:,:,:,2)
1058      aeiv(:,:,:) = zweighm1 * aeivdta(:,:,:,1) + pweigh * aeivdta(:,:,:,2)
1059      aeiw(:,:,:) = zweighm1 * aeiwdta(:,:,:,1) + pweigh * aeiwdta(:,:,:,2)
1060#  endif
1061#endif
1062     
1063#if defined key_trcbbl_dif   ||   defined key_trcbbl_adv
1064      bblx(:,:) = zweighm1 * bblxdta(:,:,1) + pweigh * bblxdta(:,:,2)
1065      bbly(:,:) = zweighm1 * bblydta(:,:,1) + pweigh * bblydta(:,:,2)
1066#endif
1067
1068   END SUBROUTINE linear_interp_dyn_data
1069
1070END MODULE dtadyn
Note: See TracBrowser for help on using the repository browser.