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.
zdfmxl.F90 in branches/UKMO/AMM15_v3_6_STABLE_package/NEMOGCM/NEMO/OPA_SRC/ZDF – NEMO

source: branches/UKMO/AMM15_v3_6_STABLE_package/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfmxl.F90 @ 8059

Last change on this file since 8059 was 8059, checked in by jgraham, 7 years ago

Merged branches required for AMM15 simulations, see ticket #1904.
Merged branches include:
branches/UKMO/CO6_KD490
branches/UKMO/CO6_Restartable_Tidal_Analysis
branches/UKMO/AMM15_v3_6_STABLE

File size: 22.2 KB
Line 
1MODULE zdfmxl
2   !!======================================================================
3   !!                       ***  MODULE  zdfmxl  ***
4   !! Ocean physics: mixed layer depth
5   !!======================================================================
6   !! History :  1.0  ! 2003-08  (G. Madec)  original code
7   !!            3.2  ! 2009-07  (S. Masson, G. Madec)  IOM + merge of DO-loop
8   !!            3.7  ! 2012-03  (G. Madec)  make public the density criteria for trdmxl
9   !!             -   ! 2014-02  (F. Roquet)  mixed layer depth calculated using N2 instead of rhop
10   !!----------------------------------------------------------------------
11   !!   zdf_mxl      : Compute the turbocline and mixed layer depths.
12   !!----------------------------------------------------------------------
13   USE oce             ! ocean dynamics and tracers variables
14   USE dom_oce         ! ocean space and time domain variables
15   USE zdf_oce         ! ocean vertical physics
16   USE in_out_manager  ! I/O manager
17   USE prtctl          ! Print control
18   USE phycst          ! physical constants
19   USE iom             ! I/O library
20   USE eosbn2          ! for zdf_mxl_zint
21   USE lib_mpp         ! MPP library
22   USE wrk_nemo        ! work arrays
23   USE timing          ! Timing
24   USE trc_oce, ONLY : lk_offline ! offline flag
25
26   IMPLICIT NONE
27   PRIVATE
28
29   PUBLIC   zdf_mxl       ! called by step.F90
30
31   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   nmln    !: number of level in the mixed layer (used by TOP)
32   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmld    !: mixing layer depth (turbocline)      [m]
33   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlp    !: mixed layer depth  (rho=rho0+zdcrit) [m]
34   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlpt   !: mixed layer depth at t-points        [m]
35   REAL(wp), PUBLIC, ALLOCATABLE,       DIMENSION(:,:) ::   hmld_zint  !: vertically-interpolated mixed layer depth   [m]
36   REAL(wp), PUBLIC, ALLOCATABLE,       DIMENSION(:,:) ::   htc_mld    ! Heat content of hmld_zint
37   LOGICAL, PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)    :: ll_found   ! Is T_b to be found by interpolation ?
38   LOGICAL, PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)  :: ll_belowml ! Flag points below mixed layer when ll_found=F
39
40   REAL(wp), PUBLIC ::   rho_c = 0.01_wp    !: density criterion for mixed layer depth
41   REAL(wp)         ::   avt_c = 5.e-4_wp   ! Kz criterion for the turbocline depth
42
43   TYPE, PUBLIC :: MXL_ZINT   !: Structure for MLD defs
44      INTEGER   :: mld_type   ! mixed layer type     
45      REAL(wp)  :: zref       ! depth of initial T_ref
46      REAL(wp)  :: dT_crit    ! Critical temp diff
47      REAL(wp)  :: iso_frac   ! Fraction of rn_dT_crit used
48   END TYPE MXL_ZINT
49
50   !! * Substitutions
51#  include "domzgr_substitute.h90"
52   !!----------------------------------------------------------------------
53   !! NEMO/OPA 4.0 , NEMO Consortium (2011)
54   !! $Id$
55   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
56   !!----------------------------------------------------------------------
57CONTAINS
58
59   INTEGER FUNCTION zdf_mxl_alloc()
60      !!----------------------------------------------------------------------
61      !!               ***  FUNCTION zdf_mxl_alloc  ***
62      !!----------------------------------------------------------------------
63      zdf_mxl_alloc = 0      ! set to zero if no array to be allocated
64      IF( .NOT. ALLOCATED( nmln ) ) THEN
65         ALLOCATE( nmln(jpi,jpj), hmld(jpi,jpj), hmlp(jpi,jpj), hmlpt(jpi,jpj), hmld_zint(jpi,jpj),       &
66        &          htc_mld(jpi,jpj),                                                                    &
67        &          ll_found(jpi,jpj), ll_belowml(jpi,jpj,jpk), STAT= zdf_mxl_alloc )
68         !
69         IF( lk_mpp             )   CALL mpp_sum ( zdf_mxl_alloc )
70         IF( zdf_mxl_alloc /= 0 )   CALL ctl_warn('zdf_mxl_alloc: failed to allocate arrays.')
71         !
72      ENDIF
73   END FUNCTION zdf_mxl_alloc
74
75
76   SUBROUTINE zdf_mxl( kt )
77      !!----------------------------------------------------------------------
78      !!                  ***  ROUTINE zdfmxl  ***
79      !!                   
80      !! ** Purpose :   Compute the turbocline depth and the mixed layer depth
81      !!              with density criteria.
82      !!
83      !! ** Method  :   The mixed layer depth is the shallowest W depth with
84      !!      the density of the corresponding T point (just bellow) bellow a
85      !!      given value defined locally as rho(10m) + rho_c
86      !!               The turbocline depth is the depth at which the vertical
87      !!      eddy diffusivity coefficient (resulting from the vertical physics
88      !!      alone, not the isopycnal part, see trazdf.F) fall below a given
89      !!      value defined locally (avt_c here taken equal to 5 cm/s2 by default)
90      !!
91      !! ** Action  :   nmln, hmld, hmlp, hmlpt
92      !!----------------------------------------------------------------------
93      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
94      !
95      INTEGER  ::   ji, jj, jk   ! dummy loop indices
96      INTEGER  ::   iikn, iiki, ikt, imkt   ! local integer
97      REAL(wp) ::   zN2_c        ! local scalar
98      INTEGER, POINTER, DIMENSION(:,:) ::   imld   ! 2D workspace
99      !!----------------------------------------------------------------------
100      !
101      IF( nn_timing == 1 )  CALL timing_start('zdf_mxl')
102      !
103      CALL wrk_alloc( jpi,jpj, imld )
104
105      IF( kt == nit000 ) THEN
106         IF(lwp) WRITE(numout,*)
107         IF(lwp) WRITE(numout,*) 'zdf_mxl : mixed layer depth'
108         IF(lwp) WRITE(numout,*) '~~~~~~~ '
109         !                             ! allocate zdfmxl arrays
110         IF( zdf_mxl_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'zdf_mxl : unable to allocate arrays' )
111      ENDIF
112
113      ! w-level of the mixing and mixed layers
114      nmln(:,:)  = nlb10               ! Initialization to the number of w ocean point
115      hmlp(:,:)  = 0._wp               ! here hmlp used as a dummy variable, integrating vertically N^2
116      zN2_c = grav * rho_c * r1_rau0   ! convert density criteria into N^2 criteria
117      DO jk = nlb10, jpkm1
118         DO jj = 1, jpj                ! Mixed layer level: w-level
119            DO ji = 1, jpi
120               ikt = mbkt(ji,jj)
121               hmlp(ji,jj) = hmlp(ji,jj) + MAX( rn2b(ji,jj,jk) , 0._wp ) * fse3w(ji,jj,jk)
122               IF( hmlp(ji,jj) < zN2_c )   nmln(ji,jj) = MIN( jk , ikt ) + 1   ! Mixed layer level
123            END DO
124         END DO
125      END DO
126      !
127      ! w-level of the turbocline
128      imld(:,:) = mbkt(:,:) + 1        ! Initialization to the number of w ocean point
129      DO jk = jpkm1, nlb10, -1         ! from the bottom to nlb10
130         DO jj = 1, jpj
131            DO ji = 1, jpi
132               imkt = mikt(ji,jj)
133               IF( avt (ji,jj,jk) < avt_c )   imld(ji,jj) = MAX( imkt, jk )      ! Turbocline
134            END DO
135         END DO
136      END DO
137      ! depth of the mixing and mixed layers
138      DO jj = 1, jpj
139         DO ji = 1, jpi
140            iiki = imld(ji,jj)
141            iikn = nmln(ji,jj)
142            imkt = mikt(ji,jj)
143            hmld (ji,jj) = ( fsdepw(ji,jj,iiki  ) - fsdepw(ji,jj,imkt )            )   * ssmask(ji,jj)    ! Turbocline depth
144            hmlp (ji,jj) = ( fsdepw(ji,jj,iikn  ) - fsdepw(ji,jj,MAX( imkt,nla10 ) ) ) * ssmask(ji,jj)    ! Mixed layer depth
145            hmlpt(ji,jj) = ( fsdept(ji,jj,iikn-1) - fsdepw(ji,jj,imkt )            )   * ssmask(ji,jj)    ! depth of the last T-point inside the mixed layer
146         END DO
147      END DO
148      IF( .NOT.lk_offline ) THEN            ! no need to output in offline mode
149         CALL iom_put( "mldr10_1", hmlp )   ! mixed layer depth
150         CALL iom_put( "mldkz5"  , hmld )   ! turbocline depth
151      ENDIF
152     
153      ! Vertically-interpolated mixed-layer depth diagnostic
154      CALL zdf_mxl_zint( kt )
155
156      IF(ln_ctl)   CALL prt_ctl( tab2d_1=REAL(nmln,wp), clinfo1=' nmln : ', tab2d_2=hmlp, clinfo2=' hmlp : ', ovlap=1 )
157      !
158      CALL wrk_dealloc( jpi,jpj, imld )
159      !
160      IF( nn_timing == 1 )  CALL timing_stop('zdf_mxl')
161      !
162   END SUBROUTINE zdf_mxl
163
164   SUBROUTINE zdf_mxl_zint_mld( sf ) 
165      !!----------------------------------------------------------------------------------
166      !!                    ***  ROUTINE zdf_mxl_zint_mld  ***
167      !                                                                       
168      !   Calculate vertically-interpolated mixed layer depth diagnostic.
169      !           
170      !   This routine can calculate the mixed layer depth diagnostic suggested by
171      !   Kara et al, 2000, JGR, 105, 16803, but is more general and can calculate
172      !   vertically-interpolated mixed-layer depth diagnostics with other parameter
173      !   settings set in the namzdf_mldzint namelist. 
174      !
175      !   If mld_type=1 the mixed layer depth is calculated as the depth at which the 
176      !   density has increased by an amount equivalent to a temperature difference of 
177      !   0.8C at the surface.
178      !
179      !   For other values of mld_type the mixed layer is calculated as the depth at 
180      !   which the temperature differs by 0.8C from the surface temperature. 
181      !                                                                       
182      !   David Acreman, Daley Calvert                                     
183      !
184      !!-----------------------------------------------------------------------------------
185
186      TYPE(MXL_ZINT), INTENT(in)  :: sf
187
188      ! Diagnostic criteria
189      INTEGER   :: nn_mld_type   ! mixed layer type     
190      REAL(wp)  :: rn_zref       ! depth of initial T_ref
191      REAL(wp)  :: rn_dT_crit    ! Critical temp diff
192      REAL(wp)  :: rn_iso_frac   ! Fraction of rn_dT_crit used
193
194      ! Local variables
195      REAL(wp), PARAMETER :: zepsilon = 1.e-30          ! local small value
196      INTEGER, POINTER, DIMENSION(:,:) :: ikmt          ! number of active tracer levels
197      INTEGER, POINTER, DIMENSION(:,:) :: ik_ref        ! index of reference level
198      INTEGER, POINTER, DIMENSION(:,:) :: ik_iso        ! index of last uniform temp level
199      REAL, POINTER, DIMENSION(:,:,:)  :: zT            ! Temperature or density
200      REAL, POINTER, DIMENSION(:,:)    :: ppzdep        ! depth for use in calculating d(rho)
201      REAL, POINTER, DIMENSION(:,:)    :: zT_ref        ! reference temperature
202      REAL    :: zT_b                                   ! base temperature
203      REAL, POINTER, DIMENSION(:,:,:)  :: zdTdz         ! gradient of zT
204      REAL, POINTER, DIMENSION(:,:,:)  :: zmoddT        ! Absolute temperature difference
205      REAL    :: zdz                                    ! depth difference
206      REAL    :: zdT                                    ! temperature difference
207      REAL, POINTER, DIMENSION(:,:)    :: zdelta_T      ! difference critereon
208      REAL, POINTER, DIMENSION(:,:)    :: zRHO1, zRHO2  ! Densities
209      INTEGER :: ji, jj, jk                             ! loop counter
210
211      !!-------------------------------------------------------------------------------------
212     
213      CALL wrk_alloc( jpi, jpj, ikmt, ik_ref, ik_iso) 
214      CALL wrk_alloc( jpi, jpj, ppzdep, zT_ref, zdelta_T, zRHO1, zRHO2 ) 
215      CALL wrk_alloc( jpi, jpj, jpk, zT, zdTdz, zmoddT ) 
216
217      ! Unpack structure
218      nn_mld_type = sf%mld_type
219      rn_zref     = sf%zref
220      rn_dT_crit  = sf%dT_crit
221      rn_iso_frac = sf%iso_frac
222
223      ! Set the mixed layer depth criterion at each grid point
224      IF( nn_mld_type == 0 ) THEN
225         zdelta_T(:,:) = rn_dT_crit
226         zT(:,:,:) = rhop(:,:,:)
227      ELSE IF( nn_mld_type == 1 ) THEN
228         ppzdep(:,:)=0.0 
229         call eos ( tsn(:,:,1,:), ppzdep(:,:), zRHO1(:,:) ) 
230! Use zT temporarily as a copy of tsn with rn_dT_crit added to SST
231! [assumes number of tracers less than number of vertical levels]
232         zT(:,:,1:jpts)=tsn(:,:,1,1:jpts) 
233         zT(:,:,jp_tem)=zT(:,:,1)+rn_dT_crit 
234         CALL eos( zT(:,:,1:jpts), ppzdep(:,:), zRHO2(:,:) ) 
235         zdelta_T(:,:) = abs( zRHO1(:,:) - zRHO2(:,:) ) * rau0 
236         ! RHO from eos (2d version) doesn't calculate north or east halo:
237         CALL lbc_lnk( zdelta_T, 'T', 1. ) 
238         zT(:,:,:) = rhop(:,:,:) 
239      ELSE
240         zdelta_T(:,:) = rn_dT_crit                     
241         zT(:,:,:) = tsn(:,:,:,jp_tem)                           
242      END IF 
243
244      ! Calculate the gradient of zT and absolute difference for use later
245      DO jk = 1 ,jpk-2 
246         zdTdz(:,:,jk)  =    ( zT(:,:,jk+1) - zT(:,:,jk) ) / fse3w(:,:,jk+1) 
247         zmoddT(:,:,jk) = abs( zT(:,:,jk+1) - zT(:,:,jk) ) 
248      END DO 
249
250      ! Find density/temperature at the reference level (Kara et al use 10m).         
251      ! ik_ref is the index of the box centre immediately above or at the reference level
252      ! Find rn_zref in the array of model level depths and find the ref   
253      ! density/temperature by linear interpolation.                                   
254      DO jk = jpkm1, 2, -1 
255         WHERE ( fsdept(:,:,jk) > rn_zref ) 
256           ik_ref(:,:) = jk - 1 
257           zT_ref(:,:) = zT(:,:,jk-1) + zdTdz(:,:,jk-1) * ( rn_zref - fsdept(:,:,jk-1) ) 
258         END WHERE
259      END DO 
260
261      ! If the first grid box centre is below the reference level then use the
262      ! top model level to get zT_ref
263      WHERE ( fsdept(:,:,1) > rn_zref ) 
264         zT_ref = zT(:,:,1) 
265         ik_ref = 1 
266      END WHERE 
267
268      ! The number of active tracer levels is 1 less than the number of active w levels
269      ikmt(:,:) = mbathy(:,:) - 1 
270
271      ! Initialize / reset
272      ll_found(:,:) = .false.
273
274      IF ( rn_iso_frac - zepsilon > 0. ) THEN
275         ! Search for a uniform density/temperature region where adjacent levels         
276         ! differ by less than rn_iso_frac * deltaT.                                     
277         ! ik_iso is the index of the last level in the uniform layer 
278         ! ll_found indicates whether the mixed layer depth can be found by interpolation
279         ik_iso(:,:)   = ik_ref(:,:) 
280         DO jj = 1, nlcj 
281            DO ji = 1, nlci 
282!CDIR NOVECTOR
283               DO jk = ik_ref(ji,jj), ikmt(ji,jj)-1 
284                  IF ( zmoddT(ji,jj,jk) > ( rn_iso_frac * zdelta_T(ji,jj) ) ) THEN
285                     ik_iso(ji,jj)   = jk 
286                     ll_found(ji,jj) = ( zmoddT(ji,jj,jk) > zdelta_T(ji,jj) ) 
287                     EXIT
288                  END IF
289               END DO
290            END DO
291         END DO 
292
293         ! Use linear interpolation to find depth of mixed layer base where possible
294         hmld_zint(:,:) = rn_zref 
295         DO jj = 1, jpj 
296            DO ji = 1, jpi 
297               IF (ll_found(ji,jj) .and. tmask(ji,jj,1) == 1.0) THEN
298                  zdz =  abs( zdelta_T(ji,jj) / zdTdz(ji,jj,ik_iso(ji,jj)) ) 
299                  hmld_zint(ji,jj) = fsdept(ji,jj,ik_iso(ji,jj)) + zdz 
300               END IF
301            END DO
302         END DO
303      END IF
304
305      ! If ll_found = .false. then calculate MLD using difference of zdelta_T   
306      ! from the reference density/temperature
307 
308! Prevent this section from working on land points
309      WHERE ( tmask(:,:,1) /= 1.0 ) 
310         ll_found = .true. 
311      END WHERE
312 
313      DO jk=1, jpk 
314         ll_belowml(:,:,jk) = abs( zT(:,:,jk) - zT_ref(:,:) ) >= zdelta_T(:,:) 
315      END DO 
316 
317! Set default value where interpolation cannot be used (ll_found=false) 
318      DO jj = 1, jpj 
319         DO ji = 1, jpi 
320            IF ( .not. ll_found(ji,jj) )  hmld_zint(ji,jj) = fsdept(ji,jj,ikmt(ji,jj)) 
321         END DO
322      END DO
323
324      DO jj = 1, jpj 
325         DO ji = 1, jpi 
326!CDIR NOVECTOR
327            DO jk = ik_ref(ji,jj)+1, ikmt(ji,jj) 
328               IF ( ll_found(ji,jj) ) EXIT
329               IF ( ll_belowml(ji,jj,jk) ) THEN               
330                  zT_b = zT_ref(ji,jj) + zdelta_T(ji,jj) * SIGN(1.0, zdTdz(ji,jj,jk-1) ) 
331                  zdT  = zT_b - zT(ji,jj,jk-1)                                     
332                  zdz  = zdT / zdTdz(ji,jj,jk-1)                                       
333                  hmld_zint(ji,jj) = fsdept(ji,jj,jk-1) + zdz 
334                  EXIT                                                   
335               END IF
336            END DO
337         END DO
338      END DO
339
340      hmld_zint(:,:) = hmld_zint(:,:)*tmask(:,:,1) 
341     
342      CALL wrk_dealloc( jpi, jpj, ikmt, ik_ref, ik_iso) 
343      CALL wrk_dealloc( jpi, jpj, ppzdep, zT_ref, zdelta_T, zRHO1, zRHO2 ) 
344      CALL wrk_dealloc( jpi,jpj, jpk, zT, zdTdz, zmoddT ) 
345      !
346   END SUBROUTINE zdf_mxl_zint_mld
347
348   SUBROUTINE zdf_mxl_zint_htc( kt )
349      !!----------------------------------------------------------------------
350      !!                  ***  ROUTINE zdf_mxl_zint_htc  ***
351      !!
352      !! ** Purpose :   
353      !!
354      !! ** Method  :   
355      !!----------------------------------------------------------------------
356
357      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
358
359      INTEGER :: ji, jj, jk
360      INTEGER :: ikmax
361      REAL(wp) :: zc, zcoef
362      !
363      INTEGER,  ALLOCATABLE, DIMENSION(:,:) ::   ilevel
364      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zthick_0, zthick
365
366      !!----------------------------------------------------------------------
367
368      IF( .NOT. ALLOCATED(ilevel) ) THEN
369         ALLOCATE( ilevel(jpi,jpj), zthick_0(jpi,jpj), &
370         &         zthick(jpi,jpj), STAT=ji )
371         IF( lk_mpp  )   CALL mpp_sum(ji)
372         IF( ji /= 0 )   CALL ctl_stop( 'STOP', 'zdf_mxl_zint_htc : unable to allocate arrays' )
373      ENDIF
374
375      ! Find last whole model T level above the MLD
376      ilevel(:,:)   = 0
377      zthick_0(:,:) = 0._wp
378
379      DO jk = 1, jpkm1 
380         DO jj = 1, jpj
381            DO ji = 1, jpi                   
382               zthick_0(ji,jj) = zthick_0(ji,jj) + fse3t(ji,jj,jk)
383               IF( zthick_0(ji,jj) < hmld_zint(ji,jj) )   ilevel(ji,jj) = jk
384            END DO
385         END DO
386         WRITE(numout,*) 'zthick_0(jk =',jk,') =',zthick_0(2,2)
387         WRITE(numout,*) 'fsdepw(jk+1 =',jk+1,') =',fsdepw(2,2,jk+1)
388      END DO
389
390      ! Surface boundary condition
391      IF( lk_vvl ) THEN   ;   zthick(:,:) = 0._wp       ;   htc_mld(:,:) = 0._wp                                   
392      ELSE                ;   zthick(:,:) = sshn(:,:)   ;   htc_mld(:,:) = tsn(:,:,1,jp_tem) * sshn(:,:) * tmask(:,:,1)   
393      ENDIF
394
395      ! Deepest whole T level above the MLD
396      ikmax = MIN( MAXVAL( ilevel(:,:) ), jpkm1 )
397
398      ! Integration down to last whole model T level
399      DO jk = 1, ikmax
400         DO jj = 1, jpj
401            DO ji = 1, jpi
402               zc = fse3t(ji,jj,jk) * REAL( MIN( MAX( 0, ilevel(ji,jj) - jk + 1 ) , 1  )  )    ! 0 below ilevel
403               zthick(ji,jj) = zthick(ji,jj) + zc
404               htc_mld(ji,jj) = htc_mld(ji,jj) + zc * tsn(ji,jj,jk,jp_tem) * tmask(ji,jj,jk)
405            END DO
406         END DO
407      END DO
408
409      ! Subsequent partial T level
410      zthick(:,:) = hmld_zint(:,:) - zthick(:,:)   !   remaining thickness to reach MLD
411
412      DO jj = 1, jpj
413         DO ji = 1, jpi
414            htc_mld(ji,jj) = htc_mld(ji,jj) + tsn(ji,jj,ilevel(ji,jj)+1,jp_tem)  & 
415      &                      * MIN( fse3t(ji,jj,ilevel(ji,jj)+1), zthick(ji,jj) ) * tmask(ji,jj,ilevel(ji,jj)+1)
416         END DO
417      END DO
418
419      WRITE(numout,*) 'htc_mld(after) =',htc_mld(2,2)
420
421      ! Convert to heat content
422      zcoef = rau0 * rcp
423      htc_mld(:,:) = zcoef * htc_mld(:,:)
424
425   END SUBROUTINE zdf_mxl_zint_htc
426
427   SUBROUTINE zdf_mxl_zint( kt )
428      !!----------------------------------------------------------------------
429      !!                  ***  ROUTINE zdf_mxl_zint  ***
430      !!
431      !! ** Purpose :   
432      !!
433      !! ** Method  :   
434      !!----------------------------------------------------------------------
435
436      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
437
438      INTEGER :: ios
439      INTEGER :: jn
440
441      INTEGER :: nn_mld_diag = 0    ! number of diagnostics
442
443      CHARACTER(len=1) :: cmld
444
445      TYPE(MXL_ZINT) :: sn_mld1, sn_mld2, sn_mld3, sn_mld4, sn_mld5
446      TYPE(MXL_ZINT), SAVE, DIMENSION(5) ::   mld_diags
447
448      NAMELIST/namzdf_mldzint/ nn_mld_diag, sn_mld1, sn_mld2, sn_mld3, sn_mld4, sn_mld5
449
450      !!----------------------------------------------------------------------
451     
452      IF( kt == nit000 ) THEN
453         REWIND( numnam_ref )              ! Namelist namzdf_mldzint in reference namelist
454         READ  ( numnam_ref, namzdf_mldzint, IOSTAT = ios, ERR = 901)
455901      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_mldzint in reference namelist', lwp )
456
457         REWIND( numnam_cfg )              ! Namelist namzdf_mldzint in configuration namelist
458         READ  ( numnam_cfg, namzdf_mldzint, IOSTAT = ios, ERR = 902 )
459902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_mldzint in configuration namelist', lwp )
460         IF(lwm) WRITE ( numond, namzdf_mldzint )
461
462         IF( nn_mld_diag > 5 )   CALL ctl_stop( 'STOP', 'zdf_mxl_ini: Specify no more than 5 MLD definitions' )
463
464         mld_diags(1) = sn_mld1
465         mld_diags(2) = sn_mld2
466         mld_diags(3) = sn_mld3
467         mld_diags(4) = sn_mld4
468         mld_diags(5) = sn_mld5
469
470         IF( nn_mld_diag > 0 ) THEN
471            WRITE(numout,*) '=============== Vertically-interpolated mixed layer ================'
472            WRITE(numout,*) '(Diagnostic number, nn_mld_type, rn_zref, rn_dT_crit, rn_iso_frac)'
473            DO jn = 1, nn_mld_diag
474               WRITE(numout,*) 'MLD criterion',jn,':'
475               WRITE(numout,*) '    nn_mld_type =', mld_diags(jn)%mld_type
476               WRITE(numout,*) '    rn_zref ='    , mld_diags(jn)%zref
477               WRITE(numout,*) '    rn_dT_crit =' , mld_diags(jn)%dT_crit
478               WRITE(numout,*) '    rn_iso_frac =', mld_diags(jn)%iso_frac
479            END DO
480            WRITE(numout,*) '===================================================================='
481         ENDIF
482      ENDIF
483
484      IF( nn_mld_diag > 0 ) THEN
485         DO jn = 1, nn_mld_diag
486            WRITE(cmld,'(I1)') jn
487            IF( iom_use( "mldzint_"//cmld ) .OR. iom_use( "mldhtc_"//cmld ) ) THEN
488               CALL zdf_mxl_zint_mld( mld_diags(jn) )
489
490               IF( iom_use( "mldzint_"//cmld ) ) THEN
491                  CALL iom_put( "mldzint_"//cmld, hmld_zint(:,:) )
492               ENDIF
493
494               IF( iom_use( "mldhtc_"//cmld ) )  THEN
495                  CALL zdf_mxl_zint_htc( kt )
496                  CALL iom_put( "mldhtc_"//cmld , htc_mld(:,:)   )
497               ENDIF
498            ENDIF
499         END DO
500      ENDIF
501
502   END SUBROUTINE zdf_mxl_zint
503
504   !!======================================================================
505END MODULE zdfmxl
Note: See TracBrowser for help on using the repository browser.