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.

zdfddm.F90 in branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/ZDF – NEMO

Context Navigation

source: branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfddm.F90 @ 4424

Last change on this file since 4424 was 4424, checked in by trackstand2, 10 years ago
Bug fixes - use jpkorig in allocates instead of jpk
Property svn:keywords set to `Id`
File size: 12.5 KB

Rev	Line
[3]	1	MODULE zdfddm
	2	!!======================================================================
	3	!! * MODULE zdfddm *
	4	!! Ocean physics : double diffusion mixing parameterization
	5	!!======================================================================
[1601]	6	!! History : OPA ! 2000-08 (G. Madec) double diffusive mixing
	7	!! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module
[2528]	8	!! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase
[1601]	9	!!----------------------------------------------------------------------
[3]	10	#if defined key_zdfddm \|\| defined key_esopa
	11	!!----------------------------------------------------------------------
	12	!! 'key_zdfddm' : double diffusion
	13	!!----------------------------------------------------------------------
	14	!! zdf_ddm : compute the Ks for salinity
	15	!! zdf_ddm_init : read namelist and control the parameters
	16	!!----------------------------------------------------------------------
	17	USE oce ! ocean dynamics and tracers variables
	18	USE dom_oce ! ocean space and time domain variables
	19	USE zdf_oce ! ocean vertical physics variables
	20	USE in_out_manager ! I/O manager
	21	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
[258]	22	USE prtctl ! Print control
[2715]	23	USE lib_mpp ! MPP library
[3]	24
	25	IMPLICIT NONE
	26	PRIVATE
	27
[2528]	28	PUBLIC zdf_ddm ! called by step.F90
	29	PUBLIC zdf_ddm_init ! called by opa.F90
[2715]	30	PUBLIC zdf_ddm_alloc ! called by nemogcm.F90
[3]	31
[1537]	32	LOGICAL , PUBLIC, PARAMETER :: lk_zdfddm = .TRUE. !: double diffusive mixing flag
[3]	33
[2715]	34	REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: avs !: salinity vertical diffusivity coeff. at w-point
	35	REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: rrau !: heat/salt buoyancy flux ratio
[3]	36
[1601]	37	! !!* Namelist namzdf_ddm : double diffusive mixing *
[1537]	38	REAL(wp) :: rn_avts = 1.e-4_wp ! maximum value of avs for salt fingering
	39	REAL(wp) :: rn_hsbfr = 1.6_wp ! heat/salt buoyancy flux ratio
	40
[3211]	41	!! * Control permutation of array indices
	42	# include "zdfddm_ftrans.h90"
	43	# include "oce_ftrans.h90"
	44	# include "dom_oce_ftrans.h90"
	45	# include "zdf_oce_ftrans.h90"
	46
[3]	47	!! * Substitutions
	48	# include "vectopt_loop_substitute.h90"
	49	!!----------------------------------------------------------------------
[2715]	50	!! NEMO/OPA 4.0 , NEMO Consortium (2011)
[2528]	51	!! $Id$
[2715]	52	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[3]	53	!!----------------------------------------------------------------------
	54	CONTAINS
	55
[2715]	56	INTEGER FUNCTION zdf_ddm_alloc()
	57	!!----------------------------------------------------------------------
	58	!! * ROUTINE zdf_ddm_alloc *
	59	!!----------------------------------------------------------------------
[4424]	60	ALLOCATE( avs(jpi,jpj,jpkorig), rrau(jpi,jpj,jpkorig), &
	61	STAT= zdf_ddm_alloc )
[2715]	62	!
	63	IF( lk_mpp ) CALL mpp_sum ( zdf_ddm_alloc )
	64	IF( zdf_ddm_alloc /= 0 ) CALL ctl_warn('zdf_ddm_alloc: failed to allocate arrays')
	65	END FUNCTION zdf_ddm_alloc
	66
	67
[3]	68	SUBROUTINE zdf_ddm( kt )
	69	!!----------------------------------------------------------------------
	70	!! * ROUTINE zdf_ddm *
	71	!!
	72	!! ** Purpose : Add to the vertical eddy diffusivity coefficient the
[1601]	73	!! effect of salt fingering and diffusive convection.
[3]	74	!!
	75	!! ** Method : Diapycnal mixing is increased in case of double
	76	!! diffusive mixing (i.e. salt fingering and diffusive layering)
	77	!! following Merryfield et al. (1999). The rate of double diffusive
	78	!! mixing depend on the buoyancy ratio: Rrau=alpha/beta dk[T]/dk[S]
	79	!! which is computed in rn2.F
	80	!! * salt fingering (Schmitt 1981):
[1537]	81	!! for Rrau > 1 and rn2 > 0 : zavfs = rn_avts / ( 1 + (Rrau/rn_hsbfr)^6 )
[3]	82	!! for Rrau > 1 and rn2 > 0 : zavfs = O
	83	!! otherwise : zavft = 0.7 zavs / Rrau
	84	!! * diffusive layering (Federov 1988):
	85	!! for 0< Rrau < 1 and rn2 > 0 : zavdt = 1.3635e-6
	86	!! * exp( 4.6 exp(-0.54 (1/Rrau-1) ) )
	87	!! otherwise : zavdt = 0
	88	!! for .5 < Rrau < 1 and rn2 > 0 : zavds = zavdt (1.885 Rrau -0.85)
	89	!! for 0 < Rrau <.5 and rn2 > 0 : zavds = zavdt 0.15 Rrau
	90	!! otherwise : zavds = 0
	91	!! * update the eddy diffusivity:
	92	!! avt = avt + zavft + zavdt
	93	!! avs = avs + zavfs + zavds
	94	!! avmu, avmv are required to remain at least above avt and avs.
	95	!!
[1601]	96	!! ** Action : avt, avs : updated vertical eddy diffusivity coef. for T & S
[3]	97	!!
[1601]	98	!! References : Merryfield et al., JPO, 29, 1124-1142, 1999.
[3]	99	!!----------------------------------------------------------------------
[2715]	100	USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released
	101	USE wrk_nemo, ONLY: zmsks => wrk_2d_1 , zmskf => wrk_2d_2 , zmskd1 => wrk_2d_3 ! 2D workspace
	102	USE wrk_nemo, ONLY: zmskd2 => wrk_2d_4 , zmskd3 => wrk_2d_5 ! - -
	103	!
[1601]	104	INTEGER, INTENT(in) :: kt ! ocean time-step indexocean time step
[2715]	105	!
[1601]	106	INTEGER :: ji, jj , jk ! dummy loop indices
	107	REAL(wp) :: zinr, zrr ! temporary scalars
	108	REAL(wp) :: zavft, zavfs ! - -
	109	REAL(wp) :: zavdt, zavds ! - -
[3]	110	!!----------------------------------------------------------------------
	111
[2715]	112	IF( wrk_in_use(2, 1,2,3,4,5) ) THEN
	113	CALL ctl_stop('zdf_ddm: Requested workspace arrays already in use') ; RETURN
	114	ENDIF
	115
[3]	116	! ! ===============
	117	DO jk = 2, jpkm1 ! Horizontal slab
	118	! ! ===============
	119	! Define the mask
	120	! ---------------
[1601]	121	rrau(:,:,jk) = MAX( 1.e-20, rrau(:,:,jk) ) ! only retains positive value of rrau
[3]	122
[1601]	123	DO jj = 1, jpj ! indicators:
[3]	124	DO ji = 1, jpi
	125	! stability indicator: msks=1 if rn2>0; 0 elsewhere
[2715]	126	IF( rn2(ji,jj,jk) + 1.e-12 <= 0. ) THEN ; zmsks(ji,jj) = 0._wp
	127	ELSE ; zmsks(ji,jj) = 1._wp
[3]	128	ENDIF
	129	! salt fingering indicator: msksf=1 if rrau>1; 0 elsewhere
