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/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM/NEMO/OPA_SRC/ZDF – NEMO

source: branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfmxl.F90 @ 9497

Last change on this file since 9497 was 9497, checked in by davestorkey, 4 years ago

branches/UKMO/dev_merge_2017_restart_datestamp_GO6_mixing : recommit science changes.

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