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.

zdfbfr.F90 in branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/OPA_SRC/ZDF – NEMO

Context Navigation

source: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfbfr.F90 @ 4616

Last change on this file since 4616 was 4616, checked in by gm, 10 years ago
#1260 : see the associated wiki page for explanation
Property svn:keywords set to `Id`
File size: 15.9 KB

Rev	Line
[3]	1	MODULE zdfbfr
	2	!!======================================================================
	3	!! * MODULE zdfbfr *
	4	!! Ocean physics: Bottom friction
	5	!!======================================================================
[1662]	6	!! History : OPA ! 1997-06 (G. Madec, A.-M. Treguier) Original code
[1601]	7	!! NEMO 1.0 ! 2002-06 (G. Madec) F90: Free form and module
[1708]	8	!! 3.2 ! 2009-09 (A.C.Coward) Correction to include barotropic contribution
[2528]	9	!! 3.3 ! 2010-10 (C. Ethe, G. Madec) reorganisation of initialisation phase
[3633]	10	!! 3.4 ! 2011-11 (H. Liu) implementation of semi-implicit bottom friction option
	11	!! ! 2012-06 (H. Liu) implementation of Log Layer bottom friction option
[1601]	12	!!----------------------------------------------------------------------
[3]	13
	14	!!----------------------------------------------------------------------
[4381]	15	!! zdf_bfr : update bottom friction coefficient (non-linear bottom friction only)
[1601]	16	!! zdf_bfr_init : read in namelist and control the bottom friction parameters.
[3]	17	!!----------------------------------------------------------------------
	18	USE oce ! ocean dynamics and tracers variables
[3764]	19	USE dom_oce ! ocean space and time domain variables
[3]	20	USE zdf_oce ! ocean vertical physics variables
	21	USE in_out_manager ! I/O manager
	22	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
[1671]	23	USE lib_mpp ! distributed memory computing
[258]	24	USE prtctl ! Print control
[3294]	25	USE timing ! Timing
[4381]	26	USE wrk_nemo ! Memory Allocation
[3633]	27	USE phycst, ONLY: vkarmn
	28
[3]	29	IMPLICIT NONE
	30	PRIVATE
	31
[2528]	32	PUBLIC zdf_bfr ! called by step.F90
[4381]	33	PUBLIC zdf_bfr_init ! called by nemogcm.F90
[2715]	34
[4147]	35	! !!* Namelist nambfr: bottom friction namelist *
	36	INTEGER , PUBLIC :: nn_bfr ! = 0/1/2/3 type of bottom friction (PUBLIC for TAM)
	37	REAL(wp), PUBLIC :: rn_bfri1 ! bottom drag coefficient (linear case) (PUBLIC for TAM)
	38	REAL(wp), PUBLIC :: rn_bfri2 ! bottom drag coefficient (non linear case) (PUBLIC for TAM)
[4381]	39	REAL(wp), PUBLIC :: rn_bfri2_max ! Maximum bottom drag coefficient (non linear case and ln_loglayer=T) (PUBLIC for TAM)
[4147]	40	REAL(wp), PUBLIC :: rn_bfeb2 ! background bottom turbulent kinetic energy [m2/s2] (PUBLIC for TAM)
	41	REAL(wp), PUBLIC :: rn_bfrien ! local factor to enhance coefficient bfri (PUBLIC for TAM)
	42	LOGICAL , PUBLIC :: ln_bfr2d ! logical switch for 2D enhancement (PUBLIC for TAM)
	43	LOGICAL , PUBLIC :: ln_loglayer ! switch for log layer bfr coeff. (PUBLIC for TAM)
	44	REAL(wp), PUBLIC :: rn_bfrz0 ! bottom roughness for loglayer bfr coeff (PUBLIC for TAM)
[4381]	45	LOGICAL , PUBLIC :: ln_bfrimp ! logical switch for implicit bottom friction
[3764]	46	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: bfrcoef2d ! 2D bottom drag coefficient (PUBLIC for TAM)
[3]	47
	48	!! * Substitutions
[2528]	49	# include "vectopt_loop_substitute.h90"
[3]	50	# include "domzgr_substitute.h90"
	51	!!----------------------------------------------------------------------
[2715]	52	!! NEMO/OPA 4.0 , NEMO Consortium (2011)
[2528]	53	!! $Id$
	54	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[3]	55	!!----------------------------------------------------------------------
	56	CONTAINS
	57
[2715]	58	INTEGER FUNCTION zdf_bfr_alloc()
	59	!!----------------------------------------------------------------------
	60	!! * FUNCTION zdf_bfr_alloc *
	61	!!----------------------------------------------------------------------
	62	ALLOCATE( bfrcoef2d(jpi,jpj), STAT=zdf_bfr_alloc )
	63	!
	64	IF( lk_mpp ) CALL mpp_sum ( zdf_bfr_alloc )
	65	IF( zdf_bfr_alloc /= 0 ) CALL ctl_warn('zdf_bfr_alloc: failed to allocate arrays.')
	66	END FUNCTION zdf_bfr_alloc
	67
	68
[3]	69	SUBROUTINE zdf_bfr( kt )
	70	!!----------------------------------------------------------------------
	71	!! * ROUTINE zdf_bfr *
[3764]	72	!!
[1708]	73	!! ** Purpose : compute the bottom friction coefficient.
[3]	74	!!
[3764]	75	!! ** Method : Calculate and store part of the momentum trend due
	76	!! to bottom friction following the chosen friction type
[1662]	77	!! (free-slip, linear, or quadratic). The component
	78	!! calculated here is multiplied by the bottom velocity in
	79	!! dyn_bfr to provide the trend term.
[1708]	80	!! The coefficients are updated at each time step only
	81	!! in the quadratic case.
[1662]	82	!!
	83	!! ** Action : bfrua , bfrva bottom friction coefficients
[3]	84	!!----------------------------------------------------------------------
	85	INTEGER, INTENT( in ) :: kt ! ocean time-step index
[1601]	86	!!
[3633]	87	INTEGER :: ji, jj ! dummy loop indices
[4381]	88	INTEGER :: ikbt, ikbu, ikbv ! local integers
	89	REAL(wp) :: zvu, zuv, zecu, zecv, ztmp ! temporary scalars
	90	REAL(wp), POINTER, DIMENSION(:,:) :: zbfrt
[3]	91	!!----------------------------------------------------------------------
[3294]	92	!
	93	IF( nn_timing == 1 ) CALL timing_start('zdf_bfr')
	94	!
[4381]	95	IF( kt == nit000 .AND. lwp ) THEN
	96	WRITE(numout,*)
	97	WRITE(numout,*) 'zdf_bfr : Set bottom friction coefficient (non-linear case)'
	98	WRITE(numout,*) '~~~~~~~~'
	99	ENDIF
	100	!
	101	IF( nn_bfr == 2 ) THEN ! quadratic bottom friction only
[1662]	102	!
[4381]	103	CALL wrk_alloc( jpi, jpj, zbfrt )
[3633]	104
[4381]	105	IF ( ln_loglayer.AND.lk_vvl ) THEN ! "log layer" bottom friction coefficient
[3820]	106
[4381]	107	DO jj = 1, jpj
	108	DO ji = 1, jpi
	109	ikbt = mbkt(ji,jj)
	110	! JC: possible WAD implementation should modify line below if layers vanish
	111	ztmp = tmask(ji,jj,ikbt) * ( vkarmn / LOG( 0.5_wp * fse3t_n(ji,jj,ikbt) / rn_bfrz0 ))**2._wp
	112	zbfrt(ji,jj) = MAX(bfrcoef2d(ji,jj), ztmp)
	113	zbfrt(ji,jj) = MIN(zbfrt(ji,jj), rn_bfri2_max)
	114	END DO
	115	END DO
	116	!
	117	ELSE
	118	zbfrt(:,:) = bfrcoef2d(:,:)
[3633]	119	ENDIF
	120
[3]	121	DO jj = 2, jpjm1
	122	DO ji = 2, jpim1
[3764]	123	ikbu = mbku(ji,jj) ! ocean bottom level at u- and v-points
[2528]	124	ikbv = mbkv(ji,jj) ! (deepest ocean u- and v-points)
[1662]	125	!
[2528]	126	zvu = 0.25 * ( vn(ji,jj ,ikbu) + vn(ji+1,jj ,ikbu) &
	127	& + vn(ji,jj-1,ikbu) + vn(ji+1,jj-1,ikbu) )
	128	zuv = 0.25 * ( un(ji,jj ,ikbv) + un(ji-1,jj ,ikbv) &
	129	& + un(ji,jj+1,ikbv) + un(ji-1,jj+1,ikbv) )
[1662]	130	!
[2528]	131	zecu = SQRT( un(ji,jj,ikbu) * un(ji,jj,ikbu) + zvu*zvu + rn_bfeb2 )
	132	zecv = SQRT( vn(ji,jj,ikbv) * vn(ji,jj,ikbv) + zuv*zuv + rn_bfeb2 )
[1662]	133	!
[4381]	134	bfrua(ji,jj) = - 0.5_wp * ( zbfrt(ji,jj) + zbfrt(ji+1,jj ) ) * zecu
	135	bfrva(ji,jj) = - 0.5_wp * ( zbfrt(ji,jj) + zbfrt(ji ,jj+1) ) * zecv
[3]	136	END DO
	137	END DO
[3633]	138
[1601]	139	!
[1708]	140	CALL lbc_lnk( bfrua, 'U', 1. ) ; CALL lbc_lnk( bfrva, 'V', 1. ) ! Lateral boundary condition
	141	!
	142	IF(ln_ctl) CALL prt_ctl( tab2d_1=bfrua, clinfo1=' bfr - u: ', mask1=umask, &
	143	& tab2d_2=bfrva, clinfo2= ' v: ', mask2=vmask,ovlap=1 )
[4381]	144	CALL wrk_dealloc( jpi,jpj,zbfrt )
[1662]	145	ENDIF
	146	!
[3294]	147	IF( nn_timing == 1 ) CALL timing_stop('zdf_bfr')
	148	!
[3]	149	END SUBROUTINE zdf_bfr
	150
	151
	152	SUBROUTINE zdf_bfr_init
	153	!!----------------------------------------------------------------------
	154	!! * ROUTINE zdf_bfr_init *
[3764]	155	!!
[3]	156	!! ** Purpose : Initialization of the bottom friction
	157	!!
[4381]	158	!! ** Method : Read the nambfr namelist and check their consistency
	159	!! called at the first timestep (nit000)
[3]	160	!!----------------------------------------------------------------------
[1662]	161	USE iom ! I/O module for ehanced bottom friction file
	162	!!
[3633]	163	INTEGER :: inum ! logical unit for enhanced bottom friction file
	164	INTEGER :: ji, jj ! dummy loop indexes
[4381]	165	INTEGER :: ikbt, ikbu, ikbv ! temporary integers
	166	INTEGER :: ictu, ictv ! - -
[4147]	167	INTEGER :: ios
[3633]	168	REAL(wp) :: zminbfr, zmaxbfr ! temporary scalars
[4381]	169	REAL(wp) :: ztmp, zfru, zfrv ! - -
[1708]	170	!!
[4381]	171	NAMELIST/nambfr/ nn_bfr, rn_bfri1, rn_bfri2, rn_bfri2_max, rn_bfeb2, rn_bfrz0, ln_bfr2d, &
[3633]	172	& rn_bfrien, ln_bfrimp, ln_loglayer
[3]	173	!!----------------------------------------------------------------------
[3294]	174	!
	175	IF( nn_timing == 1 ) CALL timing_start('zdf_bfr_init')
	176	!
[4381]	177	! !* Allocate zdfbfr arrays
	178	IF( zdf_bfr_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'zdf_bfr_init : unable to allocate arrays' )
	179	!
	180	! !* Parameter control and print
	181	!
[4147]	182	REWIND( numnam_ref ) ! Namelist nambfr in reference namelist : Bottom momentum boundary condition
	183	READ ( numnam_ref, nambfr, IOSTAT = ios, ERR = 901)
	184	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambfr in reference namelist', lwp )
[3]	185
[4147]	186	REWIND( numnam_cfg ) ! Namelist nambfr in configuration namelist : Bottom momentum boundary condition
	187	READ ( numnam_cfg, nambfr, IOSTAT = ios, ERR = 902 )
	188	902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambfr in configuration namelist', lwp )
	189	WRITE ( numond, nambfr )
[1662]	190	IF(lwp) WRITE(numout,*)
[4381]	191	IF(lwp) WRITE(numout,*) 'zdf_bfr_init : momentum bottom friction'
	192	IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~'
[1708]	193	IF(lwp) WRITE(numout,*) ' Namelist nam_bfr : set bottom friction parameters'
[4381]	194	!
[1662]	195	SELECT CASE (nn_bfr)
[2715]	196	!
[3]	197	CASE( 0 )
[1708]	198	IF(lwp) WRITE(numout,*) ' free-slip '
[1662]	199	bfrua(:,:) = 0.e0
	200	bfrva(:,:) = 0.e0
	201	!
[3]	202	CASE( 1 )
[1708]	203	IF(lwp) WRITE(numout,*) ' linear botton friction'
	204	IF(lwp) WRITE(numout,*) ' friction coef. rn_bfri1 = ', rn_bfri1
	205	IF( ln_bfr2d ) THEN
	206	IF(lwp) WRITE(numout,*) ' read coef. enhancement distribution from file ln_bfr2d = ', ln_bfr2d
	207	IF(lwp) WRITE(numout,*) ' coef rn_bfri2 enhancement factor rn_bfrien = ',rn_bfrien
[1662]	208	ENDIF
[1601]	209	!
[3764]	210	IF(ln_bfr2d) THEN
[1662]	211	! bfr_coef is a coefficient in [0,1] giving the mask where to apply the bfr enhancement
	212	CALL iom_open('bfr_coef.nc',inum)
	213	CALL iom_get (inum, jpdom_data, 'bfr_coef',bfrcoef2d,1) ! bfrcoef2d is used as tmp array
	214	CALL iom_close(inum)
[2715]	215	bfrcoef2d(:,:) = rn_bfri1 * ( 1 + rn_bfrien * bfrcoef2d(:,:) )
[4381]	216	ELSE
	217	bfrcoef2d(:,:) = rn_bfri1 ! initialize bfrcoef2d to the namelist variable
[1662]	218	ENDIF
[4381]	219	!
[1662]	220	bfrua(:,:) = - bfrcoef2d(:,:)
	221	bfrva(:,:) = - bfrcoef2d(:,:)
	222	!
[3]	223	CASE( 2 )
[4381]	224	IF(lwp) WRITE(numout,*) ' quadratic bottom friction'
[1708]	225	IF(lwp) WRITE(numout,*) ' friction coef. rn_bfri2 = ', rn_bfri2
[4381]	226	IF(lwp) WRITE(numout,*) ' Max. coef. (log case) rn_bfri2_max = ', rn_bfri2_max
[1708]	227	IF(lwp) WRITE(numout,*) ' background tke rn_bfeb2 = ', rn_bfeb2
[4381]	228	IF(lwp) WRITE(numout,*) ' log formulation ln_bfr2d = ', ln_loglayer
	229	IF(lwp) WRITE(numout,*) ' bottom roughness rn_bfrz0 [m] = ', rn_bfrz0
	230	IF( rn_bfrz0<=0.e0 ) THEN
	231	WRITE(ctmp1,*) ' bottom roughness must be strictly positive'
	232	CALL ctl_stop( ctmp1 )
	233	ENDIF
[1708]	234	IF( ln_bfr2d ) THEN
	235	IF(lwp) WRITE(numout,*) ' read coef. enhancement distribution from file ln_bfr2d = ', ln_bfr2d
	236	IF(lwp) WRITE(numout,*) ' coef rn_bfri2 enhancement factor rn_bfrien = ',rn_bfrien
[1662]	237	ENDIF
[1708]	238	!
[3764]	239	IF(ln_bfr2d) THEN
[1662]	240	! bfr_coef is a coefficient in [0,1] giving the mask where to apply the bfr enhancement
	241	CALL iom_open('bfr_coef.nc',inum)
	242	CALL iom_get (inum, jpdom_data, 'bfr_coef',bfrcoef2d,1) ! bfrcoef2d is used as tmp array
	243	CALL iom_close(inum)
[4381]	244	!
	245	bfrcoef2d(:,:) = rn_bfri2 * ( 1 + rn_bfrien * bfrcoef2d(:,:) )
	246	ELSE
	247	bfrcoef2d(:,:) = rn_bfri2 ! initialize bfrcoef2d to the namelist variable
[1662]	248	ENDIF
[1601]	249	!
[4381]	250	IF ( ln_loglayer.AND.(.NOT.lk_vvl) ) THEN ! set "log layer" bottom friction once for all
	251	DO jj = 1, jpj
	252	DO ji = 1, jpi
	253	ikbt = mbkt(ji,jj)
	254	ztmp = tmask(ji,jj,ikbt) * ( vkarmn / LOG( 0.5_wp * fse3t_n(ji,jj,ikbt) / rn_bfrz0 ))**2._wp
	255	bfrcoef2d(ji,jj) = MAX(bfrcoef2d(ji,jj), ztmp)
	256	bfrcoef2d(ji,jj) = MIN(bfrcoef2d(ji,jj), rn_bfri2_max)
	257	END DO
	258	END DO
	259	ENDIF
	260	!
[3]	261	CASE DEFAULT
[1708]	262	IF(lwp) WRITE(ctmp1,*) ' bad flag value for nn_bfr = ', nn_bfr
[474]	263	CALL ctl_stop( ctmp1 )
[1601]	264	!
[3]	265	END SELECT
[4381]	266	!
	267	IF(lwp) WRITE(numout,*) ' implicit bottom friction switch ln_bfrimp = ', ln_bfrimp
	268	!
	269	! ! Make sure ln_zdfexp=.false. when use implicit bfr
	270	IF( ln_bfrimp .AND. ln_zdfexp ) THEN
[3820]	271	IF(lwp) THEN
	272	WRITE(numout,*)
[4381]	273	WRITE(numout,*) 'Implicit bottom friction can only be used when ln_zdfexp=.false.'
	274	WRITE(numout,*) ' but you set: ln_bfrimp=.true. and ln_zdfexp=.true.!!!!'
	275	WRITE(ctmp1,*) ' set either ln_zdfexp = .false or ln_bfrimp = .false.'
[3820]	276	CALL ctl_stop( ctmp1 )
	277	END IF
	278	END IF
[1601]	279	!
[1662]	280	! Basic stability check on bottom friction coefficient
	281	!
	282	ictu = 0 ! counter for stability criterion breaches at U-pts
	283	ictv = 0 ! counter for stability criterion breaches at V-pts
[1708]	284	zminbfr = 1.e10_wp ! initialise tracker for minimum of bottom friction coefficient
	285	zmaxbfr = -1.e10_wp ! initialise tracker for maximum of bottom friction coefficient
[1662]	286	!
	287	DO jj = 2, jpjm1
	288	DO ji = 2, jpim1
[2528]	289	ikbu = mbku(ji,jj) ! deepest ocean level at u- and v-points
	290	ikbv = mbkv(ji,jj)
	291	zfru = 0.5 * fse3u(ji,jj,ikbu) / rdt
	292	zfrv = 0.5 * fse3v(ji,jj,ikbv) / rdt
[1708]	293	IF( ABS( bfrcoef2d(ji,jj) ) > zfru ) THEN
	294	IF( ln_ctl ) THEN
[2528]	295	WRITE(numout,*) 'BFR ', narea, nimpp+ji, njmpp+jj, ikbu
	296	WRITE(numout,*) 'BFR ', ABS( bfrcoef2d(ji,jj) ), zfru
[1662]	297	ENDIF
	298	ictu = ictu + 1
	299	ENDIF
[1708]	300	IF( ABS( bfrcoef2d(ji,jj) ) > zfrv ) THEN
	301	IF( ln_ctl ) THEN
	302	WRITE(numout,*) 'BFR ', narea, nimpp+ji, njmpp+jj, ikbv
	303	WRITE(numout,*) 'BFR ', bfrcoef2d(ji,jj), zfrv
[1662]	304	ENDIF
	305	ictv = ictv + 1
	306	ENDIF
[1708]	307	zminbfr = MIN( zminbfr, MIN( zfru, ABS( bfrcoef2d(ji,jj) ) ) )
	308	zmaxbfr = MAX( zmaxbfr, MIN( zfrv, ABS( bfrcoef2d(ji,jj) ) ) )
[1662]	309	END DO
	310	END DO
[1708]	311	IF( lk_mpp ) THEN
[1662]	312	CALL mpp_sum( ictu )
	313	CALL mpp_sum( ictv )
[1708]	314	CALL mpp_min( zminbfr )
	315	CALL mpp_max( zmaxbfr )
[1662]	316	ENDIF
[3294]	317	IF( .NOT.ln_bfrimp) THEN
[2528]	318	IF( lwp .AND. ictu + ictv > 0 ) THEN
[1662]	319	WRITE(numout,*) ' Bottom friction stability check failed at ', ictu, ' U-points '
	320	WRITE(numout,*) ' Bottom friction stability check failed at ', ictv, ' V-points '
[1708]	321	WRITE(numout,*) ' Bottom friction coefficient now ranges from: ', zminbfr, ' to ', zmaxbfr
	322	WRITE(numout,*) ' Bottom friction coefficient will be reduced where necessary'
[1662]	323	ENDIF
[3294]	324	ENDIF
[1662]	325	!
[3294]	326	IF( nn_timing == 1 ) CALL timing_stop('zdf_bfr_init')
	327	!
[3]	328	END SUBROUTINE zdf_bfr_init
	329
	330	!!======================================================================
	331	END MODULE zdfbfr

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

Download in other formats: