Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

trabbl.F90 @ 9157

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

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

New URL for NEMO forge! http://forge.nemo-ocean.eu

Context Navigation

source: branches/NERC/dev_r5518_GO6_Carb_Fail_from_GO6_8356/NEMOGCM/NEMO/OPA_SRC/TRA/trabbl.F90 @ 9157

Download in other formats: