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 @ 3837

Last change on this file since 3837 was 3837, checked in by trackstand2, 11 years ago
Merge of finiss
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	!!----------------------------------------------------------------------
	60	ALLOCATE( avs(jpi,jpj,jpk), rrau(jpi,jpj,jpk), STAT= zdf_ddm_alloc )
	61	!
	62	IF( lk_mpp ) CALL mpp_sum ( zdf_ddm_alloc )
	63	IF( zdf_ddm_alloc /= 0 ) CALL ctl_warn('zdf_ddm_alloc: failed to allocate arrays')
	64	END FUNCTION zdf_ddm_alloc
	65
	66
[3]	67	SUBROUTINE zdf_ddm( kt )
	68	!!----------------------------------------------------------------------
	69	!! * ROUTINE zdf_ddm *
	70	!!
	71	!! ** Purpose : Add to the vertical eddy diffusivity coefficient the
[1601]	72	!! effect of salt fingering and diffusive convection.
[3]	73	!!
	74	!! ** Method : Diapycnal mixing is increased in case of double
	75	!! diffusive mixing (i.e. salt fingering and diffusive layering)
	76	!! following Merryfield et al. (1999). The rate of double diffusive
	77	!! mixing depend on the buoyancy ratio: Rrau=alpha/beta dk[T]/dk[S]
	78	!! which is computed in rn2.F
	79	!! * salt fingering (Schmitt 1981):
[1537]	80	!! for Rrau > 1 and rn2 > 0 : zavfs = rn_avts / ( 1 + (Rrau/rn_hsbfr)^6 )
[3]	81	!! for Rrau > 1 and rn2 > 0 : zavfs = O
	82	!! otherwise : zavft = 0.7 zavs / Rrau
	83	!! * diffusive layering (Federov 1988):
	84	!! for 0< Rrau < 1 and rn2 > 0 : zavdt = 1.3635e-6
	85	!! * exp( 4.6 exp(-0.54 (1/Rrau-1) ) )
	86	!! otherwise : zavdt = 0
	87	!! for .5 < Rrau < 1 and rn2 > 0 : zavds = zavdt (1.885 Rrau -0.85)
	88	!! for 0 < Rrau <.5 and rn2 > 0 : zavds = zavdt 0.15 Rrau
	89	!! otherwise : zavds = 0
	90	!! * update the eddy diffusivity:
	91	!! avt = avt + zavft + zavdt
	92	!! avs = avs + zavfs + zavds
	93	!! avmu, avmv are required to remain at least above avt and avs.
	94	!!
[1601]	95	!! ** Action : avt, avs : updated vertical eddy diffusivity coef. for T & S
[3]	96	!!
[1601]	97	!! References : Merryfield et al., JPO, 29, 1124-1142, 1999.
[3]	98	!!----------------------------------------------------------------------
[2715]	99	USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released
	100	USE wrk_nemo, ONLY: zmsks => wrk_2d_1 , zmskf => wrk_2d_2 , zmskd1 => wrk_2d_3 ! 2D workspace
	101	USE wrk_nemo, ONLY: zmskd2 => wrk_2d_4 , zmskd3 => wrk_2d_5 ! - -
	102	!
[1601]	103	INTEGER, INTENT(in) :: kt ! ocean time-step indexocean time step
[2715]	104	!
[1601]	105	INTEGER :: ji, jj , jk ! dummy loop indices
	106	REAL(wp) :: zinr, zrr ! temporary scalars
	107	REAL(wp) :: zavft, zavfs ! - -
	108	REAL(wp) :: zavdt, zavds ! - -
[3]	109	!!----------------------------------------------------------------------
	110
[2715]	111	IF( wrk_in_use(2, 1,2,3,4,5) ) THEN
	112	CALL ctl_stop('zdf_ddm: Requested workspace arrays already in use') ; RETURN
	113	ENDIF
	114
[3]	115	! ! ===============
	116	DO jk = 2, jpkm1 ! Horizontal slab
	117	! ! ===============
	118	! Define the mask
	119	! ---------------
[1601]	120	rrau(:,:,jk) = MAX( 1.e-20, rrau(:,:,jk) ) ! only retains positive value of rrau
[3]	121
[1601]	122	DO jj = 1, jpj ! indicators:
[3]	123	DO ji = 1, jpi
	124	! stability indicator: msks=1 if rn2>0; 0 elsewhere
[2715]	125	IF( rn2(ji,jj,jk) + 1.e-12 <= 0. ) THEN ; zmsks(ji,jj) = 0._wp
	126	ELSE ; zmsks(ji,jj) = 1._wp
[3]	127	ENDIF
	128	! salt fingering indicator: msksf=1 if rrau>1; 0 elsewhere
[2715]	129	IF( rrau(ji,jj,jk) <= 1. ) THEN ; zmskf(ji,jj) = 0._wp
	130	ELSE ; zmskf(ji,jj) = 1._wp
[3]	131	ENDIF
	132	! diffusive layering indicators:
[2528]	133	! ! mskdl1=1 if 0<rrau<1; 0 elsewhere
[2715]	134	IF( rrau(ji,jj,jk) >= 1. ) THEN ; zmskd1(ji,jj) = 0._wp
	135	ELSE ; zmskd1(ji,jj) = 1._wp
[3]	136	ENDIF
[2528]	137	! ! mskdl2=1 if 0<rrau<0.5; 0 elsewhere
[2715]	138	IF( rrau(ji,jj,jk) >= 0.5 ) THEN ; zmskd2(ji,jj) = 0._wp
	139	ELSE ; zmskd2(ji,jj) = 1._wp
[3]	140	ENDIF
	141	! mskdl3=1 if 0.5<rrau<1; 0 elsewhere
[2715]	142	IF( rrau(ji,jj,jk) <= 0.5 .OR. rrau(ji,jj,jk) >= 1. ) THEN ; zmskd3(ji,jj) = 0._wp
	143	ELSE ; zmskd3(ji,jj) = 1._wp
[3]	144	ENDIF
	145	END DO
	146	END DO
	147	! mask zmsk in order to have avt and avs masked
	148	zmsks(:,:) = zmsks(:,:) * tmask(:,:,jk)
	149
	150
	151	! Update avt and avs
	152	! ------------------
	153	! Constant eddy coefficient: reset to the background value
	154	!CDIR NOVERRCHK
	155	DO jj = 1, jpj
	156	!CDIR NOVERRCHK
	157	DO ji = 1, jpi
	158	zinr = 1./rrau(ji,jj,jk)
	159	! salt fingering
[1537]	160	zrr = rrau(ji,jj,jk)/rn_hsbfr
[3]	161	zrr = zrr * zrr
[1601]	162	zavfs = rn_avts / ( 1 + zrrzrrzrr ) * zmsks(ji,jj) * zmskf(ji,jj)
[1163]	163	zavft = 0.7 * zavfs * zinr
[3]	164	! diffusive layering
[1601]	165	zavdt = 1.3635e-6 * EXP( 4.6 * EXP( -0.54(zinr-1.) ) ) zmsks(ji,jj) * zmskd1(ji,jj)
	166	zavds = zavdt * zmsks(ji,jj) * ( (1.85 * rrau(ji,jj,jk) - 0.85 ) * zmskd3(ji,jj) &
	167	& + 0.15 * rrau(ji,jj,jk) * zmskd2(ji,jj) )
[3]	168	! add to the eddy viscosity coef. previously computed
	169	avs (ji,jj,jk) = avt(ji,jj,jk) + zavfs + zavds
	170	avt (ji,jj,jk) = avt(ji,jj,jk) + zavft + zavdt
[1527]	171	avm (ji,jj,jk) = avm(ji,jj,jk) + MAX( zavft + zavdt, zavfs + zavds )
[3]	172	END DO
	173	END DO
	174
	175
	176	! Increase avmu, avmv if necessary
	177	! --------------------------------
[1527]	178	!!gm to be changed following the definition of avm.
[3]	179	DO jj = 1, jpjm1
	180	DO ji = 1, fs_jpim1 ! vector opt.
	181	avmu(ji,jj,jk) = MAX( avmu(ji,jj,jk), &
[1601]	182	& avt(ji,jj,jk), avt(ji+1,jj,jk), &
	183	& avs(ji,jj,jk), avs(ji+1,jj,jk) ) * umask(ji,jj,jk)
[3]	184	avmv(ji,jj,jk) = MAX( avmv(ji,jj,jk), &
[1601]	185	& avt(ji,jj,jk), avt(ji,jj+1,jk), &
	186	& avs(ji,jj,jk), avs(ji,jj+1,jk) ) * vmask(ji,jj,jk)
[3]	187	END DO
	188	END DO
	189	! ! ===============
	190	END DO ! End of slab
	191	! ! ===============
[1601]	192	!
[2715]	193	CALL lbc_lnk( avt , 'W', 1._wp ) ! Lateral boundary conditions (unchanged sign)
	194	CALL lbc_lnk( avs , 'W', 1._wp )
	195	CALL lbc_lnk( avm , 'W', 1._wp )
	196	CALL lbc_lnk( avmu, 'U', 1._wp )
	197	CALL lbc_lnk( avmv, 'V', 1._wp )
[49]	198
[258]	199	IF(ln_ctl) THEN
	200	CALL prt_ctl(tab3d_1=avt , clinfo1=' ddm - t: ', tab3d_2=avs , clinfo2=' s: ', ovlap=1, kdim=jpk)
[516]	201	CALL prt_ctl(tab3d_1=avmu, clinfo1=' ddm - u: ', mask1=umask, &
	202	& tab3d_2=avmv, clinfo2= ' v: ', mask2=vmask, ovlap=1, kdim=jpk)
[49]	203	ENDIF
[1601]	204	!
[2715]	205	IF( wrk_not_released(2, 1,2,3,4,5) ) CALL ctl_stop('zdf_ddm: Release of workspace arrays failed')
	206	!
[3]	207	END SUBROUTINE zdf_ddm
	208
	209
	210	SUBROUTINE zdf_ddm_init
	211	!!----------------------------------------------------------------------
	212	!! * ROUTINE zdf_ddm_init *
	213	!!
	214	!! ** Purpose : Initialization of double diffusion mixing scheme
	215	!!
[1601]	216	!! ** Method : Read the namzdf_ddm namelist and check the parameter values
[1537]	217	!! called by zdf_ddm at the first timestep (nit000)
[3]	218	!!----------------------------------------------------------------------
[1601]	219	NAMELIST/namzdf_ddm/ rn_avts, rn_hsbfr
[3]	220	!!----------------------------------------------------------------------
[1537]	221	!
[2715]	222	REWIND( numnam ) ! Read Namelist namzdf_ddm : double diffusion mixing scheme
	223	READ ( numnam, namzdf_ddm )
[1537]	224	!
	225	IF(lwp) THEN ! Parameter print
[3]	226	WRITE(numout,*)
	227	WRITE(numout,*) 'zdf_ddm : double diffusive mixing'
	228	WRITE(numout,*) '~~~~~~~'
[1601]	229	WRITE(numout,*) ' Namelist namzdf_ddm : set dd mixing parameter'
[1537]	230	WRITE(numout,*) ' maximum avs for dd mixing rn_avts = ', rn_avts
	231	WRITE(numout,*) ' heat/salt buoyancy flux ratio rn_hsbfr = ', rn_hsbfr
[3]	232	ENDIF
[1537]	233	!
[2715]	234	! ! allocate zdfddm arrays
	235	IF( zdf_ddm_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'zdf_ddm_init : unable to allocate arrays' )
	236	!
[3837]	237	! Initialise new array
	238	avs(:,:,:) = 0.0_wp
	239
[3]	240	END SUBROUTINE zdf_ddm_init
	241
	242	#else
	243	!!----------------------------------------------------------------------
	244	!! Default option : Dummy module No double diffusion
	245	!!----------------------------------------------------------------------
[16]	246	LOGICAL, PUBLIC, PARAMETER :: lk_zdfddm = .FALSE. !: double diffusion flag
[3]	247	CONTAINS
	248	SUBROUTINE zdf_ddm( kt ) ! Dummy routine
[16]	249	WRITE(,) 'zdf_ddm: You should not have seen this print! error?', kt
[3]	250	END SUBROUTINE zdf_ddm
[2528]	251	SUBROUTINE zdf_ddm_init ! Dummy routine
	252	END SUBROUTINE zdf_ddm_init
[3]	253	#endif
	254
	255	!!======================================================================
	256	END MODULE zdfddm

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

Download in other formats: