source: NEMO/branches/UKMO/NEMO_4.0.1_GO8_package/src/OCE/ZDF/zdfmxl.F90 @ 11717

Last change on this file since 11717 was 11717, checked in by davestorkey, 12 months ago

UKMO/NEMO_4.0.1_GO8_package: copy over changes from NEMO_4.0_GO8_package branch.

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