[2715]	130	IF( rrau(ji,jj,jk) <= 1. ) THEN ; zmskf(ji,jj) = 0._wp
	131	ELSE ; zmskf(ji,jj) = 1._wp
[3]	132	ENDIF
	133	! diffusive layering indicators:
[2528]	134	! ! mskdl1=1 if 0<rrau<1; 0 elsewhere
[2715]	135	IF( rrau(ji,jj,jk) >= 1. ) THEN ; zmskd1(ji,jj) = 0._wp
	136	ELSE ; zmskd1(ji,jj) = 1._wp
[3]	137	ENDIF
[2528]	138	! ! mskdl2=1 if 0<rrau<0.5; 0 elsewhere
[2715]	139	IF( rrau(ji,jj,jk) >= 0.5 ) THEN ; zmskd2(ji,jj) = 0._wp
	140	ELSE ; zmskd2(ji,jj) = 1._wp
[3]	141	ENDIF
	142	! mskdl3=1 if 0.5<rrau<1; 0 elsewhere
[2715]	143	IF( rrau(ji,jj,jk) <= 0.5 .OR. rrau(ji,jj,jk) >= 1. ) THEN ; zmskd3(ji,jj) = 0._wp
	144	ELSE ; zmskd3(ji,jj) = 1._wp
[3]	145	ENDIF
	146	END DO
	147	END DO
	148	! mask zmsk in order to have avt and avs masked
	149	zmsks(:,:) = zmsks(:,:) * tmask(:,:,jk)
	150
	151
	152	! Update avt and avs
	153	! ------------------
	154	! Constant eddy coefficient: reset to the background value
	155	!CDIR NOVERRCHK
	156	DO jj = 1, jpj
	157	!CDIR NOVERRCHK
	158	DO ji = 1, jpi
	159	zinr = 1./rrau(ji,jj,jk)
	160	! salt fingering
[1537]	161	zrr = rrau(ji,jj,jk)/rn_hsbfr
[3]	162	zrr = zrr * zrr
[1601]	163	zavfs = rn_avts / ( 1 + zrrzrrzrr ) * zmsks(ji,jj) * zmskf(ji,jj)
[1163]	164	zavft = 0.7 * zavfs * zinr
[3]	165	! diffusive layering
[1601]	166	zavdt = 1.3635e-6 * EXP( 4.6 * EXP( -0.54(zinr-1.) ) ) zmsks(ji,jj) * zmskd1(ji,jj)
	167	zavds = zavdt * zmsks(ji,jj) * ( (1.85 * rrau(ji,jj,jk) - 0.85 ) * zmskd3(ji,jj) &
	168	& + 0.15 * rrau(ji,jj,jk) * zmskd2(ji,jj) )
[3]	169	! add to the eddy viscosity coef. previously computed
	170	avs (ji,jj,jk) = avt(ji,jj,jk) + zavfs + zavds
	171	avt (ji,jj,jk) = avt(ji,jj,jk) + zavft + zavdt
[1527]	172	avm (ji,jj,jk) = avm(ji,jj,jk) + MAX( zavft + zavdt, zavfs + zavds )
[3]	173	END DO
	174	END DO
	175
	176
	177	! Increase avmu, avmv if necessary
	178	! --------------------------------
[1527]	179	!!gm to be changed following the definition of avm.
[3]	180	DO jj = 1, jpjm1
	181	DO ji = 1, fs_jpim1 ! vector opt.
	182	avmu(ji,jj,jk) = MAX( avmu(ji,jj,jk), &
[1601]	183	& avt(ji,jj,jk), avt(ji+1,jj,jk), &
	184	& avs(ji,jj,jk), avs(ji+1,jj,jk) ) * umask(ji,jj,jk)
[3]	185	avmv(ji,jj,jk) = MAX( avmv(ji,jj,jk), &
[1601]	186	& avt(ji,jj,jk), avt(ji,jj+1,jk), &
	187	& avs(ji,jj,jk), avs(ji,jj+1,jk) ) * vmask(ji,jj,jk)
[3]	188	END DO
	189	END DO
	190	! ! ===============
	191	END DO ! End of slab
	192	! ! ===============
[1601]	193	!
[2715]	194	CALL lbc_lnk( avt , 'W', 1._wp ) ! Lateral boundary conditions (unchanged sign)
	195	CALL lbc_lnk( avs , 'W', 1._wp )
	196	CALL lbc_lnk( avm , 'W', 1._wp )
	197	CALL lbc_lnk( avmu, 'U', 1._wp )
	198	CALL lbc_lnk( avmv, 'V', 1._wp )
[49]	199
[258]	200	IF(ln_ctl) THEN
	201	CALL prt_ctl(tab3d_1=avt , clinfo1=' ddm - t: ', tab3d_2=avs , clinfo2=' s: ', ovlap=1, kdim=jpk)
[516]	202	CALL prt_ctl(tab3d_1=avmu, clinfo1=' ddm - u: ', mask1=umask, &
	203	& tab3d_2=avmv, clinfo2= ' v: ', mask2=vmask, ovlap=1, kdim=jpk)
[49]	204	ENDIF
[1601]	205	!
[2715]	206	IF( wrk_not_released(2, 1,2,3,4,5) ) CALL ctl_stop('zdf_ddm: Release of workspace arrays failed')
	207	!
[3]	208	END SUBROUTINE zdf_ddm
	209
	210
	211	SUBROUTINE zdf_ddm_init
	212	!!----------------------------------------------------------------------
	213	!! * ROUTINE zdf_ddm_init *
	214	!!
	215	!! ** Purpose : Initialization of double diffusion mixing scheme
	216	!!
[1601]	217	!! ** Method : Read the namzdf_ddm namelist and check the parameter values
[1537]	218	!! called by zdf_ddm at the first timestep (nit000)
[3]	219	!!----------------------------------------------------------------------
[1601]	220	NAMELIST/namzdf_ddm/ rn_avts, rn_hsbfr
[3]	221	!!----------------------------------------------------------------------
[1537]	222	!
[2715]	223	REWIND( numnam ) ! Read Namelist namzdf_ddm : double diffusion mixing scheme
	224	READ ( numnam, namzdf_ddm )
[1537]	225	!
	226	IF(lwp) THEN ! Parameter print
[3]	227	WRITE(numout,*)
	228	WRITE(numout,*) 'zdf_ddm : double diffusive mixing'
	229	WRITE(numout,*) '~~~~~~~'
[1601]	230	WRITE(numout,*) ' Namelist namzdf_ddm : set dd mixing parameter'
[1537]	231	WRITE(numout,*) ' maximum avs for dd mixing rn_avts = ', rn_avts
	232	WRITE(numout,*) ' heat/salt buoyancy flux ratio rn_hsbfr = ', rn_hsbfr
[3]	233	ENDIF
[1537]	234	!
[2715]	235	! ! allocate zdfddm arrays
	236	IF( zdf_ddm_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'zdf_ddm_init : unable to allocate arrays' )
	237	!
[3837]	238	! Initialise new array
	239	avs(:,:,:) = 0.0_wp
	240
[3]	241	END SUBROUTINE zdf_ddm_init
	242
	243	#else
	244	!!----------------------------------------------------------------------
	245	!! Default option : Dummy module No double diffusion
	246	!!----------------------------------------------------------------------
[16]	247	LOGICAL, PUBLIC, PARAMETER :: lk_zdfddm = .FALSE. !: double diffusion flag
[3]	248	CONTAINS
	249	SUBROUTINE zdf_ddm( kt ) ! Dummy routine
[16]	250	WRITE(,) 'zdf_ddm: You should not have seen this print! error?', kt
[3]	251	END SUBROUTINE zdf_ddm
[2528]	252	SUBROUTINE zdf_ddm_init ! Dummy routine
	253	END SUBROUTINE zdf_ddm_init
[3]	254	#endif
	255
	256	!!======================================================================
	257	END MODULE zdfddm

Note: See TracBrowser for help on using the repository browser.

Download in other formats: