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.

trabbl.F90 in NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/TRA – NEMO

Context Navigation

source: NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/TRA/trabbl.F90 @ 10402

Last change on this file since 10402 was 10402, checked in by smasson, 5 years ago
dev_r10164_HPC09_ESIWACE_PREP_MERGE: no more need of lk_mpp for mpp_sum/max/min, see #2133
Property svn:keywords set to `Id`
File size: 31.9 KB

Rev	Line
[3]	1	MODULE trabbl
	2	!!==============================================================================
	3	!! * MODULE trabbl *
	4	!! Ocean physics : advective and/or diffusive bottom boundary layer scheme
	5	!!==============================================================================
[2528]	6	!! History : OPA ! 1996-06 (L. Mortier) Original code
	7	!! 8.0 ! 1997-11 (G. Madec) Optimization
	8	!! NEMO 1.0 ! 2002-08 (G. Madec) free form + modules
	9	!! - ! 2004-01 (A. de Miranda, G. Madec, J.M. Molines ) add advective bbl
[3764]	10	!! 3.3 ! 2009-11 (G. Madec) merge trabbl and trabbl_adv + style + optimization
	11	!! - ! 2010-04 (G. Madec) Campin & Goosse advective bbl
[2528]	12	!! - ! 2010-06 (C. Ethe, G. Madec) merge TRA-TRC
	13	!! - ! 2010-11 (G. Madec) add mbk. arrays associated to the deepest ocean level
[4990]	14	!! - ! 2013-04 (F. Roquet, G. Madec) use of eosbn2 instead of local hard coded alpha and beta
[9019]	15	!! 4.0 ! 2017-04 (G. Madec) ln_trabbl namelist variable instead of a CPP key
[503]	16	!!----------------------------------------------------------------------
[9019]	17
[3]	18	!!----------------------------------------------------------------------
[2715]	19	!! tra_bbl_alloc : allocate trabbl arrays
[2528]	20	!! tra_bbl : update the tracer trends due to the bottom boundary layer (advective and/or diffusive)
	21	!! tra_bbl_dif : generic routine to compute bbl diffusive trend
	22	!! tra_bbl_adv : generic routine to compute bbl advective trend
	23	!! bbl : computation of bbl diffu. flux coef. & transport in bottom boundary layer
	24	!! tra_bbl_init : initialization, namelist read, parameters control
[503]	25	!!----------------------------------------------------------------------
[2528]	26	USE oce ! ocean dynamics and active tracers
	27	USE dom_oce ! ocean space and time domain
	28	USE phycst ! physical constant
	29	USE eosbn2 ! equation of state
[5836]	30	USE trd_oce ! trends: ocean variables
[2528]	31	USE trdtra ! trends: active tracers
[4990]	32	!
	33	USE iom ! IOM library
[2528]	34	USE in_out_manager ! I/O manager
	35	USE lbclnk ! ocean lateral boundary conditions
	36	USE prtctl ! Print control
[3294]	37	USE timing ! Timing
[4990]	38	USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
[3]	39
	40	IMPLICIT NONE
	41	PRIVATE
	42
[2528]	43	PUBLIC tra_bbl ! routine called by step.F90
[9124]	44	PUBLIC tra_bbl_init ! routine called by nemogcm.F90
[2528]	45	PUBLIC tra_bbl_dif ! routine called by trcbbl.F90
[9124]	46	PUBLIC tra_bbl_adv ! - - -
[2528]	47	PUBLIC bbl ! routine called by trcbbl.F90 and dtadyn.F90
[3]	48
[4147]	49	! !!* Namelist nambbl *
[9019]	50	LOGICAL , PUBLIC :: ln_trabbl !: bottom boundary layer flag
[4147]	51	INTEGER , PUBLIC :: nn_bbl_ldf !: =1 : diffusive bbl or not (=0)
	52	INTEGER , PUBLIC :: nn_bbl_adv !: =1/2 : advective bbl or not (=0)
[2528]	53	! ! =1 : advective bbl using the bottom ocean velocity
	54	! ! =2 : - - using utr_bbl proportional to grad(rho)
[4147]	55	REAL(wp), PUBLIC :: rn_ahtbbl !: along slope bbl diffusive coefficient [m2/s]
	56	REAL(wp), PUBLIC :: rn_gambbl !: lateral coeff. for bottom boundary layer scheme [s]
[409]	57
[4990]	58	LOGICAL , PUBLIC :: l_bbl !: flag to compute bbl diffu. flux coef and transport
[3764]	59
[2715]	60	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: utr_bbl , vtr_bbl ! u- (v-) transport in the bottom boundary layer
	61	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: ahu_bbl , ahv_bbl ! masked diffusive bbl coeff. at u & v-pts
[3]	62
[3764]	63	INTEGER , ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: mbku_d , mbkv_d ! vertical index of the "lower" bottom ocean U/V-level (PUBLIC for TAM)
	64	INTEGER , ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: mgrhu , mgrhv ! = +/-1, sign of grad(H) in u-(v-)direction (PUBLIC for TAM)
	65	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ahu_bbl_0, ahv_bbl_0 ! diffusive bbl flux coefficients at u and v-points
	66	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: e3u_bbl_0, e3v_bbl_0 ! thichness of the bbl (e3) at u and v-points (PUBLIC for TAM)
[3]	67
	68	!! * Substitutions
	69	# include "vectopt_loop_substitute.h90"
	70	!!----------------------------------------------------------------------
[9598]	71	!! NEMO/OCE 4.0 , NEMO Consortium (2018)
[2528]	72	!! $Id$
[10068]	73	!! Software governed by the CeCILL license (see ./LICENSE)
[3]	74	!!----------------------------------------------------------------------
	75	CONTAINS
	76
[2715]	77	INTEGER FUNCTION tra_bbl_alloc()
	78	!!----------------------------------------------------------------------
	79	!! * FUNCTION tra_bbl_alloc *
	80	!!----------------------------------------------------------------------
[9019]	81	ALLOCATE( utr_bbl (jpi,jpj) , ahu_bbl (jpi,jpj) , mbku_d(jpi,jpj) , mgrhu(jpi,jpj) , &
	82	& vtr_bbl (jpi,jpj) , ahv_bbl (jpi,jpj) , mbkv_d(jpi,jpj) , mgrhv(jpi,jpj) , &
	83	& ahu_bbl_0(jpi,jpj) , ahv_bbl_0(jpi,jpj) , &
	84	& e3u_bbl_0(jpi,jpj) , e3v_bbl_0(jpi,jpj) , STAT=tra_bbl_alloc )
[2715]	85	!
[10402]	86	CALL mpp_sum ( 'trabbl', tra_bbl_alloc )
[2715]	87	IF( tra_bbl_alloc > 0 ) CALL ctl_warn('tra_bbl_alloc: allocation of arrays failed.')
	88	END FUNCTION tra_bbl_alloc
	89
	90
[2528]	91	SUBROUTINE tra_bbl( kt )
[3]	92	!!----------------------------------------------------------------------
[2528]	93	!! * ROUTINE bbl *
[3764]	94	!!
	95	!! ** Purpose : Compute the before tracer (t & s) trend associated
[2528]	96	!! with the bottom boundary layer and add it to the general
	97	!! trend of tracer equations.
[3]	98	!!
[2528]	99	!! ** Method : Depending on namtra_bbl namelist parameters the bbl
	100	!! diffusive and/or advective contribution to the tracer trend
	101	!! is added to the general tracer trend
[3764]	102	!!----------------------------------------------------------------------
	103	INTEGER, INTENT( in ) :: kt ! ocean time-step
[4990]	104	!
[9019]	105	REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdt, ztrds
[2528]	106	!!----------------------------------------------------------------------
[3294]	107	!
[9019]	108	IF( ln_timing ) CALL timing_start( 'tra_bbl')
[3294]	109	!
[9019]	110	IF( l_trdtra ) THEN !* Save the T-S input trends
	111	ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) )
[7753]	112	ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
	113	ztrds(:,:,:) = tsa(:,:,:,jp_sal)
[2528]	114	ENDIF
	115
[4990]	116	IF( l_bbl ) CALL bbl( kt, nit000, 'TRA' ) !* bbl coef. and transport (only if not already done in trcbbl)
[3764]	117
[4990]	118	IF( nn_bbl_ldf == 1 ) THEN !* Diffusive bbl
[3294]	119	!
[2528]	120	CALL tra_bbl_dif( tsb, tsa, jpts )
	121	IF( ln_ctl ) &
	122	CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' bbl_ldf - Ta: ', mask1=tmask, &
[4990]	123	& tab3d_2=tsa(:,:,:,jp_sal), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' )
[3764]	124	! lateral boundary conditions ; just need for outputs
[10170]	125	CALL lbc_lnk_multi( 'trabbl', ahu_bbl, 'U', 1. , ahv_bbl, 'V', 1. )
[3764]	126	CALL iom_put( "ahu_bbl", ahu_bbl ) ! bbl diffusive flux i-coef
[2528]	127	CALL iom_put( "ahv_bbl", ahv_bbl ) ! bbl diffusive flux j-coef
[3294]	128	!
[9168]	129	ENDIF
[6140]	130	!
[4990]	131	IF( nn_bbl_adv /= 0 ) THEN !* Advective bbl
[3294]	132	!
[2528]	133	CALL tra_bbl_adv( tsb, tsa, jpts )
	134	IF(ln_ctl) &
	135	CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' bbl_adv - Ta: ', mask1=tmask, &
[4990]	136	& tab3d_2=tsa(:,:,:,jp_sal), clinfo2= ' Sa: ', mask2=tmask, clinfo3='tra' )
[3764]	137	! lateral boundary conditions ; just need for outputs
[10170]	138	CALL lbc_lnk_multi( 'trabbl', utr_bbl, 'U', 1. , vtr_bbl, 'V', 1. )
[3764]	139	CALL iom_put( "uoce_bbl", utr_bbl ) ! bbl i-transport
[2528]	140	CALL iom_put( "voce_bbl", vtr_bbl ) ! bbl j-transport
[3294]	141	!
[9168]	142	ENDIF
[2528]	143
[6140]	144	IF( l_trdtra ) THEN ! send the trends for further diagnostics
[7753]	145	ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
	146	ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)
[4990]	147	CALL trd_tra( kt, 'TRA', jp_tem, jptra_bbl, ztrdt )
	148	CALL trd_tra( kt, 'TRA', jp_sal, jptra_bbl, ztrds )
[9019]	149	DEALLOCATE( ztrdt, ztrds )
[2528]	150	ENDIF
	151	!
[9019]	152	IF( ln_timing ) CALL timing_stop( 'tra_bbl')
[3294]	153	!
[2528]	154	END SUBROUTINE tra_bbl
	155
	156
	157	SUBROUTINE tra_bbl_dif( ptb, pta, kjpt )
	158	!!----------------------------------------------------------------------
	159	!! * ROUTINE tra_bbl_dif *
[3764]	160	!!
[2528]	161	!! ** Purpose : Computes the bottom boundary horizontal and vertical
[3764]	162	!! advection terms.
[2528]	163	!!
[4990]	164	!! ** Method : * diffusive bbl only (nn_bbl_ldf=1) :
[2528]	165	!! When the product grad( rho) * grad(h) < 0 (where grad is an
	166	!! along bottom slope gradient) an additional lateral 2nd order
	167	!! diffusion along the bottom slope is added to the general
	168	!! tracer trend, otherwise the additional trend is set to 0.
	169	!! A typical value of ahbt is 2000 m2/s (equivalent to
[3]	170	!! a downslope velocity of 20 cm/s if the condition for slope
	171	!! convection is satified)
	172	!!
[2528]	173	!! ** Action : pta increased by the bbl diffusive trend
[3]	174	!!
[503]	175	!! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591.
[2528]	176	!! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430.
[3764]	177	!!----------------------------------------------------------------------
[2715]	178	INTEGER , INTENT(in ) :: kjpt ! number of tracers
	179	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields
[3764]	180	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend
[2715]	181	!
[2528]	182	INTEGER :: ji, jj, jn ! dummy loop indices
	183	INTEGER :: ik ! local integers
	184	REAL(wp) :: zbtr ! local scalars
[9019]	185	REAL(wp), DIMENSION(jpi,jpj) :: zptb ! workspace
[3]	186	!!----------------------------------------------------------------------
[2528]	187	!
	188	DO jn = 1, kjpt ! tracer loop
	189	! ! ===========
	190	DO jj = 1, jpj
[3764]	191	DO ji = 1, jpi
[5836]	192	ik = mbkt(ji,jj) ! bottom T-level index
	193	zptb(ji,jj) = ptb(ji,jj,ik,jn) ! bottom before T and S
[2528]	194	END DO
	195	END DO
[4990]	196	!
	197	DO jj = 2, jpjm1 ! Compute the trend
[2528]	198	DO ji = 2, jpim1
[5836]	199	ik = mbkt(ji,jj) ! bottom T-level index
	200	pta(ji,jj,ik,jn) = pta(ji,jj,ik,jn) &
	201	& + ( ahu_bbl(ji ,jj ) * ( zptb(ji+1,jj ) - zptb(ji ,jj ) ) &
	202	& - ahu_bbl(ji-1,jj ) * ( zptb(ji ,jj ) - zptb(ji-1,jj ) ) &
	203	& + ahv_bbl(ji ,jj ) * ( zptb(ji ,jj+1) - zptb(ji ,jj ) ) &
	204	& - ahv_bbl(ji ,jj-1) * ( zptb(ji ,jj ) - zptb(ji ,jj-1) ) ) &
[6140]	205	& * r1_e1e2t(ji,jj) / e3t_n(ji,jj,ik)
[2528]	206	END DO
[3]	207	END DO
[2528]	208	! ! ===========
	209	END DO ! end tracer
	210	! ! ===========
	211	END SUBROUTINE tra_bbl_dif
[3]	212
[3764]	213
[2528]	214	SUBROUTINE tra_bbl_adv( ptb, pta, kjpt )
	215	!!----------------------------------------------------------------------
	216	!! * ROUTINE trc_bbl *
	217	!!
[3764]	218	!! ** Purpose : Compute the before passive tracer trend associated
[2528]	219	!! with the bottom boundary layer and add it to the general trend
	220	!! of tracer equations.
	221	!! ** Method : advective bbl (nn_bbl_adv = 1 or 2) :
	222	!! nn_bbl_adv = 1 use of the ocean near bottom velocity as bbl velocity
[3764]	223	!! nn_bbl_adv = 2 follow Campin and Goosse (1999) implentation i.e.
	224	!! transport proportional to the along-slope density gradient
[2528]	225	!!
	226	!! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591.
	227	!! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430.
[3764]	228	!!----------------------------------------------------------------------
[2715]	229	INTEGER , INTENT(in ) :: kjpt ! number of tracers
	230	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields
[3764]	231	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend
[2715]	232	!
[2528]	233	INTEGER :: ji, jj, jk, jn ! dummy loop indices
	234	INTEGER :: iis , iid , ijs , ijd ! local integers
	235	INTEGER :: ikus, ikud, ikvs, ikvd ! - -
	236	REAL(wp) :: zbtr, ztra ! local scalars
	237	REAL(wp) :: zu_bbl, zv_bbl ! - -
	238	!!----------------------------------------------------------------------
	239	!
	240	! ! ===========
	241	DO jn = 1, kjpt ! tracer loop
[3764]	242	! ! ===========
[457]	243	DO jj = 1, jpjm1
[2528]	244	DO ji = 1, jpim1 ! CAUTION start from i=1 to update i=2 when cyclic east-west
	245	IF( utr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero i-direction bbl advection
	246	! down-slope i/k-indices (deep) & up-slope i/k indices (shelf)
	247	iid = ji + MAX( 0, mgrhu(ji,jj) ) ; iis = ji + 1 - MAX( 0, mgrhu(ji,jj) )
	248	ikud = mbku_d(ji,jj) ; ikus = mbku(ji,jj)
	249	zu_bbl = ABS( utr_bbl(ji,jj) )
	250	!
	251	! ! up -slope T-point (shelf bottom point)
[6140]	252	zbtr = r1_e1e2t(iis,jj) / e3t_n(iis,jj,ikus)
[2528]	253	ztra = zu_bbl * ( ptb(iid,jj,ikus,jn) - ptb(iis,jj,ikus,jn) ) * zbtr
	254	pta(iis,jj,ikus,jn) = pta(iis,jj,ikus,jn) + ztra
[3764]	255	!
[2528]	256	DO jk = ikus, ikud-1 ! down-slope upper to down T-point (deep column)
[6140]	257	zbtr = r1_e1e2t(iid,jj) / e3t_n(iid,jj,jk)
[2528]	258	ztra = zu_bbl * ( ptb(iid,jj,jk+1,jn) - ptb(iid,jj,jk,jn) ) * zbtr
	259	pta(iid,jj,jk,jn) = pta(iid,jj,jk,jn) + ztra
	260	END DO
[3764]	261	!
[6140]	262	zbtr = r1_e1e2t(iid,jj) / e3t_n(iid,jj,ikud)
[2528]	263	ztra = zu_bbl * ( ptb(iis,jj,ikus,jn) - ptb(iid,jj,ikud,jn) ) * zbtr
	264	pta(iid,jj,ikud,jn) = pta(iid,jj,ikud,jn) + ztra
	265	ENDIF
	266	!
	267	IF( vtr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero j-direction bbl advection
	268	! down-slope j/k-indices (deep) & up-slope j/k indices (shelf)
	269	ijd = jj + MAX( 0, mgrhv(ji,jj) ) ; ijs = jj + 1 - MAX( 0, mgrhv(ji,jj) )
	270	ikvd = mbkv_d(ji,jj) ; ikvs = mbkv(ji,jj)
	271	zv_bbl = ABS( vtr_bbl(ji,jj) )
[3764]	272	!
[2528]	273	! up -slope T-point (shelf bottom point)
[6140]	274	zbtr = r1_e1e2t(ji,ijs) / e3t_n(ji,ijs,ikvs)
[2528]	275	ztra = zv_bbl * ( ptb(ji,ijd,ikvs,jn) - ptb(ji,ijs,ikvs,jn) ) * zbtr
	276	pta(ji,ijs,ikvs,jn) = pta(ji,ijs,ikvs,jn) + ztra
[3764]	277	!
[2528]	278	DO jk = ikvs, ikvd-1 ! down-slope upper to down T-point (deep column)
[6140]	279	zbtr = r1_e1e2t(ji,ijd) / e3t_n(ji,ijd,jk)
[2528]	280	ztra = zv_bbl * ( ptb(ji,ijd,jk+1,jn) - ptb(ji,ijd,jk,jn) ) * zbtr
	281	pta(ji,ijd,jk,jn) = pta(ji,ijd,jk,jn) + ztra
	282	END DO
	283	! ! down-slope T-point (deep bottom point)
[6140]	284	zbtr = r1_e1e2t(ji,ijd) / e3t_n(ji,ijd,ikvd)
[2528]	285	ztra = zv_bbl * ( ptb(ji,ijs,ikvs,jn) - ptb(ji,ijd,ikvd,jn) ) * zbtr
	286	pta(ji,ijd,ikvd,jn) = pta(ji,ijd,ikvd,jn) + ztra
	287	ENDIF
[457]	288	END DO
[2528]	289	!
[457]	290	END DO
[9019]	291	! ! ===========
	292	END DO ! end tracer
	293	! ! ===========
[2528]	294	END SUBROUTINE tra_bbl_adv
[3]	295
	296
[3294]	297	SUBROUTINE bbl( kt, kit000, cdtype )
[2528]	298	!!----------------------------------------------------------------------
	299	!! * ROUTINE bbl *
[3764]	300	!!
[2528]	301	!! ** Purpose : Computes the bottom boundary horizontal and vertical
[3764]	302	!! advection terms.
[2528]	303	!!
[4990]	304	!! ** Method : * diffusive bbl (nn_bbl_ldf=1) :
[2528]	305	!! When the product grad( rho) * grad(h) < 0 (where grad is an
	306	!! along bottom slope gradient) an additional lateral 2nd order
	307	!! diffusion along the bottom slope is added to the general
	308	!! tracer trend, otherwise the additional trend is set to 0.
	309	!! A typical value of ahbt is 2000 m2/s (equivalent to
	310	!! a downslope velocity of 20 cm/s if the condition for slope
	311	!! convection is satified)
[4990]	312	!! * advective bbl (nn_bbl_adv=1 or 2) :
[2528]	313	!! nn_bbl_adv = 1 use of the ocean velocity as bbl velocity
	314	!! nn_bbl_adv = 2 follow Campin and Goosse (1999) implentation
	315	!! i.e. transport proportional to the along-slope density gradient
	316	!!
	317	!! NB: the along slope density gradient is evaluated using the
	318	!! local density (i.e. referenced at a common local depth).
	319	!!
	320	!! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591.
	321	!! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430.
[3764]	322	!!----------------------------------------------------------------------
[2528]	323	INTEGER , INTENT(in ) :: kt ! ocean time-step index
[4990]	324	INTEGER , INTENT(in ) :: kit000 ! first time step index
[2528]	325	CHARACTER(len=3), INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator)
[6140]	326	!
[2528]	327	INTEGER :: ji, jj ! dummy loop indices
	328	INTEGER :: ik ! local integers
[4990]	329	INTEGER :: iis, iid, ikus, ikud ! - -
	330	INTEGER :: ijs, ijd, ikvs, ikvd ! - -
	331	REAL(wp) :: za, zb, zgdrho ! local scalars
	332	REAL(wp) :: zsign, zsigna, zgbbl ! - -
	333	REAL(wp), DIMENSION(jpi,jpj,jpts) :: zts, zab ! 3D workspace
	334	REAL(wp), DIMENSION(jpi,jpj) :: zub, zvb, zdep ! 2D workspace
[2528]	335	!!----------------------------------------------------------------------
[3294]	336	!
	337	IF( kt == kit000 ) THEN
[2528]	338	IF(lwp) WRITE(numout,*)
	339	IF(lwp) WRITE(numout,*) 'trabbl:bbl : Compute bbl velocities and diffusive coefficients in ', cdtype
	340	IF(lwp) WRITE(numout,*) '~~~~~~~~~~'
	341	ENDIF
[4990]	342	! !* bottom variables (T, S, alpha, beta, depth, velocity)
[2528]	343	DO jj = 1, jpj
	344	DO ji = 1, jpi
[4990]	345	ik = mbkt(ji,jj) ! bottom T-level index
	346	zts (ji,jj,jp_tem) = tsb(ji,jj,ik,jp_tem) ! bottom before T and S
	347	zts (ji,jj,jp_sal) = tsb(ji,jj,ik,jp_sal)
[2528]	348	!
[6140]	349	zdep(ji,jj) = gdept_n(ji,jj,ik) ! bottom T-level reference depth
[4990]	350	zub (ji,jj) = un(ji,jj,mbku(ji,jj)) ! bottom velocity
	351	zvb (ji,jj) = vn(ji,jj,mbkv(ji,jj))
[2528]	352	END DO
	353	END DO
[4990]	354	!
	355	CALL eos_rab( zts, zdep, zab )
	356	!
[2528]	357	! !-------------------!
	358	IF( nn_bbl_ldf == 1 ) THEN ! diffusive bbl !
	359	! !-------------------!
	360	DO jj = 1, jpjm1 ! (criteria for non zero flux: grad(rho).grad(h) < 0 )
[4990]	361	DO ji = 1, fs_jpim1 ! vector opt.
	362	! ! i-direction
	363	za = zab(ji+1,jj,jp_tem) + zab(ji,jj,jp_tem) ! 2*(alpha,beta) at u-point
	364	zb = zab(ji+1,jj,jp_sal) + zab(ji,jj,jp_sal)
	365	! ! 2*masked bottom density gradient
	366	zgdrho = ( za * ( zts(ji+1,jj,jp_tem) - zts(ji,jj,jp_tem) ) &
	367	& - zb * ( zts(ji+1,jj,jp_sal) - zts(ji,jj,jp_sal) ) ) * umask(ji,jj,1)
[3764]	368	!
[4990]	369	zsign = SIGN( 0.5, -zgdrho * REAL( mgrhu(ji,jj) ) ) ! sign of ( i-gradient * i-slope )
	370	ahu_bbl(ji,jj) = ( 0.5 - zsign ) * ahu_bbl_0(ji,jj) ! masked diffusive flux coeff.
[1601]	371	!
[4990]	372	! ! j-direction
	373	za = zab(ji,jj+1,jp_tem) + zab(ji,jj,jp_tem) ! 2*(alpha,beta) at v-point
	374	zb = zab(ji,jj+1,jp_sal) + zab(ji,jj,jp_sal)
	375	! ! 2*masked bottom density gradient
	376	zgdrho = ( za * ( zts(ji,jj+1,jp_tem) - zts(ji,jj,jp_tem) ) &
	377	& - zb * ( zts(ji,jj+1,jp_sal) - zts(ji,jj,jp_sal) ) ) * vmask(ji,jj,1)
[3764]	378	!
[4990]	379	zsign = SIGN( 0.5, -zgdrho * REAL( mgrhv(ji,jj) ) ) ! sign of ( j-gradient * j-slope )
[2528]	380	ahv_bbl(ji,jj) = ( 0.5 - zsign ) * ahv_bbl_0(ji,jj)
[1601]	381	END DO
[3]	382	END DO
[1601]	383	!
[2528]	384	ENDIF
[6140]	385	!
[2528]	386	! !-------------------!
	387	IF( nn_bbl_adv /= 0 ) THEN ! advective bbl !
	388	! !-------------------!
	389	SELECT CASE ( nn_bbl_adv ) !* bbl transport type
[503]	390	!
[2528]	391	CASE( 1 ) != use of upper velocity
	392	DO jj = 1, jpjm1 ! criteria: grad(rho).grad(h)<0 and grad(rho).grad(h)<0
	393	DO ji = 1, fs_jpim1 ! vector opt.
[4990]	394	! ! i-direction
	395	za = zab(ji+1,jj,jp_tem) + zab(ji,jj,jp_tem) ! 2*(alpha,beta) at u-point
	396	zb = zab(ji+1,jj,jp_sal) + zab(ji,jj,jp_sal)
	397	! ! 2*masked bottom density gradient
	398	zgdrho = ( za * ( zts(ji+1,jj,jp_tem) - zts(ji,jj,jp_tem) ) &
	399	- zb * ( zts(ji+1,jj,jp_sal) - zts(ji,jj,jp_sal) ) ) * umask(ji,jj,1)
[3764]	400	!
[4990]	401	zsign = SIGN( 0.5, - zgdrho * REAL( mgrhu(ji,jj) ) ) ! sign of i-gradient * i-slope
	402	zsigna= SIGN( 0.5, zub(ji,jj) * REAL( mgrhu(ji,jj) ) ) ! sign of u * i-slope
[2528]	403	!
[4990]	404	! ! bbl velocity
[2528]	405	utr_bbl(ji,jj) = ( 0.5 + zsigna ) * ( 0.5 - zsign ) * e2u(ji,jj) * e3u_bbl_0(ji,jj) * zub(ji,jj)
	406	!
[4990]	407	! ! j-direction
	408	za = zab(ji,jj+1,jp_tem) + zab(ji,jj,jp_tem) ! 2*(alpha,beta) at v-point
	409	zb = zab(ji,jj+1,jp_sal) + zab(ji,jj,jp_sal)
	410	! ! 2*masked bottom density gradient
	411	zgdrho = ( za * ( zts(ji,jj+1,jp_tem) - zts(ji,jj,jp_tem) ) &
	412	& - zb * ( zts(ji,jj+1,jp_sal) - zts(ji,jj,jp_sal) ) ) * vmask(ji,jj,1)
	413	zsign = SIGN( 0.5, - zgdrho * REAL( mgrhv(ji,jj) ) ) ! sign of j-gradient * j-slope
	414	zsigna= SIGN( 0.5, zvb(ji,jj) * REAL( mgrhv(ji,jj) ) ) ! sign of u * i-slope
[2528]	415	!
[4990]	416	! ! bbl transport
[2528]	417	vtr_bbl(ji,jj) = ( 0.5 + zsigna ) * ( 0.5 - zsign ) * e1v(ji,jj) * e3v_bbl_0(ji,jj) * zvb(ji,jj)
	418	END DO
	419	END DO
[503]	420	!
[2528]	421	CASE( 2 ) != bbl velocity = F( delta rho )
	422	zgbbl = grav * rn_gambbl
	423	DO jj = 1, jpjm1 ! criteria: rho_up > rho_down
	424	DO ji = 1, fs_jpim1 ! vector opt.
[4990]	425	! ! i-direction
[2528]	426	! down-slope T-point i/k-index (deep) & up-slope T-point i/k-index (shelf)
[4990]	427	iid = ji + MAX( 0, mgrhu(ji,jj) )
	428	iis = ji + 1 - MAX( 0, mgrhu(ji,jj) )
[2528]	429	!
[4990]	430	ikud = mbku_d(ji,jj)
	431	ikus = mbku(ji,jj)
[2528]	432	!
[4990]	433	za = zab(ji+1,jj,jp_tem) + zab(ji,jj,jp_tem) ! 2*(alpha,beta) at u-point
	434	zb = zab(ji+1,jj,jp_sal) + zab(ji,jj,jp_sal)
	435	! ! masked bottom density gradient
	436	zgdrho = 0.5 * ( za * ( zts(iid,jj,jp_tem) - zts(iis,jj,jp_tem) ) &
	437	& - zb * ( zts(iid,jj,jp_sal) - zts(iis,jj,jp_sal) ) ) * umask(ji,jj,1)
	438	zgdrho = MAX( 0.e0, zgdrho ) ! only if shelf is denser than deep
	439	!
	440	! ! bbl transport (down-slope direction)
[2528]	441	utr_bbl(ji,jj) = e2u(ji,jj) * e3u_bbl_0(ji,jj) * zgbbl * zgdrho * REAL( mgrhu(ji,jj) )
	442	!
[4990]	443	! ! j-direction
[2528]	444	! down-slope T-point j/k-index (deep) & of the up -slope T-point j/k-index (shelf)
[4990]	445	ijd = jj + MAX( 0, mgrhv(ji,jj) )
	446	ijs = jj + 1 - MAX( 0, mgrhv(ji,jj) )
[2528]	447	!
[4990]	448	ikvd = mbkv_d(ji,jj)
	449	ikvs = mbkv(ji,jj)
[2528]	450	!
[4990]	451	za = zab(ji,jj+1,jp_tem) + zab(ji,jj,jp_tem) ! 2*(alpha,beta) at v-point
	452	zb = zab(ji,jj+1,jp_sal) + zab(ji,jj,jp_sal)
	453	! ! masked bottom density gradient
	454	zgdrho = 0.5 * ( za * ( zts(ji,ijd,jp_tem) - zts(ji,ijs,jp_tem) ) &
	455	& - zb * ( zts(ji,ijd,jp_sal) - zts(ji,ijs,jp_sal) ) ) * vmask(ji,jj,1)
	456	zgdrho = MAX( 0.e0, zgdrho ) ! only if shelf is denser than deep
	457	!
	458	! ! bbl transport (down-slope direction)
[2528]	459	vtr_bbl(ji,jj) = e1v(ji,jj) * e3v_bbl_0(ji,jj) * zgbbl * zgdrho * REAL( mgrhv(ji,jj) )
	460	END DO
	461	END DO
[3]	462	END SELECT
[2528]	463	!
[3]	464	ENDIF
[503]	465	!
[2528]	466	END SUBROUTINE bbl
[3]	467
	468
	469	SUBROUTINE tra_bbl_init
	470	!!----------------------------------------------------------------------
	471	!! * ROUTINE tra_bbl_init *
	472	!!
	473	!! ** Purpose : Initialization for the bottom boundary layer scheme.
	474	!!
	475	!! ** Method : Read the nambbl namelist and check the parameters
[3294]	476	!! called by nemo_init at the first timestep (kit000)
[3]	477	!!----------------------------------------------------------------------
[9019]	478	INTEGER :: ji, jj ! dummy loop indices
	479	INTEGER :: ii0, ii1, ij0, ij1, ios ! local integer
[9094]	480	REAL(wp), DIMENSION(jpi,jpj) :: zmbku, zmbkv ! workspace
[9019]	481	!!
	482	NAMELIST/nambbl/ ln_trabbl, nn_bbl_ldf, nn_bbl_adv, rn_ahtbbl, rn_gambbl
[3]	483	!!----------------------------------------------------------------------
[3294]	484	!
[4147]	485	REWIND( numnam_ref ) ! Namelist nambbl in reference namelist : Bottom boundary layer scheme
	486	READ ( numnam_ref, nambbl, IOSTAT = ios, ERR = 901)
[6140]	487	901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambbl in reference namelist', lwp )
	488	!
[4147]	489	REWIND( numnam_cfg ) ! Namelist nambbl in configuration namelist : Bottom boundary layer scheme
	490	READ ( numnam_cfg, nambbl, IOSTAT = ios, ERR = 902 )
[9168]	491	902 IF( ios > 0 ) CALL ctl_nam ( ios , 'nambbl in configuration namelist', lwp )
[4624]	492	IF(lwm) WRITE ( numond, nambbl )
[2528]	493	!
	494	l_bbl = .TRUE. !* flag to compute bbl coef and transport
	495	!
	496	IF(lwp) THEN !* Parameter control and print
[3]	497	WRITE(numout,*)
[2528]	498	WRITE(numout,*) 'tra_bbl_init : bottom boundary layer initialisation'
[3]	499	WRITE(numout,*) '~~~~~~~~~~~~'
[9019]	500	WRITE(numout,*) ' Namelist nambbl : set bbl parameters'
	501	WRITE(numout,*) ' bottom boundary layer flag ln_trabbl = ', ln_trabbl
[3]	502	ENDIF
[9019]	503	IF( .NOT.ln_trabbl ) RETURN
	504	!
	505	IF(lwp) THEN
	506	WRITE(numout,*) ' diffusive bbl (=1) or not (=0) nn_bbl_ldf = ', nn_bbl_ldf
	507	WRITE(numout,*) ' advective bbl (=1/2) or not (=0) nn_bbl_adv = ', nn_bbl_adv
	508	WRITE(numout,*) ' diffusive bbl coefficient rn_ahtbbl = ', rn_ahtbbl, ' m2/s'
	509	WRITE(numout,*) ' advective bbl coefficient rn_gambbl = ', rn_gambbl, ' s'
	510	ENDIF
	511	!
[2715]	512	! ! allocate trabbl arrays
	513	IF( tra_bbl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'tra_bbl_init : unable to allocate arrays' )
[9019]	514	!
[2528]	515	IF( nn_bbl_adv == 1 ) WRITE(numout,) ' Advective BBL using upper velocity'
	516	IF( nn_bbl_adv == 2 ) WRITE(numout,) ' Advective BBL using velocity = F( delta rho)'
[9019]	517	!
[2528]	518	! !* vertical index of "deep" bottom u- and v-points
	519	DO jj = 1, jpjm1 ! (the "shelf" bottom k-indices are mbku and mbkv)
	520	DO ji = 1, jpim1
	521	mbku_d(ji,jj) = MAX( mbkt(ji+1,jj ) , mbkt(ji,jj) ) ! >= 1 as mbkt=1 over land
	522	mbkv_d(ji,jj) = MAX( mbkt(ji ,jj+1) , mbkt(ji,jj) )
[3]	523	END DO
	524	END DO
[4990]	525	! converte into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk
[9094]	526	zmbku(:,:) = REAL( mbku_d(:,:), wp ) ; zmbkv(:,:) = REAL( mbkv_d(:,:), wp )
[10170]	527	CALL lbc_lnk_multi( 'trabbl', zmbku,'U',1., zmbkv,'V',1.)
[9919]	528	mbku_d(:,:) = MAX( INT( zmbku(:,:) ), 1 ) ; mbkv_d(:,:) = MAX( NINT( zmbkv(:,:) ), 1 )
[9019]	529	!
[9168]	530	! !* sign of grad(H) at u- and v-points; zero if grad(H) = 0
[8509]	531	mgrhu(:,:) = 0 ; mgrhv(:,:) = 0
[3764]	532	DO jj = 1, jpjm1
[3]	533	DO ji = 1, jpim1
[8509]	534	IF( gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) /= 0._wp ) THEN
	535	mgrhu(ji,jj) = INT( SIGN( 1.e0, gdept_0(ji+1,jj,mbkt(ji+1,jj)) - gdept_0(ji,jj,mbkt(ji,jj)) ) )
	536	ENDIF
	537	!
	538	IF( gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) /= 0._wp ) THEN
	539	mgrhv(ji,jj) = INT( SIGN( 1.e0, gdept_0(ji,jj+1,mbkt(ji,jj+1)) - gdept_0(ji,jj,mbkt(ji,jj)) ) )
	540	ENDIF
[3]	541	END DO
	542	END DO
[7646]	543	!
[3764]	544	DO jj = 1, jpjm1 !* bbl thickness at u- (v-) point
[4990]	545	DO ji = 1, jpim1 ! minimum of top & bottom e3u_0 (e3v_0)
[4292]	546	e3u_bbl_0(ji,jj) = MIN( e3u_0(ji,jj,mbkt(ji+1,jj )), e3u_0(ji,jj,mbkt(ji,jj)) )
	547	e3v_bbl_0(ji,jj) = MIN( e3v_0(ji,jj,mbkt(ji ,jj+1)), e3v_0(ji,jj,mbkt(ji,jj)) )
[3764]	548	END DO
[2528]	549	END DO
[10170]	550	CALL lbc_lnk_multi( 'trabbl', e3u_bbl_0, 'U', 1. , e3v_bbl_0, 'V', 1. ) ! lateral boundary conditions
[7646]	551	!
[3764]	552	! !* masked diffusive flux coefficients
[7753]	553	ahu_bbl_0(:,:) = rn_ahtbbl * e2_e1u(:,:) * e3u_bbl_0(:,:) * umask(:,:,1)
	554	ahv_bbl_0(:,:) = rn_ahtbbl * e1_e2v(:,:) * e3v_bbl_0(:,:) * vmask(:,:,1)
[503]	555	!
[3]	556	END SUBROUTINE tra_bbl_init
	557
	558	!!======================================================================
	559	END MODULE trabbl

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

Download in other formats: