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_tam.F90 in branches/2012/dev_r3604_LEGI8_TAM/NEMOGCM/NEMO/OPATAM_SRC/TRA – NEMO

Context Navigation

source: branches/2012/dev_r3604_LEGI8_TAM/NEMOGCM/NEMO/OPATAM_SRC/TRA/trabbl_tam.F90 @ 3627

Last change on this file since 3627 was 3627, checked in by rblod, 11 years ago
Correct long lines in TAM 4.3 see ticket #1010
Property svn:executable set to ``*
File size: 85.9 KB

Rev	Line
[3611]	1	MODULE trabbl_tam
	2	!!==============================================================================
	3	!! * MODULE trabbl *
	4	!! Ocean physics : advective and/or diffusive bottom boundary layer scheme
	5	!!==============================================================================
	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
	10	!! 3.3 ! 2009-11 (G. Madec) merge trabbl and trabbl_adv + style + optimization
	11	!! - ! 2010-04 (G. Madec) Campin & Goosse advective bbl
	12	!! - ! 2010-06 (C. Ethe, G. Madec) merge TRA-TRC
	13	!! - ! 2010-11 (G. Madec) add mbk. arrays associated to the deepest ocean level
	14	!! History of the T&A module
	15	!! NEMO 3.2 ! 2011-02 (A. Vidard) Original version
	16	!! 3.4 ! 2012-09 (A. Vidard) Update to 3.4
	17	!!
	18	!!----------------------------------------------------------------------
	19	#if defined key_trabbl \|\| defined key_esopa
	20	!!----------------------------------------------------------------------
	21	!! 'key_trabbl' or bottom boundary layer
	22	!!----------------------------------------------------------------------
	23	!! tra_bbl_alloc : allocate trabbl arrays
	24	!! tra_bbl : update the tracer trends due to the bottom boundary layer (advective and/or diffusive)
	25	!! tra_bbl_dif : generic routine to compute bbl diffusive trend
	26	!! tra_bbl_adv : generic routine to compute bbl advective trend
	27	!! bbl : computation of bbl diffu. flux coef. & transport in bottom boundary layer
	28	!! tra_bbl_init : initialization, namelist read, parameters control
	29	!!----------------------------------------------------------------------
	30	USE oce ! ocean dynamics and active tracers
	31	USE oce_tam
	32	USE dom_oce ! ocean space and time domain
	33	USE phycst ! physical constant
	34	USE eosbn2 ! equation of state
	35	USE iom ! IOM server
	36	USE in_out_manager ! I/O manager
	37	USE lbclnk ! ocean lateral boundary conditions
	38	USE prtctl ! Print control
	39	USE wrk_nemo ! Memory Allocation
	40	USE timing ! Timing
	41	USE trabbl
	42	USE gridrandom
	43	USE dotprodfld
	44	USE tstool_tam
	45
	46	IMPLICIT NONE
	47	PRIVATE
	48
	49	PUBLIC tra_bbl_tan ! routine called by step.F90
	50	PUBLIC tra_bbl_init_tam ! routine called by opa.F90
	51	PUBLIC tra_bbl_dif_tan ! routine called by trcbbl.F90
	52	PUBLIC tra_bbl_adv_tan ! - - - -
	53	PUBLIC bbl_tan ! routine called by trcbbl.F90 and dtadyn.F90
	54	PUBLIC tra_bbl_adj ! routine called by step.F90
	55	PUBLIC tra_bbl_dif_adj ! routine called by trcbbl.F90
	56	PUBLIC tra_bbl_adv_adj ! - - - -
	57	PUBLIC bbl_adj ! routine called by trcbbl.F90 and dtadyn.F90
	58	PUBLIC tra_bbl_adj_tst ! routine called by tamtst
	59
	60	REAL(WP), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: utr_bbl_tl, vtr_bbl_tl
	61	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: ahu_bbl_tl, ahv_bbl_tl
	62	REAL(WP), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: utr_bbl_ad, vtr_bbl_ad
	63	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:), PUBLIC :: ahu_bbl_ad, ahv_bbl_ad
	64
	65	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ahu_bbl_0_tl, ahv_bbl_0_tl
	66	REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ahu_bbl_0_ad, ahv_bbl_0_ad
	67
	68	LOGICAL, PRIVATE :: ll_alloctl = .FALSE., ll_allocad = .FALSE.
	69
	70	!! * Substitutions
	71	# include "domzgr_substitute.h90"
	72	# include "vectopt_loop_substitute.h90"
	73	!!----------------------------------------------------------------------
	74	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
	75	!! $Id$
	76	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
	77	!!----------------------------------------------------------------------
	78	CONTAINS
	79
	80	INTEGER FUNCTION tra_bbl_alloc_tam( kmode )
	81	!!----------------------------------------------------------------------
	82	!! * FUNCTION tra_bbl_alloc_tam *
	83	!!----------------------------------------------------------------------
	84	INTEGER, OPTIONAL :: kmode
	85	INTEGER, DIMENSION(2) :: ierr
	86	INTEGER :: jmode
	87
	88	IF ( PRESENT( kmode ) ) THEN
	89	jmode = kmode
	90	ELSE
	91	jmode = 0
	92	END IF
	93	ierr(:) = 0.0_wp
	94
	95	IF ( ( jmode == 0 ) .OR. ( jmode == 1 ) .AND. ( .NOT. ll_alloctl ) ) THEN
	96	ALLOCATE( utr_bbl_tl (jpi,jpj), vtr_bbl_tl (jpi,jpj), &
	97	& ahu_bbl_tl (jpi,jpj), ahv_bbl_tl (jpi,jpj), &
	98	& STAT= ierr(1) )
	99	ll_alloctl = .TRUE.
	100	END IF
	101	!
	102	IF ( ( jmode == 0 ) .OR. ( jmode == 2 ) .AND. ( .NOT. ll_allocad ) ) THEN
	103	ALLOCATE( utr_bbl_ad (jpi,jpj), vtr_bbl_ad (jpi,jpj), &
	104	& ahu_bbl_ad (jpi,jpj), ahv_bbl_ad (jpi,jpj), &
	105	& STAT= ierr(2) )
	106	ll_allocad = .TRUE.
	107	END IF
	108	tra_bbl_alloc_tam = SUM(ierr)
	109	!
	110	IF( lk_mpp ) CALL mpp_sum ( tra_bbl_alloc_tam )
	111	IF( tra_bbl_alloc_tam > 0 ) CALL ctl_warn('tra_bbl_alloc_tam: allocation of arrays failed.')
	112	END FUNCTION tra_bbl_alloc_tam
	113
	114
	115	INTEGER FUNCTION tra_bbl_dealloc_tam( kmode )
	116	!!----------------------------------------------------------------------
	117	!! * FUNCTION tra_bbl_dealloc *
	118	!!----------------------------------------------------------------------
	119	INTEGER, OPTIONAL :: kmode
	120	INTEGER, DIMENSION(2) :: ierr
	121
	122	IF ( .NOT. PRESENT( kmode ) ) kmode=0
	123	ierr(:) = 0.0_wp
	124
	125	IF ( ( kmode == 0 ) .OR. ( kmode == 1 ) .AND. ( ll_alloctl ) ) THEN
	126	DEALLOCATE( utr_bbl_tl, vtr_bbl_tl, &
	127	& ahu_bbl_tl, ahv_bbl_tl, &
	128	& STAT= ierr(1) )
	129	ll_alloctl = .FALSE.
	130	END IF
	131	!
	132	IF ( ( kmode == 0 ) .OR. ( kmode == 1 ) .AND. ( ll_allocad ) ) THEN
	133	DEALLOCATE( utr_bbl_ad, vtr_bbl_ad, &
	134	& ahu_bbl_ad, ahv_bbl_ad, &
	135	& STAT= ierr(2) )
	136	ll_allocad = .FALSE.
	137	END IF
	138	tra_bbl_dealloc_tam = SUM(ierr)
	139	!
	140	IF( lk_mpp ) CALL mpp_sum ( tra_bbl_dealloc_tam )
	141	IF( tra_bbl_dealloc_tam > 0 ) CALL ctl_warn('tra_bbl_dealloc_tam: allocation of arrays failed.')
	142	END FUNCTION tra_bbl_dealloc_tam
	143
	144
	145	SUBROUTINE tra_bbl_tan( kt )
	146	!!----------------------------------------------------------------------
	147	!! * ROUTINE bbl_tan *
	148	!!
	149	!! ** Purpose : Compute the before tracer (t & s) trend associated
	150	!! with the bottom boundary layer and add it to the general
	151	!! trend of tracer equations.
	152	!!
	153	!! ** Method : Depending on namtra_bbl namelist parameters the bbl
	154	!! diffusive and/or advective contribution to the tracer trend
	155	!! is added to the general tracer trend
	156	!!----------------------------------------------------------------------
	157	INTEGER, INTENT( in ) :: kt ! ocean time-step
	158	!!
	159	REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdttl, ztrdstl
	160	!!----------------------------------------------------------------------
	161	!
	162	IF( nn_timing == 1 ) CALL timing_start( 'tra_bbl_tan')
	163	!
	164	IF( l_bbl ) CALL bbl_tan( kt, nit000, 'TRA' ) !* bbl coef. and transport (only if not already done in trcbbl)
	165
	166	IF( nn_bbl_ldf == 1 ) THEN !* Diffusive bbl
	167	!
	168	CALL tra_bbl_dif_tan( tsb, tsb_tl, tsa_tl, jpts )
	169	!
	170	END IF
	171
	172	IF( nn_bbl_adv /= 0 ) THEN !* Advective bbl
	173	!
	174	CALL tra_bbl_adv_tan( tsb, tsb_tl, tsa_tl, jpts )
	175	!
	176	END IF
	177	!
	178	IF( nn_timing == 1 ) CALL timing_stop( 'tra_bbl_tan')
	179	!
	180	END SUBROUTINE tra_bbl_tan
	181
	182
	183	SUBROUTINE tra_bbl_dif_tan( ptb, ptb_tl, pta_tl, kjpt )
	184	!!----------------------------------------------------------------------
	185	!! * ROUTINE tra_bbl_dif_tan *
	186	!!
	187	!! ** Purpose : Computes the bottom boundary horizontal and vertical
	188	!! advection terms.
	189	!!
	190	!! ** Method :
	191	!! * diffusive bbl (nn_bbl_ldf=1) :
	192	!! When the product grad( rho) * grad(h) < 0 (where grad is an
	193	!! along bottom slope gradient) an additional lateral 2nd order
	194	!! diffusion along the bottom slope is added to the general
	195	!! tracer trend, otherwise the additional trend is set to 0.
	196	!! A typical value of ahbt is 2000 m2/s (equivalent to
	197	!! a downslope velocity of 20 cm/s if the condition for slope
	198	!! convection is satified)
	199	!!
	200	!! ** Action : pta increased by the bbl diffusive trend
	201	!!
	202	!! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591.
	203	!! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430.
	204	!!----------------------------------------------------------------------
	205	!
	206	INTEGER , INTENT(in ) :: kjpt ! number of tracers
	207	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields
	208	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb_tl ! before and now tracer fields
	209	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta_tl ! tracer trend
	210	!
	211	INTEGER :: ji, jj, jn ! dummy loop indices
	212	INTEGER :: ik ! local integers
	213	REAL(wp) :: zbtr ! local scalars
	214	REAL(wp), POINTER, DIMENSION(:,:) :: zptb, zptbtl
	215	!!----------------------------------------------------------------------
	216	!
	217	IF( nn_timing == 1 ) CALL timing_start('tra_bbl_dif_tan')
	218	!
	219	CALL wrk_alloc( jpi, jpj, zptb, zptbtl )
	220	!
	221	DO jn = 1, kjpt ! tracer loop
	222	! ! ===========
	223	# if defined key_vectopt_loop
	224	DO jj = 1, 1 ! vector opt. (forced unrolling)
	225	DO ji = 1, jpij
	226	#else
	227	DO jj = 1, jpj
	228	DO ji = 1, jpi
	229	#endif
	230	ik = mbkt(ji,jj) ! bottom T-level index
	231	zptbtl(ji,jj) = ptb_tl(ji,jj,ik,jn) ! bottom before T and S
	232	zptb (ji,jj) = ptb (ji,jj,ik,jn) ! bottom before T and S
	233	END DO
	234	END DO
	235	! ! Compute the trend
	236	# if defined key_vectopt_loop
	237	DO jj = 1, 1 ! vector opt. (forced unrolling)
	238	DO ji = jpi+1, jpij-jpi-1
	239	# else
	240	DO jj = 2, jpjm1
	241	DO ji = 2, jpim1
	242	# endif
	243	ik = mbkt(ji,jj) ! bottom T-level index
	244	zbtr = r1_e1e2t(ji,jj) / fse3t(ji,jj,ik)
	245	pta_tl(ji,jj,ik,jn) = pta_tl(ji,jj,ik,jn) &
	246	# if defined control_param
	247	& + ( ahu_bbl_tl(ji ,jj ) * ( zptb(ji+1,jj ) - zptb(ji ,jj ) ) &
	248	& - ahu_bbl_tl(ji-1,jj ) * ( zptb(ji ,jj ) - zptb(ji-1,jj ) ) &
	249	& + ahv_bbl_tl(ji ,jj ) * ( zptb(ji ,jj+1) - zptb(ji ,jj ) ) &
	250	& - ahv_bbl_tl(ji ,jj-1) * ( zptb(ji ,jj ) - zptb(ji ,jj-1) ) ) * zbtr &
	251	#endif
	252	& + ( ahu_bbl(ji ,jj ) * ( zptbtl(ji+1,jj ) - zptbtl(ji ,jj ) ) &
	253	& - ahu_bbl(ji-1,jj ) * ( zptbtl(ji ,jj ) - zptbtl(ji-1,jj ) ) &
	254	& + ahv_bbl(ji ,jj ) * ( zptbtl(ji ,jj+1) - zptbtl(ji ,jj ) ) &
	255	& - ahv_bbl(ji ,jj-1) * ( zptbtl(ji ,jj ) - zptbtl(ji ,jj-1) ) ) * zbtr
	256	END DO
	257	END DO
	258	! ! ===========
	259	END DO ! end tracer
	260	! ! ===========
	261	CALL wrk_dealloc( jpi, jpj, zptbtl, zptb )
	262	!
	263	IF( nn_timing == 1 ) CALL timing_stop('tra_bbl_dif_tan')
	264	!
	265	END SUBROUTINE tra_bbl_dif_tan
	266
	267
	268	SUBROUTINE tra_bbl_adv_tan( ptb, ptb_tl, pta_tl, kjpt )
	269	!!----------------------------------------------------------------------
	270	!! * ROUTINE trc_bbl *
	271	!!
	272	!! ** Purpose : Compute the before passive tracer trend associated
	273	!! with the bottom boundary layer and add it to the general trend
	274	!! of tracer equations.
	275	!! ** Method : advective bbl (nn_bbl_adv = 1 or 2) :
	276	!! nn_bbl_adv = 1 use of the ocean near bottom velocity as bbl velocity
	277	!! nn_bbl_adv = 2 follow Campin and Goosse (1999) implentation i.e.
	278	!! transport proportional to the along-slope density gradient
	279	!!
	280	!! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591.
	281	!! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430.
	282	!!----------------------------------------------------------------------
	283	INTEGER , INTENT(in ) :: kjpt ! number of tracers
	284	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields
	285	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb_tl ! before and now tangent tracer fields
	286	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta_tl ! tracer trend
	287	!
	288	INTEGER :: ji, jj, jk, jn ! dummy loop indices
	289	INTEGER :: iis , iid , ijs , ijd ! local integers
	290	INTEGER :: ikus, ikud, ikvs, ikvd ! - -
	291	REAL(wp) :: zbtr, ztra ! local scalars
	292	REAL(wp) :: ztratl ! - -
	293	REAL(wp) :: zu_bbl, zv_bbl ! - -
	294	REAL(wp) :: zu_bbltl, zv_bbltl ! - -
	295	!!----------------------------------------------------------------------
	296	!
	297	IF( nn_timing == 1 ) CALL timing_start( 'tra_bbl_adv_tan')
	298	! ! ===========
	299	DO jn = 1, kjpt ! tracer loop
	300	! ! ===========
	301	# if defined key_vectopt_loop
	302	DO jj = 1, 1
	303	DO ji = 1, jpij-jpi-1 ! vector opt. (forced unrolling)
	304	# else
	305	DO jj = 1, jpjm1
	306	DO ji = 1, jpim1 ! CAUTION start from i=1 to update i=2 when cyclic east-west
	307	# endif
	308	IF( utr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero i-direction bbl advection
	309	! down-slope i/k-indices (deep) & up-slope i/k indices (shelf)
	310	iid = ji + MAX( 0, mgrhu(ji,jj) ) ; iis = ji + 1 - MAX( 0, mgrhu(ji,jj) )
	311	ikud = mbku_d(ji,jj) ; ikus = mbku(ji,jj)
	312	zu_bbl = ABS( utr_bbl(ji,jj) )
	313	zu_bbltl = SIGN( utr_bbl_tl(ji,jj), utr_bbl(ji,jj) )
	314	!
	315	! ! up -slope T-point (shelf bottom point)
	316	zbtr = r1_e1e2t(iis,jj) / fse3t(iis,jj,ikus)
	317	ztratl = ( zu_bbltl * ( ptb (iid,jj,ikus,jn) - ptb (iis,jj,ikus,jn) ) &
	318	& + zu_bbl * ( ptb_tl(iid,jj,ikus,jn) - ptb_tl(iis,jj,ikus,jn) ) ) * zbtr
	319	pta_tl(iis,jj,ikus,jn) = pta_tl(iis,jj,ikus,jn) + ztratl
	320	!
	321	DO jk = ikus, ikud-1 ! down-slope upper to down T-point (deep column)
	322	zbtr = r1_e1e2t(iid,jj) / fse3t(iid,jj,jk)
	323	ztratl = ( zu_bbltl * ( ptb (iid,jj,jk+1,jn) - ptb (iid,jj,jk,jn) ) &
	324	& + zu_bbl * ( ptb_tl(iid,jj,jk+1,jn) - ptb_tl(iid,jj,jk,jn) ) ) * zbtr
	325	pta_tl(iid,jj,jk,jn) = pta_tl(iid,jj,jk,jn) + ztratl
	326	END DO
	327	!
	328	zbtr = r1_e1e2t(iid,jj) / fse3t(iid,jj,ikud)
	329	ztratl = ( zu_bbltl * ( ptb (iis,jj,ikus,jn) - ptb (iid,jj,ikud,jn) ) &
	330	& + zu_bbl * ( ptb_tl(iis,jj,ikus,jn) - ptb_tl(iid,jj,ikud,jn) ) ) * zbtr
	331	pta_tl(iid,jj,ikud,jn) = pta_tl(iid,jj,ikud,jn) + ztratl
	332	ENDIF
	333	!
	334	IF( vtr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero j-direction bbl advection
	335	! down-slope j/k-indices (deep) & up-slope j/k indices (shelf)
	336	ijd = jj + MAX( 0, mgrhv(ji,jj) ) ; ijs = jj + 1 - MAX( 0, mgrhv(ji,jj) )
	337	ikvd = mbkv_d(ji,jj) ; ikvs = mbkv(ji,jj)
	338	zv_bbl = ABS( vtr_bbl(ji,jj) )
	339	zv_bbltl = SIGN( vtr_bbl_tl(ji,jj), vtr_bbl(ji,jj) )
	340	!
	341	! up -slope T-point (shelf bottom point)
	342	zbtr = r1_e1e2t(ji,ijs) / fse3t(ji,ijs,ikvs)
	343	ztratl = ( zv_bbltl * ( ptb (ji,ijd,ikvs,jn) - ptb (ji,ijs,ikvs,jn) ) &
	344	& + zv_bbl * ( ptb_tl(ji,ijd,ikvs,jn) - ptb_tl(ji,ijs,ikvs,jn) ) ) * zbtr
	345	pta_tl(ji,ijs,ikvs,jn) = pta_tl(ji,ijs,ikvs,jn) + ztratl
	346	!
	347	DO jk = ikvs, ikvd-1 ! down-slope upper to down T-point (deep column)
	348	zbtr = r1_e1e2t(ji,ijd) / fse3t(ji,ijd,jk)
	349	ztratl = ( zv_bbltl * ( ptb (ji,ijd,jk+1,jn) - ptb (ji,ijd,jk,jn) ) &
	350	& + zv_bbl * ( ptb_tl(ji,ijd,jk+1,jn) - ptb_tl(ji,ijd,jk,jn) ) ) * zbtr
	351	pta_tl(ji,ijd,jk,jn) = pta_tl(ji,ijd,jk,jn) + ztratl
	352	END DO
	353	! ! down-slope T-point (deep bottom point)
	354	zbtr = r1_e1e2t(ji,ijd) / fse3t(ji,ijd,ikvd)
	355	ztratl = ( zv_bbltl * ( ptb (ji,ijs,ikvs,jn) - ptb (ji,ijd,ikvd,jn) ) &
	356	& + zv_bbl * ( ptb_tl(ji,ijs,ikvs,jn) - ptb_tl(ji,ijd,ikvd,jn) ) ) * zbtr
	357	pta_tl(ji,ijd,ikvd,jn) = pta_tl(ji,ijd,ikvd,jn) + ztratl
	358	ENDIF
	359	END DO
	360	!
	361	END DO
	362	! ! ===========
	363	END DO ! end tracer
	364	! ! ===========
	365	!
	366	IF( nn_timing == 1 ) CALL timing_stop( 'tra_bbl_adv_tan')
	367	!
	368	END SUBROUTINE tra_bbl_adv_tan
	369
	370
	371	SUBROUTINE bbl_tan( kt, kit000, cdtype )
	372	!!----------------------------------------------------------------------
	373	!! * ROUTINE bbl *
	374	!!
	375	!! ** Purpose : Computes the bottom boundary horizontal and vertical
	376	!! advection terms.
	377	!!
	378	!! ** Method :
	379	!! * diffusive bbl (nn_bbl_ldf=1) :
	380	!! When the product grad( rho) * grad(h) < 0 (where grad is an
	381	!! along bottom slope gradient) an additional lateral 2nd order
	382	!! diffusion along the bottom slope is added to the general
	383	!! tracer trend, otherwise the additional trend is set to 0.
	384	!! A typical value of ahbt is 2000 m2/s (equivalent to
	385	!! a downslope velocity of 20 cm/s if the condition for slope
	386	!! convection is satified)
	387	!! * advective bbl (nn_bbl_adv=1 or 2) :
	388	!! nn_bbl_adv = 1 use of the ocean velocity as bbl velocity
	389	!! nn_bbl_adv = 2 follow Campin and Goosse (1999) implentation
	390	!! i.e. transport proportional to the along-slope density gradient
	391	!!
	392	!! NB: the along slope density gradient is evaluated using the
	393	!! local density (i.e. referenced at a common local depth).
	394	!!
	395	!! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591.
	396	!! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430.
	397	!!----------------------------------------------------------------------
	398	!
	399	INTEGER , INTENT(in ) :: kt ! ocean time-step index
	400	INTEGER , INTENT(in ) :: kit000 ! first time step index
	401	CHARACTER(len=3), INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator)
	402	!!
	403	INTEGER :: ji, jj ! dummy loop indices
	404	INTEGER :: ik ! local integers
	405	INTEGER :: iis , iid , ijs , ijd ! - -
	406	INTEGER :: ikus, ikud, ikvs, ikvd ! - -
	407	REAL(wp) :: zsign, zsigna, zgbbl ! local scalars
	408	REAL(wp) :: zgdrho, zt, zs, zh ! - -
	409	REAL(wp) :: zgdrhotl, zttl, zstl, zhtl! - -
	410	!!
	411	REAL(wp) :: fsalbt, fsbeta, pft, pfs, pfh ! statement function
	412	REAL(wp) :: fsalbt_tan, fsbeta_tan, pfttl, pfstl, pfhtl ! statement function
	413	REAL(wp), POINTER, DIMENSION(:,:) :: zub , zvb , ztb , zsb , zdep
	414	REAL(wp), POINTER, DIMENSION(:,:) :: zubtl, zvbtl, ztbtl, zsbtl
	415	!!----------------------- zv_bbl-----------------------------------------------
	416	! ratio alpha/beta = fsalbt : ratio of thermal over saline expension coefficients
	417	! ================ pft : potential temperature in degrees celcius
	418	! pfs : salinity anomaly (s-35) in psu
	419	! pfh : depth in meters
	420	! nn_eos = 0 (Jackett and McDougall 1994 formulation)
	421	fsalbt( pft, pfs, pfh ) = & ! alpha/beta
	422	( ( ( -0.255019e-07 * pft + 0.298357e-05 ) * pft &
	423	- 0.203814e-03 ) * pft &
	424	+ 0.170907e-01 ) * pft &
	425	+ 0.665157e-01 &
	426	+(-0.678662e-05 * pfs - 0.846960e-04 * pft + 0.378110e-02 ) * pfs &
	427	+ ( ( - 0.302285e-13 * pfh &
	428	- 0.251520e-11 * pfs &
	429	+ 0.512857e-12 * pft * pft ) * pfh &
	430	- 0.164759e-06 * pfs &
	431	+( 0.791325e-08 * pft - 0.933746e-06 ) * pft &
	432	+ 0.380374e-04 ) * pfh
	433	fsbeta( pft, pfs, pfh ) = & ! beta
	434	( ( -0.415613e-09 * pft + 0.555579e-07 ) * pft &
	435	- 0.301985e-05 ) * pft &
	436	+ 0.785567e-03 &
	437	+ ( 0.515032e-08 * pfs &
	438	+ 0.788212e-08 * pft - 0.356603e-06 ) * pfs &
	439	+( ( 0.121551e-17 * pfh &
	440	- 0.602281e-15 * pfs &
	441	- 0.175379e-14 * pft + 0.176621e-12 ) * pfh &
	442	+ 0.408195e-10 * pfs &
	443	+ ( - 0.213127e-11 * pft + 0.192867e-09 ) * pft &
	444	- 0.121555e-07 ) * pfh
	445
	446	fsalbt_tan( pft, pfs, pfh, pfttl, pfstl, pfhtl ) = & ! alpha/beta
	447	& ( - 0.255019e-07 * 4 * pft * pft * pft &
	448	& + 0.298357e-05 * 3 * pft * pft &
	449	& - 0.203814e-03 * 2 * pft &
	450	& - 0.846960e-04 * pfs &
	451	& + 0.512857e-12 * 2 * pft * pfh * pfh &
	452	& + 0.791325e-08 * pft * pfh &
	453	& - 0.933746e-06 * pfh &
	454	& + 0.170907e-01 ) * pfttl &
	455	& + ( - 0.678662e-05 * 2 * pfs &
	456	& - 0.846960e-04 * pft &
	457	& - 0.251520e-11 * pfh * pfh &
	458	& - 0.164759e-06 * pfh &
	459	& + 0.378110e-02 ) * pfstl &
	460	& + ( - 0.302285e-13 * 3 * pfh * pfh &
	461	& - 0.251520e-11 * pfs * pfh &
	462	& + 0.512857e-12 * pft * pft * pfh &
	463	& - 0.164759e-06 * pfs &
	464	& + 0.791325e-08 * pft * pft &
	465	& - 0.933746e-06 * pft &
	466	& + 0.380374e-04 ) * pfhtl
	467
	468
	469	fsbeta_tan( pft, pfs, pfh, pfttl, pfstl, pfhtl ) = & ! beta
	470	& ( - 0.415613e-09 * 3 * pft * pft &
	471	& + 0.555579e-07 * 2 * pft &
	472	& - 0.301985e-05 &
	473	& + 0.788212e-08 * pfs &
	474	& - 0.213127e-11 * 2 * pfh * pft &
	475	& - 0.175379e-14 * pfh * pfh ) * pfttl &
	476	& + ( 0.788212e-08 * pft &
	477	& + 0.515032e-08 * 2 * pfs &
	478	& - 0.356603e-06 &
	479	& + 0.408195e-10 * pfh &
	480	& - 0.602281e-15 * pfh * pfh ) * pfstl &
	481	& + ( 0.121551e-17 * 3 * pfh * pfh &
	482	& - 0.602281e-15 * 2 * pfs * pfh &
	483	& - 0.175379e-14 * 2 * pft * pfh &
	484	& + 0.176621e-12 * 2 * pfh &
	485	& + 0.408195e-10 * pfs &
	486	& + 0.192867e-09 * pfh &
	487	& - 0.213127e-11 * pft * pft &
	488	& + 0.192867e-09 * pft &
	489	& - 0.121555e-07 ) * pfhtl
	490
	491
	492	!!----------------------------------------------------------------------
	493
	494	!
	495	IF( nn_timing == 1 ) CALL timing_start( 'bbl_tan')
	496	!
	497	CALL wrk_alloc( jpi, jpj, zub , zvb , ztb , zsb , zdep, &
	498	& zubtl, zvbtl, ztbtl, zsbtl )
	499	!
	500
	501	IF( kt == kit000 ) THEN
	502	IF(lwp) WRITE(numout,*)
	503	IF(lwp) WRITE(numout,*) 'trabbl_tam:bbl_tan : Compute bbl velocities and diffusive coefficients in ', cdtype
	504	IF(lwp) WRITE(numout,*) '~~~~~~~~~~'
	505	ENDIF
	506
	507	! !* bottom temperature, salinity, velocity and depth
	508	#if defined key_vectopt_loop
	509	DO jj = 1, 1 ! vector opt. (forced unrolling)
	510	DO ji = 1, jpij
	511	#else
	512	DO jj = 1, jpj
	513	DO ji = 1, jpi
	514	#endif
	515	ik = mbkt(ji,jj) ! bottom T-level index
	516	ztb (ji,jj) = tsb(ji,jj,ik,jp_tem) * tmask(ji,jj,1) ! bottom before T and S
	517	zsb (ji,jj) = tsb(ji,jj,ik,jp_sal) * tmask(ji,jj,1)
	518	ztbtl(ji,jj) = tsb_tl(ji,jj,ik,jp_tem) * tmask(ji,jj,1) ! bottom before T and S
	519	zsbtl(ji,jj) = tsb_tl(ji,jj,ik,jp_sal) * tmask(ji,jj,1)
	520	zdep(ji,jj) = fsdept_0(ji,jj,ik) ! bottom T-level reference depth
	521	!
	522	zub(ji,jj) = un(ji,jj,mbku(ji,jj)) ! bottom velocity
	523	zvb(ji,jj) = vn(ji,jj,mbkv(ji,jj))
	524	END DO
	525	END DO
	526
	527	! !-------------------!
	528	IF( nn_bbl_ldf == 1 ) THEN ! diffusive bbl !
	529	! !-------------------!
	530	! AV NOTE : while rn_ahtbbl remains a passive variable, the code below will only yield ah_bbl_tl=0, so i put it under key
	531	#if defined key_control_param
	532	DO jj = 1, jpjm1 ! (criteria for non zero flux: grad(rho).grad(h) < 0 )
	533	DO ji = 1, jpim1
	534	! ! i-direction
	535	zt = 0.5 * ( ztb (ji,jj) + ztb (ji+1,jj) ) ! T, S anomalie, and depth
	536	zs = 0.5 * ( zsb (ji,jj) + zsb (ji+1,jj) ) - 35.0
	537	zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) )
	538	! ! masked bbl i-gradient of density
	539	zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji+1,jj) - ztb(ji,jj) ) &
	540	& - ( zsb(ji+1,jj) - zsb(ji,jj) ) ) * umask(ji,jj,1)
	541	!
	542	zsign = SIGN( 0.5, - zgdrho * REAL( mgrhu(ji,jj) ) ) ! sign of ( i-gradient * i-slope )
	543	ahu_bbl_tl(ji,jj) = ( 0.5 - zsign ) * ahu_bbl_0_tl(ji,jj) ! masked diffusive flux coeff.
	544	!
	545	! ! j-direction
	546	zt = 0.5 * ( ztb (ji,jj+1) + ztb (ji,jj) ) ! T, S anomalie, and depth
	547	zs = 0.5 * ( zsb (ji,jj+1) + zsb (ji,jj) ) - 35.0
	548	zh = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) )
	549	! ! masked bbl j-gradient of density
	550	zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji,jj+1) - ztb(ji,jj) ) &
	551	& - ( zsb(ji,jj+1) - zsb(ji,jj) ) ) * vmask(ji,jj,1)
	552	!
	553	zsign = SIGN( 0.5, -zgdrho * REAL( mgrhv(ji,jj) ) ) ! sign of ( j-gradient * j-slope )
	554	ahv_bbl_tl(ji,jj) = ( 0.5 - zsign ) * ahv_bbl_0_tl(ji,jj)
	555	!
	556	END DO
	557	END DO
	558	#else
	559	DO jj = 1, jpjm1
	560	DO ji = 1, jpim1
	561	ahu_bbl_tl(ji,jj)=0.0_wp
	562	ahv_bbl_tl(ji,jj)=0.0_wp
	563	END DO
	564	END DO
	565	#endif
	566	!
	567	ENDIF
	568
	569	! !-------------------!
	570	IF( nn_bbl_adv /= 0 ) THEN ! advective bbl !
	571	! !-------------------!
	572	SELECT CASE ( nn_bbl_adv ) !* bbl transport type
	573	!
	574	CASE( 1 ) != use of upper velocity
	575	! AV NOTE: not much needed for deriving, almost all the computations are for the SIGN, which is kept identical as in the NL
	576	DO jj = 1, jpjm1 ! criteria: grad(rho).grad(h)<0 and grad(rho).grad(h)<0
	577	DO ji = 1, fs_jpim1 ! vector opt.
	578	! ! i-direction
	579	zt = 0.5 * ( ztb (ji,jj) + ztb (ji+1,jj) ) ! T, S anomalie, and depth
	580	zs = 0.5 * ( zsb (ji,jj) + zsb (ji+1,jj) ) - 35.0
	581	zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) )
	582	! ! masked bbl i-gradient of density
	583	zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji+1,jj) - ztb(ji,jj) ) &
	584	& - ( zsb(ji+1,jj) - zsb(ji,jj) ) ) * umask(ji,jj,1)
	585	!
	586	zsign = SIGN( 0.5, - zgdrho * REAL( mgrhu(ji,jj) ) ) ! sign of i-gradient * i-slope
	587	zsigna= SIGN( 0.5, zub(ji,jj) * REAL( mgrhu(ji,jj) ) ) ! sign of u * i-slope
	588	!
	589	! ! bbl velocity
	590	utr_bbl_tl(ji,jj) = ( 0.5 + zsigna ) * ( 0.5 - zsign ) * e2u(ji,jj) * e3u_bbl_0(ji,jj) * zubtl(ji,jj)
	591	!
	592	! ! j-direction
	593	zt = 0.5 * ( ztb (ji,jj+1) + ztb (ji,jj) ) ! T, S anomalie, and depth
	594	zs = 0.5 * ( zsb (ji,jj+1) + zsb (ji,jj) ) - 35.0
	595	zh = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) )
	596	! ! masked bbl j-gradient of density
	597	zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji,jj+1) - ztb(ji,jj) ) &
	598	& - ( zsb(ji,jj+1) - zsb(ji,jj) ) ) * vmask(ji,jj,1)
	599	zsign = SIGN( 0.5, - zgdrho * REAL( mgrhv(ji,jj) ) ) ! sign of j-gradient * j-slope
	600	zsigna= SIGN( 0.5, zvb(ji,jj) * REAL( mgrhv(ji,jj) ) ) ! sign of u * i-slope
	601	!
	602	! ! bbl velocity
	603	vtr_bbl_tl(ji,jj) = ( 0.5 + zsigna ) * ( 0.5 - zsign ) * e1v(ji,jj) * e3v_bbl_0(ji,jj) * zvbtl(ji,jj)
	604	END DO
	605	END DO
	606	!
	607	CASE( 2 ) != bbl velocity = F( delta rho )
	608	! AV NOTE: this one is nastier
	609	zgbbl = grav * rn_gambbl
	610	DO jj = 1, jpjm1 ! criteria: rho_up > rho_down
	611	DO ji = 1, fs_jpim1 ! vector opt.
	612	! ! i-direction
	613	! down-slope T-point i/k-index (deep) & up-slope T-point i/k-index (shelf)
	614	iid = ji + MAX( 0, mgrhu(ji,jj) ) ; iis = ji + 1 - MAX( 0, mgrhu(ji,jj) )
	615	ikud = mbku_d(ji,jj) ; ikus = mbku(ji,jj)
	616	!
	617	! ! mid-depth density anomalie (up-slope minus down-slope)
	618	zt = 0.5 * ( ztb (ji,jj) + ztb (ji+1,jj) ) ! mid slope depth of T, S, and depth
	619	zs = 0.5 * ( zsb (ji,jj) + zsb (ji+1,jj) ) - 35.0
	620	zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) )
	621	zgdrho = fsbeta( zt, zs, zh ) &
	622	& * ( fsalbt( zt, zs, zh ) * ( ztb(iid,jj) - ztb(iis,jj) ) &
	623	& - ( zsb(iid,jj) - zsb(iis,jj) ) ) * umask(ji,jj,1)
	624	zttl = 0.5 * ( ztbtl (ji,jj) + ztbtl (ji+1,jj) ) ! mid slope depth of T, S, and depth
	625	zstl = 0.5 * ( zsbtl (ji,jj) + zsbtl (ji+1,jj) )
	626	zhtl = 0.0_wp
	627	zgdrhotl = ( fsbeta_tan( zt, zs, zh, zttl, zstl, zhtl ) &
	628	& * ( fsalbt( zt, zs, zh ) * ( ztb(iid,jj) - ztb(iis,jj) ) &
	629	& - ( zsb(iid,jj) - zsb(iis,jj) ) ) &
	630	& + fsbeta( zt, zs, zh ) &
	631	& * ( fsalbt_tan( zt, zs, zh, zttl, zstl, zhtl ) &
	632	& * ( ztb (iid,jj) - ztb (iis,jj) ) &
	633	& + fsalbt ( zt, zs, zh ) * ( ztbtl(iid,jj) - ztbtl(iis,jj) ) &
	634	& - ( zsbtl(iid,jj) - zsbtl(iis,jj) ) ) ) * umask(ji,jj,1)
	635
	636	zsign = SIGN( 0.5_wp, zgdrho ) ! tangent of zgdrho = MAX( 0.e0, zgdrho )
	637	! ! bbl transport (down-slope direction)
	638	utr_bbl_tl(ji,jj) = zsign * e2u(ji,jj) * e3u_bbl_0(ji,jj) * zgbbl * zgdrhotl * REAL( mgrhu(ji,jj) )
	639	!
	640	! ! j-direction
	641	! down-slope T-point j/k-index (deep) & of the up -slope T-point j/k-index (shelf)
	642	ijd = jj + MAX( 0, mgrhv(ji,jj) ) ; ijs = jj + 1 - MAX( 0, mgrhv(ji,jj) )
	643	ikvd = mbkv_d(ji,jj) ; ikvs = mbkv(ji,jj)
	644	!
	645	! ! mid-depth density anomalie (up-slope minus down-slope)
	646	zt = 0.5 * ( ztb (ji,jj) + ztb (ji,jj+1) ) ! mid slope depth of T, S, and depth
	647	zs = 0.5 * ( zsb (ji,jj) + zsb (ji,jj+1) ) - 35.0
	648	zh = 0.5 * ( zdep(ji,jj) + zdep(ji,jj+1) )
	649	zgdrho = fsbeta( zt, zs, zh ) &
	650	& * ( fsalbt( zt, zs, zh ) * ( ztb(ji,ijd) - ztb(ji,ijs) ) &
	651	& - ( zsb(ji,ijd) - zsb(ji,ijs) ) ) * vmask(ji,jj,1)
	652	zttl = 0.5 * ( ztbtl (ji,jj) + ztbtl (ji,jj+1) ) ! mid slope depth of T, S, and depth
	653	zstl = 0.5 * ( zsbtl (ji,jj) + zsbtl (ji,jj+1) )
	654	zhtl = 0.0_wp
	655	zgdrhotl = ( fsbeta_tan( zt, zs, zh, zttl, zstl, zhtl ) &
	656	& * ( fsalbt( zt, zs, zh ) * ( ztb(ji,ijd) - ztb(ji,ijs) ) &
	657	& - ( zsb(ji,ijd) - zsb(ji,ijs) ) ) &
	658	& + fsbeta( zt, zs, zh ) &
	659	& * ( fsalbt_tan( zt, zs, zh, zttl, zstl, zhtl ) &
	660	& * ( ztb (ji,ijd) - ztb (ji,ijs) ) &
	661	& + fsalbt ( zt, zs, zh ) * ( ztbtl(ji,ijd) - ztbtl(ji,ijs) ) &
	662	& - ( zsbtl(ji,ijd) - zsbtl(ji,ijs) ) ) ) * vmask(ji,jj,1)
	663	!
	664	zsign = SIGN( 0.5_wp, zgdrho ) ! tangent of zgdrho = MAX( 0.e0, zgdrho )
	665	! ! bbl transport (down-slope direction)
	666	vtr_bbl_tl(ji,jj) = zsign * e1v(ji,jj) * e3v_bbl_0(ji,jj) * zgbbl * zgdrhotl * REAL( mgrhv(ji,jj) )
	667	END DO
	668	END DO
	669	END SELECT
	670	!
	671	ENDIF
	672	!
	673	CALL wrk_dealloc( jpi, jpj, zub , zvb , ztb , zsb , zdep, &
	674	& zubtl, zvbtl, ztbtl, zsbtl )
	675	!
	676	IF( nn_timing == 1 ) CALL timing_stop( 'bbl_tan')
	677	!
	678	END SUBROUTINE bbl_tan
	679
	680
	681	SUBROUTINE tra_bbl_adj( kt )
	682	!!----------------------------------------------------------------------
	683	!! * ROUTINE bbl_adj *
	684	!!
	685	!! ** Purpose : Compute the before tracer (t & s) trend associated
	686	!! with the bottom boundary layer and add it to the general
	687	!! trend of tracer equations.
	688	!!
	689	!! ** Method : Depending on namtra_bbl namelist parameters the bbl
	690	!! diffusive and/or advective contribution to the tracer trend
	691	!! is added to the general tracer trend
	692	!!----------------------------------------------------------------------
	693	INTEGER, INTENT( in ) :: kt ! ocean time-step
	694	!!----------------------------------------------------------------------
	695	!
	696	IF( nn_timing == 1 ) CALL timing_start( 'tra_bbl_adj')
	697	!
	698	IF( l_bbl ) CALL bbl_adj( kt, nit000, 'TRA' ) !* bbl coef. and transport (only if not already done in trcbbl)
	699
	700	IF( nn_bbl_ldf == 1 ) THEN !* Diffusive bbl
	701	!
	702	CALL tra_bbl_dif_adj( tsb, tsb_ad, tsa_ad, jpts )
	703	!
	704	END IF
	705
	706	IF( nn_bbl_adv /= 0 ) THEN !* Advective bbl
	707	!
	708	CALL tra_bbl_adv_adj( tsb, tsb_ad, tsa_ad, jpts )
	709	!
	710	END IF
	711	!
	712	IF( nn_timing == 1 ) CALL timing_stop( 'tra_bbl_adj')
	713	!
	714	END SUBROUTINE tra_bbl_adj
	715
	716
	717	SUBROUTINE tra_bbl_dif_adj( ptb, ptb_ad, pta_ad, kjpt )
	718	!!----------------------------------------------------------------------
	719	!! * ROUTINE tra_bbl_dif_adj *
	720	!!
	721	!! ** Purpose : Computes the bottom boundary horizontal and vertical
	722	!! advection terms.
	723	!!
	724	!! ** Method :
	725	!! * diffusive bbl (nn_bbl_ldf=1) :
	726	!! When the product grad( rho) * grad(h) < 0 (where grad is an
	727	!! along bottom slope gradient) an additional lateral 2nd order
	728	!! diffusion along the bottom slope is added to the general
	729	!! tracer trend, otherwise the additional trend is set to 0.
	730	!! A typical value of ahbt is 2000 m2/s (equivalent to
	731	!! a downslope velocity of 20 cm/s if the condition for slope
	732	!! convection is satified)
	733	!!
	734	!! ** Action : pta increased by the bbl diffusive trend
	735	!!
	736	!! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591.
	737	!! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430.
	738	!!----------------------------------------------------------------------
	739	!
	740	INTEGER , INTENT(in ) :: kjpt ! number of tracers
	741	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields
	742	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: ptb_ad ! before and now tracer fields
	743	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta_ad ! tracer trend
	744	!
	745	INTEGER :: ji, jj, jn ! dummy loop indices
	746	INTEGER :: ik ! local integers
	747	REAL(wp) :: zbtr ! local scalars
	748	REAL(wp), POINTER, DIMENSION(:,:) :: zptb, zptbad
	749	!!----------------------------------------------------------------------
	750	!
	751	IF( nn_timing == 1 ) CALL timing_start('tra_bbl_dif_adj')
	752	!
	753	CALL wrk_alloc( jpi, jpj, zptb, zptbad )
	754	zptbad(:,:) = 0.0_wp
	755	!
	756	DO jn = 1, kjpt ! tracer loop
	757	! ! ===========
	758	# if defined key_vectopt_loop
	759	DO jj = 1, 1 ! vector opt. (forced unrolling)
	760	DO ji = 1, jpij
	761	# else
	762	DO jj = 1, jpj
	763	DO ji = 1, jpi
	764	#endif
	765	ik = mbkt(ji,jj) ! bottom T-level index
	766	zptb (ji,jj) = ptb (ji,jj,ik,jn) ! bottom before T and S
	767	END DO
	768	END DO
	769	! ! ===========
	770	! ! Compute the trend
	771	# if defined key_vectopt_loop
	772	DO jj = 1, 1 ! vector opt. (forced unrolling)
	773	DO ji = jpi+1, jpij-jpi-1
	774	# else
	775	DO jj = jpjm1, 2, -1
	776	DO ji = jpim1, 2, -1
	777	# endif
	778	ik = mbkt(ji,jj) ! bottom T-level index
	779	zbtr = r1_e1e2t(ji,jj) / fse3t(ji,jj,ik)
	780	# if defined control_param
	781	ahu_bbl_ad(ji ,jj ) = ahu_bbl_ad(ji ,jj ) + pta_ad(ji,jj,ik,jn) * ( zptb(ji+1,jj ) - zptb(ji ,jj ) ) * zbtr
	782	ahu_bbl_ad(ji-1,jj ) = ahu_bbl_ad(ji-1,jj ) - pta_ad(ji,jj,ik,jn) * ( zptb(ji ,jj ) - zptb(ji-1,jj ) ) * zbtr
	783	ahv_bbl_ad(ji ,jj ) = ahv_bbl_ad(ji ,jj ) + pta_ad(ji,jj,ik,jn) * ( zptb(ji ,jj+1) - zptb(ji ,jj ) ) * zbtr
	784	ahv_bbl_ad(ji ,jj-1) = ahv_bbl_ad(ji ,jj-1) - pta_ad(ji,jj,ik,jn) * ( zptb(ji ,jj ) - zptb(ji ,jj-1) ) * zbtr
	785	# endif
	786	zptbad(ji ,jj ) = zptbad(ji ,jj ) - pta_ad(ji,jj,ik,jn) * ( ahu_bbl(ji ,jj ) + ahu_bbl(ji-1,jj ) &
	787	& + ahv_bbl(ji ,jj ) + ahv_bbl(ji ,jj-1) ) * zbtr
	788	zptbad(ji+1,jj ) = zptbad(ji+1,jj ) + pta_ad(ji,jj,ik,jn) * ahu_bbl(ji ,jj ) * zbtr
	789	zptbad(ji-1,jj ) = zptbad(ji-1,jj ) + pta_ad(ji,jj,ik,jn) * ahu_bbl(ji-1,jj ) * zbtr
	790	zptbad(ji ,jj+1) = zptbad(ji ,jj+1) + pta_ad(ji,jj,ik,jn) * ahv_bbl(ji ,jj ) * zbtr
	791	zptbad(ji ,jj-1) = zptbad(ji ,jj-1) + pta_ad(ji,jj,ik,jn) * ahv_bbl(ji ,jj-1) * zbtr
	792
	793	pta_ad(ji,jj,ik,jn) = pta_ad(ji,jj,ik,jn) &
	794	& + ( ahu_bbl(ji ,jj ) * ( zptbad(ji+1,jj ) - zptbad(ji ,jj ) ) &
	795	& - ahu_bbl(ji-1,jj ) * ( zptbad(ji ,jj ) - zptbad(ji-1,jj ) ) &
	796	& + ahv_bbl(ji ,jj ) * ( zptbad(ji ,jj+1) - zptbad(ji ,jj ) ) &
	797	& - ahv_bbl(ji ,jj-1) * ( zptbad(ji ,jj ) - zptbad(ji ,jj-1) ) ) * zbtr
	798	END DO
	799	END DO
	800	# if defined key_vectopt_loop
	801	DO jj = 1, 1 ! vector opt. (forced unrolling)
	802	DO ji = 1, jpij
	803	#else
	804	DO jj = 1, jpj
	805	DO ji = 1, jpi
	806	#endif
	807	ik = mbkt(ji,jj) ! bottom T-level index
	808	ptb_ad(ji,jj,ik,jn) = ptb_ad(ji,jj,ik,jn) + zptbad(ji,jj)
	809	zptbad(ji,jj) = 0.0_wp ! bottom before T and S
	810	END DO
	811	END DO
	812	! ! ===========
	813	END DO ! end tracer
	814	! ! ===========
	815	CALL wrk_dealloc( jpi, jpj, zptbad, zptb )
	816	!
	817	IF( nn_timing == 1 ) CALL timing_stop('tra_bbl_dif_adj')
	818	!
	819	END SUBROUTINE tra_bbl_dif_adj
	820
	821
	822	SUBROUTINE tra_bbl_adv_adj( ptb, ptb_ad, pta_ad, kjpt )
	823	!!----------------------------------------------------------------------
	824	!! * ROUTINE trc_bbl *
	825	!!
	826	!! ** Purpose : Compute the before passive tracer trend associated
	827	!! with the bottom boundary layer and add it to the general trend
	828	!! of tracer equations.
	829	!! ** Method : advective bbl (nn_bbl_adv = 1 or 2) :
	830	!! nn_bbl_adv = 1 use of the ocean near bottom velocity as bbl velocity
	831	!! nn_bbl_adv = 2 follow Campin and Goosse (1999) implentation i.e.
	832	!! transport proportional to the along-slope density gradient
	833	!!
	834	!! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591.
	835	!! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430.
	836	!!----------------------------------------------------------------------
	837	INTEGER , INTENT(in ) :: kjpt ! number of tracers
	838	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields
	839	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: ptb_ad ! before and now adjoint tracer fields
	840	REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta_ad ! tracer trend
	841	!
	842	INTEGER :: ji, jj, jk, jn ! dummy loop indices
	843	INTEGER :: iis , iid , ijs , ijd ! local integers
	844	INTEGER :: ikus, ikud, ikvs, ikvd ! - -
	845	REAL(wp) :: zbtr, ztra ! local scalars
	846	REAL(wp) :: ztraad ! - -
	847	REAL(wp) :: zu_bbl, zv_bbl ! - -
	848	REAL(wp) :: zu_bblad, zv_bblad ! - -
	849	!!----------------------------------------------------------------------
	850	!
	851	IF( nn_timing == 1 ) CALL timing_start( 'tra_bbl_adv_adj')
	852	!
	853	zu_bblad = 0.0_wp ; zv_bblad = 0.0_wp
	854	! ! ===========
	855	DO jn = 1, kjpt ! tracer loop
	856	! ! ===========
	857	# if defined key_vectopt_loop
	858	DO jj = 1, 1
	859	DO ji = jpij-jpi-1, 1, -1 ! vector opt. (forced unrolling)
	860	# else
	861	DO jj = jpjm1, 1, -1
	862	DO ji = jpim1, 1, -1 ! CAUTION start from i=1 to update i=2 when cyclic east-west
	863	# endif
	864	IF( vtr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero j-direction bbl advection
	865	! down-slope j/k-indices (deep) & up-slope j/k indices (shelf)
	866	ijd = jj + MAX( 0, mgrhv(ji,jj) ) ; ijs = jj + 1 - MAX( 0, mgrhv(ji,jj) )
	867	ikvd = mbkv_d(ji,jj) ; ikvs = mbkv(ji,jj)
	868	zv_bbl = ABS ( vtr_bbl(ji,jj) )
	869	! ! down-slope T-point (deep bottom point)
	870	zbtr = r1_e1e2t(ji,ijd) / fse3t(ji,ijd,ikvd)
	871	ztraad = pta_ad(ji,ijd,ikvd,jn)
	872	zv_bblad = zv_bblad + ztraad * ( ptb(ji,ijs,ikvs,jn) - ptb(ji,ijd,ikvd,jn) ) * zbtr
	873	ptb_ad(ji,ijs,ikvs,jn) = ptb_ad(ji,ijs,ikvs,jn) + ztraad * zv_bbl * zbtr
	874	ptb_ad(ji,ijd,ikvd,jn) = ptb_ad(ji,ijd,ikvd,jn) - ztraad * zv_bbl * zbtr
	875	!
	876	DO jk = ikvd-1, ikvs, -1 ! down-slope upper to down T-point (deep column)
	877	zbtr = r1_e1e2t(ji,ijd) / fse3t(ji,ijd,jk)
	878	ztraad = pta_ad(ji,ijd,jk,jn)
	879	zv_bblad = zv_bblad + ztraad * ( ptb(ji,ijd,jk+1,jn) - ptb(ji,ijd,jk,jn) ) * zbtr
	880	ptb_ad(ji,ijd,jk+1,jn) = ptb_ad(ji,ijd,jk+1,jn) + ztraad * zv_bbl * zbtr
	881	ptb_ad(ji,ijd,jk ,jn) = ptb_ad(ji,ijd,jk ,jn) - ztraad * zv_bbl * zbtr
	882	END DO
	883	! up -slope T-point (shelf bottom point)
	884	zbtr = r1_e1e2t(ji,ijs) / fse3t(ji,ijs,ikvs)
	885	ztraad = pta_ad(ji,ijs,ikvs,jn)
	886	zv_bblad = zv_bblad + ztraad * ( ptb(ji,ijd,ikvs,jn) - ptb(ji,ijs,ikvs,jn) ) * zbtr
	887	ptb_ad(ji,ijd,ikvs,jn) = ptb_ad(ji,ijd,ikvs,jn) + ztraad * zv_bbl * zbtr
	888	ptb_ad(ji,ijs,ikvs,jn) = ptb_ad(ji,ijs,ikvs,jn) - ztraad * zv_bbl * zbtr
	889
	890	!
	891	vtr_bbl_ad(ji,jj) = vtr_bbl_ad(ji,jj) + SIGN( zv_bblad, vtr_bbl(ji,jj) )
	892	zv_bblad = 0.0_wp
	893	!
	894	ENDIF
	895
	896
	897
	898	IF( utr_bbl(ji,jj) /= 0.e0 ) THEN ! non-zero i-direction bbl advection
	899	! down-slope i/k-indices (deep) & up-slope i/k indices (shelf)
	900	iid = ji + MAX( 0, mgrhu(ji,jj) ) ; iis = ji + 1 - MAX( 0, mgrhu(ji,jj) )
	901	ikud = mbku_d(ji,jj) ; ikus = mbku(ji,jj)
	902	zu_bbl = ABS( utr_bbl(ji,jj) )
	903	!
	904	zbtr = r1_e1e2t(iid,jj) / fse3t(iid,jj,ikud)
	905	ztraad = pta_ad(iid,jj,ikud,jn)
	906	zu_bblad = zu_bblad + ztraad * ( ptb(iis,jj,ikus,jn) - ptb(iid,jj,ikud,jn) ) * zbtr
	907	ptb_ad(iis,jj,ikus,jn) = ptb_ad(iis,jj,ikus,jn) + ztraad * zu_bbl * zbtr
	908	ptb_ad(iid,jj,ikud,jn) = ptb_ad(iid,jj,ikud,jn) + ztraad * zu_bbl * zbtr
	909	!
	910	DO jk = ikus, ikud-1 ! down-slope upper to down T-point (deep column)
	911	zbtr = r1_e1e2t(iid,jj) / fse3t(iid,jj,jk)
	912	ztraad = pta_ad(iid,jj,jk,jn)
	913	zu_bblad = zu_bblad + ztraad * ( ptb(iid,jj,jk+1,jn) - ptb(iid,jj,jk,jn) ) * zbtr
	914	ptb_ad(iid,jj,jk+1,jn) = ptb_ad(iid,jj,jk+1,jn) + ztraad * zu_bbl * zbtr
	915	ptb_ad(iid,jj,jk ,jn) = ptb_ad(iid,jj,jk ,jn) - ztraad * zu_bbl * zbtr
	916	END DO
	917	! ! up -slope T-point (shelf bottom point)
	918	zbtr = r1_e1e2t(iis,jj) / fse3t(iis,jj,ikus)
	919	ztraad = pta_ad(iis,jj,ikus,jn)
	920	zu_bblad = zu_bblad + ztraad * ( ptb(iid,jj,ikus,jn) - ptb(iis,jj,ikus,jn) ) * zbtr
	921	ptb_ad(iid,jj,ikus,jn) = ptb_ad(iid,jj,ikus,jn) + ztraad * zu_bbl * zbtr
	922	ptb_ad(iis,jj,ikus,jn) = ptb_ad(iis,jj,ikus,jn) - ztraad * zu_bbl * zbtr
	923	!
	924	utr_bbl_ad(ji,jj) = utr_bbl_ad(ji,jj) + SIGN( zu_bblad, utr_bbl(ji,jj) )
	925	zu_bblad = 0.0_wp
	926	!
	927	ENDIF
	928	!
	929	END DO
	930	!
	931	END DO
	932	! ! ===========
	933	END DO ! end tracer
	934	! ! ===========
	935	!
	936	IF( nn_timing == 1 ) CALL timing_stop( 'tra_bbl_adv_adj')
	937	!
	938	END SUBROUTINE tra_bbl_adv_adj
	939
	940
	941	SUBROUTINE bbl_adj( kt, kit000, cdtype )
	942	!!----------------------------------------------------------------------
	943	!! * ROUTINE bbl *
	944	!!
	945	!! ** Purpose : Computes the bottom boundary horizontal and vertical
	946	!! advection terms.
	947	!!
	948	!! ** Method :
	949	!! * diffusive bbl (nn_bbl_ldf=1) :
	950	!! When the product grad( rho) * grad(h) < 0 (where grad is an
	951	!! along bottom slope gradient) an additional lateral 2nd order
	952	!! diffusion along the bottom slope is added to the general
	953	!! tracer trend, otherwise the additional trend is set to 0.
	954	!! A typical value of ahbt is 2000 m2/s (equivalent to
	955	!! a downslope velocity of 20 cm/s if the condition for slope
	956	!! convection is satified)
	957	!! * advective bbl (nn_bbl_adv=1 or 2) :
	958	!! nn_bbl_adv = 1 use of the ocean velocity as bbl velocity
	959	!! nn_bbl_adv = 2 follow Campin and Goosse (1999) implentation
	960	!! i.e. transport proportional to the along-slope density gradient
	961	!!
	962	!! NB: the along slope density gradient is evaluated using the
	963	!! local density (i.e. referenced at a common local depth).
	964	!!
	965	!! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591.
	966	!! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430.
	967	!!----------------------------------------------------------------------
	968	!
	969	INTEGER , INTENT(in ) :: kt ! ocean time-step index
	970	INTEGER , INTENT(in ) :: kit000 ! first time step index
	971	CHARACTER(len=3), INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator)
	972	!!
	973	INTEGER :: ji, jj ! dummy loop indices
	974	INTEGER :: ik ! local integers
	975	INTEGER :: iis , iid , ijs , ijd ! - -
	976	INTEGER :: ikus, ikud, ikvs, ikvd ! - -
	977	REAL(wp) :: zsign, zsigna, zgbbl ! local scalars
	978	REAL(wp) :: zgdrho, zt, zs, zh ! - -
	979	REAL(wp) :: zgdrhoad, ztad, zsad, zhad! - -
	980	!!
	981	REAL(wp) :: fsalbt, fsbeta, pft, pfs, pfh ! statement function
	982	REAL(wp) :: pftad, pfsad, pfhad
	983	REAL(wp) :: fsalbt_adj_t, fsbeta_adj_t
	984	REAL(wp) :: fsalbt_adj_s, fsbeta_adj_s
	985	REAL(wp) :: fsalbt_adj_h, fsbeta_adj_h
	986	REAL(wp), POINTER, DIMENSION(:,:) :: zub, zvb, ztb, zsb, zdep
	987	REAL(wp), POINTER, DIMENSION(:,:) :: zubad, zvbad, ztbad, zsbad
	988	!!----------------------- zv_bbl-----------------------------------------------
	989	! ratio alpha/beta = fsalbt : ratio of thermal over saline expension coefficients
	990	! ================ pft : potential temperature in degrees celcius
	991	! pfs : salinity anomaly (s-35) in psu
	992	! pfh : depth in meters
	993	! nn_eos = 0 (Jackett and McDougall 1994 formulation)
	994	fsalbt( pft, pfs, pfh ) = & ! alpha/beta
	995	( ( ( -0.255019e-07 * pft + 0.298357e-05 ) * pft &
	996	- 0.203814e-03 ) * pft &
	997	+ 0.170907e-01 ) * pft &
	998	+ 0.665157e-01 &
	999	+(-0.678662e-05 * pfs - 0.846960e-04 * pft + 0.378110e-02 ) * pfs &
	1000	+ ( ( - 0.302285e-13 * pfh &
	1001	- 0.251520e-11 * pfs &
	1002	+ 0.512857e-12 * pft * pft ) * pfh &
	1003	- 0.164759e-06 * pfs &
	1004	+( 0.791325e-08 * pft - 0.933746e-06 ) * pft &
	1005	+ 0.380374e-04 ) * pfh
	1006	fsbeta( pft, pfs, pfh ) = & ! beta
	1007	( ( -0.415613e-09 * pft + 0.555579e-07 ) * pft &
	1008	- 0.301985e-05 ) * pft &
	1009	+ 0.785567e-03 &
	1010	+ ( 0.515032e-08 * pfs &
	1011	+ 0.788212e-08 * pft - 0.356603e-06 ) * pfs &
	1012	+( ( 0.121551e-17 * pfh &
	1013	- 0.602281e-15 * pfs &
	1014	- 0.175379e-14 * pft + 0.176621e-12 ) * pfh &
	1015	+ 0.408195e-10 * pfs &
	1016	+ ( - 0.213127e-11 * pft + 0.192867e-09 ) * pft &
	1017	- 0.121555e-07 ) * pfh
	1018
	1019	fsalbt_adj_t( pft, pfs, pfh, pftad ) = & ! alpha/beta
	1020	& ( - 0.255019e-07 * 4 * pft * pft * pft &
	1021	& + 0.298357e-05 * 3 * pft * pft &
	1022	& - 0.203814e-03 * 2 * pft &
	1023	& - 0.846960e-04 * pfs &
	1024	& + 0.512857e-12 * 2 * pft * pfh * pfh &
	1025	& + 0.791325e-08 * pft * pfh &
	1026	& - 0.933746e-06 * pfh &
	1027	& + 0.170907e-01 ) * pftad
	1028
	1029	fsalbt_adj_s( pft, pfs, pfh, pfsad ) = & ! alpha/beta
	1030	& + ( - 0.678662e-05 * 2 * pfs &
	1031	& - 0.846960e-04 * pft &
	1032	& - 0.251520e-11 * pfh * pfh &
	1033	& - 0.164759e-06 * pfh &
	1034	& + 0.378110e-02 ) * pfsad
	1035
	1036	fsalbt_adj_h( pft, pfs, pfh, pfhad ) = & ! alpha/beta
	1037	& + ( - 0.302285e-13 * 3 * pfh * pfh &
	1038	& - 0.251520e-11 * pfs * pfh &
	1039	& + 0.512857e-12 * pft * pft * pfh &
	1040	& - 0.164759e-06 * pfs &
	1041	& + 0.791325e-08 * pft * pft &
	1042	& - 0.933746e-06 * pft &
	1043	& + 0.380374e-04 ) * pfhad
	1044
	1045
	1046	fsbeta_adj_t( pft, pfs, pfh, pftad ) = & ! beta
	1047	& ( - 0.415613e-09 * 3 * pft * pft &
	1048	& + 0.555579e-07 * 2 * pft &
	1049	& - 0.301985e-05 &
	1050	& + 0.788212e-08 * pfs &
	1051	& - 0.213127e-11 * 2 * pfh * pft &
	1052	& - 0.175379e-14 * pfh * pfh ) * pftad
	1053	fsbeta_adj_s( pft, pfs, pfh, pfsad ) = & ! beta
	1054	& ( 0.788212e-08 * pft &
	1055	& + 0.515032e-08 * 2 * pfs &
	1056	& - 0.356603e-06 &
	1057	& + 0.408195e-10 * pfh &
	1058	& - 0.602281e-15 * pfh * pfh ) * pfsad
	1059	fsbeta_adj_h( pft, pfs, pfh, pfhad ) = & ! beta
	1060	& ( 0.121551e-17 * 3 * pfh * pfh &
	1061	& - 0.602281e-15 * 2 * pfs * pfh &
	1062	& - 0.175379e-14 * 2 * pft * pfh &
	1063	& + 0.176621e-12 * 2 * pfh &
	1064	& + 0.408195e-10 * pfs &
	1065	& + 0.192867e-09 * pfh &
	1066	& - 0.213127e-11 * pft * pft &
	1067	& + 0.192867e-09 * pft &
	1068	& - 0.121555e-07 ) * pfhad
	1069	!!----------------------------------------------------------------------
	1070
	1071	!
	1072	IF( nn_timing == 1 ) CALL timing_start( 'bbl_adj')
	1073	!
	1074	CALL wrk_alloc( jpi, jpj, zub , zvb , ztb , zsb, zdep, &
	1075	& zubad, zvbad, ztbad, zsbad )
	1076	!
	1077	zubad(:,:) = 0.0_wp ; zvbad(:,:) = 0.0_wp ; ztbad(:,:) = 0.0_wp ; zsbad(:,:) = 0.0_wp
	1078
	1079	IF( kt == kit000 ) THEN
	1080	IF(lwp) WRITE(numout,*)
	1081	IF(lwp) WRITE(numout,*) 'trabbl_tam:bbl_adj : Compute bbl velocities and diffusive coefficients in ', cdtype
	1082	IF(lwp) WRITE(numout,*) '~~~~~~~~~~'
	1083	ENDIF
	1084	! !* bottom temperature, salinity, velocity and depth
	1085	#if defined key_vectopt_loop
	1086	DO jj = 1, 1 ! vector opt. (forced unrolling)
	1087	DO ji = 1, jpij
	1088	#else
	1089	DO jj = 1, jpj
	1090	DO ji = 1, jpi
	1091	#endif
	1092	ik = mbkt(ji,jj) ! bottom T-level index
	1093	ztb (ji,jj) = tsb(ji,jj,ik,jp_tem) * tmask(ji,jj,1) ! bottom before T and S
	1094	zsb (ji,jj) = tsb(ji,jj,ik,jp_sal) * tmask(ji,jj,1)
	1095	zdep(ji,jj) = fsdept_0(ji,jj,ik) ! bottom T-level reference depth
	1096	!
	1097	zub(ji,jj) = un(ji,jj,mbku(ji,jj)) ! bottom velocity
	1098	zvb(ji,jj) = vn(ji,jj,mbkv(ji,jj))
	1099	END DO
	1100	END DO
	1101	! !-------------------!
	1102	IF( nn_bbl_adv /= 0 ) THEN ! advective bbl !
	1103	! !-------------------!
	1104	SELECT CASE ( nn_bbl_adv ) !* bbl transport type
	1105	!
	1106	CASE( 1 ) != use of upper velocity
	1107	! NOTE: not much needed for deriving, almost all the computations are for the SIGN, which is kept as in the NL
	1108	DO jj = 1, jpjm1 ! criteria: grad(rho).grad(h)<0 and grad(rho).grad(h)<0
	1109	DO ji = 1, fs_jpim1 ! vector opt.
	1110	! ! j-direction
	1111	zt = 0.5 * ( ztb (ji,jj+1) + ztb (ji,jj) ) ! T, S anomalie, and depth
	1112	zs = 0.5 * ( zsb (ji,jj+1) + zsb (ji,jj) ) - 35.0
	1113	zh = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) )
	1114	! ! masked bbl j-gradient of density
	1115	zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji,jj+1) - ztb(ji,jj) ) &
	1116	& - ( zsb(ji,jj+1) - zsb(ji,jj) ) ) * vmask(ji,jj,1)
	1117	zsign = SIGN( 0.5, - zgdrho * REAL( mgrhv(ji,jj) ) ) ! sign of j-gradient * j-slope
	1118	zsigna= SIGN( 0.5, zvb(ji,jj) * REAL( mgrhv(ji,jj) ) ) ! sign of u * i-slope
	1119	!
	1120	! ! bbl velocity
[3627]	1121	zvbad(ji,jj) = zvbad(ji,jj) + vtr_bbl_ad(ji,jj) * ( 0.5 + zsigna ) * ( 0.5 - zsign ) &
	1122	& * e1v(ji,jj) * e3v_bbl_0(ji,jj)
[3611]	1123	vtr_bbl_ad(ji,jj) = 0.0_wp
	1124	! ! i-direction
	1125	zt = 0.5 * ( ztb (ji,jj) + ztb (ji+1,jj) ) ! T, S anomalie, and depth
	1126	zs = 0.5 * ( zsb (ji,jj) + zsb (ji+1,jj) ) - 35.0
	1127	zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) )
	1128	! ! masked bbl i-gradient of density
	1129	zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji+1,jj) - ztb(ji,jj) ) &
	1130	& - ( zsb(ji+1,jj) - zsb(ji,jj) ) ) * umask(ji,jj,1)
	1131	!
	1132	zsign = SIGN( 0.5, - zgdrho * REAL( mgrhu(ji,jj) ) ) ! sign of i-gradient * i-slope
	1133	zsigna= SIGN( 0.5, zub(ji,jj) * REAL( mgrhu(ji,jj) ) ) ! sign of u * i-slope
	1134	!
	1135	! ! bbl velocity
[3627]	1136	zubad(ji,jj) = zubad(ji,jj) + utr_bbl_ad(ji,jj) * ( 0.5 + zsigna ) * ( 0.5 - zsign ) &
	1137	& * e2u(ji,jj) * e3u_bbl_0(ji,jj)
[3611]	1138	utr_bbl_ad(ji,jj) = 0.0_wp
	1139	!
	1140	END DO
	1141	END DO
	1142	!
	1143	CASE( 2 ) != bbl velocity = F( delta rho )
	1144	! NOTE: this one is nastier
	1145	zgbbl = grav * rn_gambbl
	1146	DO jj = jpjm1, 1, -1 ! criteria: rho_up > rho_down
	1147	DO ji = fs_jpim1, 1, -1 ! vector opt.
	1148	! ! j-direction
	1149	! down-slope T-point j/k-index (deep) & of the up -slope T-point j/k-index (shelf)
	1150	ijd = jj + MAX( 0, mgrhv(ji,jj) ) ; ijs = jj + 1 - MAX( 0, mgrhv(ji,jj) )
	1151	ikvd = mbkv_d(ji,jj) ; ikvs = mbkv(ji,jj)
	1152	!
	1153	! ! mid-depth density anomalie (up-slope minus down-slope)
	1154	zt = 0.5 * ( ztb (ji,jj) + ztb (ji,jj+1) ) ! mid slope depth of T, S, and depth
	1155	zs = 0.5 * ( zsb (ji,jj) + zsb (ji,jj+1) ) - 35.0
	1156	zh = 0.5 * ( zdep(ji,jj) + zdep(ji,jj+1) )
	1157	zgdrho = fsbeta( zt, zs, zh ) &
	1158	& * ( fsalbt( zt, zs, zh ) * ( ztb(ji,ijd) - ztb(ji,ijs) ) &
	1159	& - ( zsb(ji,ijd) - zsb(ji,ijs) ) ) * vmask(ji,jj,1)
	1160	!
	1161	zsign = SIGN( 0.5_wp, zgdrho ) ! adjoint of zgdrho = MAX( 0.e0, zgdrho )
	1162	! ! bbl transport (down-slope direction)
	1163	zgdrhoad = zsign * e1v(ji,jj) * e3v_bbl_0(ji,jj) * zgbbl * vtr_bbl_ad(ji,jj) * REAL( mgrhv(ji,jj) )
	1164	vtr_bbl_ad(ji,jj) = 0.0_wp
	1165
	1166	ztad = ( fsbeta_adj_t( zt, zs, zh, zgdrhoad ) &
	1167	& * ( fsalbt( zt, zs, zh ) * ( ztb(ji,ijd) - ztb(ji,ijs) ) &
	1168	& - ( zsb(ji,ijd) - zsb(ji,ijs) ) ) &
	1169	& + fsbeta( zt, zs, zh ) * fsalbt_adj_t( zt, zs, zh, zgdrhoad ) &
	1170	& * ( ztb(ji,ijd) - ztb(ji,ijs) ) &
	1171	& ) * vmask(ji,jj,1)
	1172	zsad = ( fsbeta_adj_s( zt, zs, zh, zgdrhoad ) &
	1173	& * ( fsalbt( zt, zs, zh ) * ( ztb(ji,ijd) - ztb(ji,ijs) ) &
	1174	& - ( zsb(ji,ijd) - zsb(ji,ijs) ) ) &
	1175	& + fsbeta( zt, zs, zh ) * fsalbt_adj_s( zt, zs, zh, zgdrhoad ) &
	1176	& * ( ztb(ji,ijd) - ztb(ji,ijs) ) &
	1177	& ) * vmask(ji,jj,1)
	1178	zhad = ( fsbeta_adj_h( zt, zs, zh, zgdrhoad ) &
	1179	& * ( fsalbt( zt, zs, zh ) * ( ztb(ji,ijd) - ztb(ji,ijs) ) &
	1180	& - ( zsb(ji,ijd) - zsb(ji,ijs) ) ) &
	1181	& + fsbeta( zt, zs, zh ) * fsalbt_adj_h( zt, zs, zh, zgdrhoad ) &
	1182	& * ( ztb(ji,ijd) - ztb(ji,ijs) ) &
	1183	& ) * vmask(ji,jj,1)
	1184
	1185	ztbad(ji,ijd) = ztbad(ji,ijd) + zgdrhoad * fsbeta( zt, zs, zh ) * fsalbt( zt, zs, zh ) * vmask(ji,jj,1)
	1186	ztbad(ji,ijs) = ztbad(ji,ijs) - zgdrhoad * fsbeta( zt, zs, zh ) * fsalbt( zt, zs, zh ) * vmask(ji,jj,1)
	1187	zsbad(ji,ijd) = zsbad(ji,ijd) - zgdrhoad * fsbeta( zt, zs, zh ) * vmask(ji,jj,1)
	1188	zsbad(ji,ijs) = zsbad(ji,ijs) + zgdrhoad * fsbeta( zt, zs, zh ) * vmask(ji,jj,1)
	1189
	1190	ztbad (ji,jj ) = ztbad (ji,jj ) + 0.5 * ztad
	1191	ztbad (ji,jj+1) = ztbad (ji,jj+1) + 0.5 * ztad
	1192	zsbad (ji,jj ) = zsbad (ji,jj ) + 0.5 * zsad
	1193	zsbad (ji,jj+1) = zsbad (ji,jj+1) + 0.5 * zsad
	1194	ztad = 0.0_wp ; zsad = 0.0_wp ; zhad = 0.0_wp
	1195
	1196	! ! i-direction
	1197	! down-slope T-point i/k-index (deep) & up-slope T-point i/k-index (shelf)
	1198	iid = ji + MAX( 0, mgrhu(ji,jj) ) ; iis = ji + 1 - MAX( 0, mgrhu(ji,jj) )
	1199	ikud = mbku_d(ji,jj) ; ikus = mbku(ji,jj)
	1200	!
	1201	! ! mid-depth density anomalie (up-slope minus down-slope)
	1202	zt = 0.5 * ( ztb (ji,jj) + ztb (ji+1,jj) ) ! mid slope depth of T, S, and depth
	1203	zs = 0.5 * ( zsb (ji,jj) + zsb (ji+1,jj) ) - 35.0
	1204	zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) )
	1205	zgdrho = fsbeta( zt, zs, zh ) &
	1206	& * ( fsalbt( zt, zs, zh ) * ( ztb(iid,jj) - ztb(iis,jj) ) &
	1207	& - ( zsb(iid,jj) - zsb(iis,jj) ) ) * umask(ji,jj,1)
	1208	zsign = SIGN( 0.5_wp, zgdrho ) ! adjoint of zgdrho = MAX( 0.e0, zgdrho )
	1209	! ! bbl transport (down-slope direction)
	1210	zgdrhoad = zsign * e2u(ji,jj) * e3u_bbl_0(ji,jj) * zgbbl * utr_bbl_ad(ji,jj) * REAL( mgrhu(ji,jj) )
	1211	utr_bbl_ad(ji,jj) = 0.0_wp
	1212	!
	1213	ztad = ( fsbeta_adj_t( zt, zs, zh, zgdrhoad ) &
	1214	& * ( fsalbt( zt, zs, zh ) * ( ztb(iid,jj) - ztb(iis, jj) ) &
	1215	& - ( zsb(iid,jj) - zsb(iis, jj) ) ) &
	1216	& + fsbeta( zt, zs, zh ) * fsalbt_adj_t( zt, zs, zh, zgdrhoad ) &
	1217	& * ( ztb(iid,jj) - ztb(iis, jj) ) &
	1218	& ) * umask(ji,jj,1)
	1219	zsad = ( fsbeta_adj_s( zt, zs, zh, zgdrhoad ) &
	1220	& * ( fsalbt( zt, zs, zh ) * ( ztb(iid,jj) - ztb(iis, jj) ) &
	1221	& - ( zsb(iid,jj) - zsb(iis, jj) ) ) &
	1222	& + fsbeta( zt, zs, zh ) * fsalbt_adj_s( zt, zs, zh, zgdrhoad ) &
	1223	& * ( ztb(iid,jj) - ztb(iis, jj) ) &
	1224	& ) * umask(ji,jj,1)
	1225	zhad = ( fsbeta_adj_h( zt, zs, zh, zgdrhoad ) &
	1226	& * ( fsalbt( zt, zs, zh ) * ( ztb(iid,jj) - ztb(iis, jj) ) &
	1227	& - ( zsb(iid,jj) - zsb(iis, jj) ) ) &
	1228	& + fsbeta( zt, zs, zh ) * fsalbt_adj_h( zt, zs, zh, zgdrhoad ) &
	1229	& * ( ztb(iid,jj) - ztb(iis, jj) ) &
	1230	& ) * umask(ji,jj,1)
	1231
	1232	ztbad(iid,jj) = ztbad(iid,jj) + zgdrhoad * fsbeta( zt, zs, zh ) * fsalbt( zt, zs, zh ) * umask(ji,jj,1)
	1233	ztbad(iis,jj) = ztbad(iis,jj) - zgdrhoad * fsbeta( zt, zs, zh ) * fsalbt( zt, zs, zh ) * umask(ji,jj,1)
	1234	zsbad(iid,jj) = zsbad(iid,jj) - zgdrhoad * fsbeta( zt, zs, zh ) * umask(ji,jj,1)
	1235	zsbad(iis,jj) = zsbad(iis,jj) + zgdrhoad * fsbeta( zt, zs, zh ) * umask(ji,jj,1)
	1236	zgdrhoad = 0.0_wp
	1237
	1238	ztbad (ji,jj ) = ztbad (ji,jj ) + 0.5 * ztad
	1239	ztbad (ji+1,jj) = ztbad (ji+1,jj) + 0.5 * ztad
	1240	zsbad (ji,jj ) = zsbad (ji,jj ) + 0.5 * zsad
	1241	zsbad (ji+1,jj) = zsbad (ji+1,jj) + 0.5 * zsad
	1242	ztad = 0.0_wp ; zsad = 0.0_wp ; zhad = 0.0_wp
	1243	!
	1244	END DO
	1245	END DO
	1246	END SELECT
	1247	!
	1248	ENDIF
	1249	IF( nn_bbl_ldf == 1 ) THEN ! diffusive bbl !
	1250	! !-------------------!
	1251	! NOTE : while rn_ahtbbl remains a passive variable, the code below will only yield ah_bbl_ad=0
	1252	#if defined key_control_param
	1253	DO jj = 1, jpjm1 ! (criteria for non zero flux: grad(rho).grad(h) < 0 )
	1254	DO ji = 1, jpim1
	1255	! ! j-direction
	1256	zt = 0.5 * ( ztb (ji,jj+1) + ztb (ji,jj) ) ! T, S anomalie, and depth
	1257	zs = 0.5 * ( zsb (ji,jj+1) + zsb (ji,jj) ) - 35.0
	1258	zh = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) )
	1259	! ! masked bbl j-gradient of density
	1260	zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji,jj+1) - ztb(ji,jj) ) &
	1261	& - ( zsb(ji,jj+1) - zsb(ji,jj) ) ) * vmask(ji,jj,1)
	1262	!
	1263	zsign = SIGN( 0.5, -zgdrho * REAL( mgrhv(ji,jj) ) ) ! sign of ( j-gradient * j-slope )
	1264	ahv_bbl_0_ad(ji,jj) = ( 0.5 - zsign ) * ahv_bbl_ad(ji,jj)
	1265	ahv_bbl_ad(ji,jj) = 0.0_wp
	1266	!
	1267	! ! i-direction
	1268	zt = 0.5 * ( ztb (ji,jj) + ztb (ji+1,jj) ) ! T, S anomalie, and depth
	1269	zs = 0.5 * ( zsb (ji,jj) + zsb (ji+1,jj) ) - 35.0
	1270	zh = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) )
	1271	! ! masked bbl i-gradient of density
	1272	zgdrho = ( fsalbt( zt, zs, zh ) * ( ztb(ji+1,jj) - ztb(ji,jj) ) &
	1273	& - ( zsb(ji+1,jj) - zsb(ji,jj) ) ) * umask(ji,jj,1)
	1274	!
	1275	zsign = SIGN( 0.5, - zgdrho * REAL( mgrhu(ji,jj) ) ) ! sign of ( i-gradient * i-slope )
	1276	ahu_bbl_0_ad(ji,jj) = ( 0.5 - zsign ) * ahu_bbl_ad(ji,jj) ! masked diffusive flux coeff.
	1277	ahu_bbl_ad(ji,jj) = 0.0_wp
	1278	!
	1279	END DO
	1280	END DO
	1281	#else
	1282	DO jj = 1, jpjm1
	1283	DO ji = 1, jpim1
	1284	ahu_bbl_ad(ji,jj)=0.0_wp
	1285	ahv_bbl_ad(ji,jj)=0.0_wp
	1286	END DO
	1287	END DO
	1288	#endif
	1289	!
	1290	ENDIF
	1291	! !* bottom temperature, salinity, velocity and depth
	1292	#if defined key_vectopt_loop
	1293	DO jj = 1, 1 ! vector opt. (forced unrolling)
	1294	DO ji = 1, jpij
	1295	#else
	1296	DO jj = 1, jpj
	1297	DO ji = 1, jpi
	1298	#endif
	1299	ik = mbkt(ji,jj) ! bottom T-level index
	1300	tsb_ad(ji,jj,ik,jp_tem) = tsb_ad(ji,jj,ik,jp_tem) + ztbad(ji,jj) * tmask(ji,jj,1)
	1301	tsb_ad(ji,jj,ik,jp_sal) = tsb_ad(ji,jj,ik,jp_sal) + zsbad(ji,jj) * tmask(ji,jj,1)
	1302	ztbad (ji,jj) = 0.0_wp
	1303	zsbad (ji,jj) = 0.0_wp
	1304	END DO
	1305	END DO
	1306	! !-------------------!
	1307	!
	1308	CALL wrk_dealloc( jpi, jpj, zub , zvb , ztb , zsb, zdep, &
	1309	& ztbad, zsbad, ztbad, zsbad )
	1310	!
	1311	IF( nn_timing == 1 ) CALL timing_stop( 'bbl_adj')
	1312	!
	1313	END SUBROUTINE bbl_adj
	1314
	1315
	1316	SUBROUTINE tra_bbl_init_tam
	1317	!!----------------------------------------------------------------------
	1318	!! * ROUTINE tra_bbl_init *
	1319	!!
	1320	!! ** Purpose : Initialization for the bottom boundary layer scheme.
	1321	!!
	1322	!! ** Method :
	1323	!!----------------------------------------------------------------------
	1324	!
	1325	integer :: ierr
	1326	IF( nn_timing == 1 ) CALL timing_start( 'tra_bbl_init_tam')
	1327	!
	1328	ierr = tra_bbl_alloc_tam( 0 )
	1329
	1330	ahu_bbl_0_tl = 0.0_wp
	1331	ahv_bbl_0_tl = 0.0_wp
	1332	ahu_bbl_0_ad = 0.0_wp
	1333	ahv_bbl_0_ad = 0.0_wp
	1334	!
	1335	IF( nn_timing == 1 ) CALL timing_stop( 'tra_bbl_init_tam')
	1336	!
	1337	END SUBROUTINE tra_bbl_init_tam
	1338
	1339	SUBROUTINE tra_bbl_adj_tst( kumadt )
	1340	!!-----------------------------------------------------------------------
	1341	!!
	1342	!! * ROUTINE tra_bbl_adj_tst *
	1343	!!
	1344	!! ** Purpose : Test the adjoint routine.
	1345	!!
	1346	!! ** Method : Verify the scalar product
	1347	!!
	1348	!! ( L dx )^T W dy = dx^T L^T W dy
	1349	!!
	1350	!! where L = tangent routine
	1351	!! L^T = adjoint routine
	1352	!! W = diagonal matrix of scale factors
	1353	!! dx = input perturbation (random field)
	1354	!! dy = L dx
	1355	!!
	1356	!!
	1357	!! History :
	1358	!! ! 08-08 (A. Vidard)
	1359	!!-----------------------------------------------------------------------
	1360	!! * Modules used
	1361
	1362	!! * Arguments
	1363	INTEGER, INTENT(IN) :: &
	1364	& kumadt ! Output unit
	1365
	1366	!! * Local declarations
	1367	INTEGER :: &
	1368	& ji, & ! dummy loop indices
	1369	& jj, &
	1370	& jk, &
	1371	& jtst
	1372	INTEGER :: &
	1373	& jsav1, &
	1374	& jsav2
	1375	REAL(KIND=wp) :: &
	1376	& zsp1, & ! scalar product involving the tangent routine
	1377	& zsp2 ! scalar product involving the adjoint routine
	1378	REAL(KIND=wp), POINTER, DIMENSION(:,:,:,:) :: &
	1379	& ztsa_tlin , & ! Tangent input
	1380	& ztsa_tlout, &
	1381	& ztsb_tlin , &
	1382	& ztsa_adout, & ! Adjoint input
	1383	& ztsa_adin , &
	1384	& ztsb_adout, &
	1385	& zrts ! 2*3D random field
	1386	REAL(KIND=wp), POINTER, DIMENSION(:,:) :: &
	1387	& zutr_tlin , &
	1388	& zutr_tlout, &
	1389	& zvtr_tlin , &
	1390	& zvtr_tlout, &
	1391	& zutr_adout, &
	1392	& zutr_adin , &
	1393	& zvtr_adout, &
	1394	& zvtr_adin , &
	1395	& zr2 ! 2D random field
	1396	CHARACTER(LEN=14) :: &
	1397	& cl_name
	1398	! Allocate memory
	1399
	1400	CALL wrk_alloc( jpi, jpj, jpk, jpts, ztsa_tlin , ztsa_tlout, ztsb_tlin , &
	1401	& ztsa_adout, ztsa_adin , ztsb_adout, &
	1402	& zrts )
	1403	CALL wrk_alloc( jpi, jpj, zutr_tlin , zutr_tlout, zvtr_tlin , zvtr_tlout, &
	1404	& zutr_adout, zutr_adin , zvtr_adout, zvtr_adin , &
	1405	& zr2 )
	1406
	1407	CALL grid_random( utr_bbl(:,:), 'U', 0.0_wp, stdu )
	1408	CALL grid_random( vtr_bbl(:,:), 'V', 0.0_wp, stdv )
	1409
	1410	jsav1 = nn_bbl_ldf
	1411	jsav2 = nn_bbl_adv
	1412
	1413	DO jtst = 1, 2
	1414	!==================================================================
	1415	! 1) dx = ( un_tl, vn_tl, hdivn_tl ) and
	1416	! dy = ( hdivb_tl, hdivn_tl )
	1417	!==================================================================
	1418
	1419	SELECT CASE( jtst)
	1420	CASE ( 1 )
	1421	nn_bbl_ldf = 1
	1422	nn_bbl_adv = 0
	1423	CASE ( 2 )
	1424	nn_bbl_ldf = 0
	1425	nn_bbl_adv = 1
	1426	CASE ( 3 )
	1427	nn_bbl_ldf = 0
	1428	nn_bbl_adv = 2
	1429	END SELECT
	1430	!--------------------------------------------------------------------
	1431	! Reset the tangent and adjoint variables
	1432	!--------------------------------------------------------------------
	1433	ztsa_tlin (:,:,:,:) = 0.0_wp
	1434	ztsa_tlout(:,:,:,:) = 0.0_wp
	1435	ztsb_tlin (:,:,:,:) = 0.0_wp
	1436	ztsa_adout(:,:,:,:) = 0.0_wp
	1437	ztsa_adin (:,:,:,:) = 0.0_wp
	1438	ztsb_adout(:,:,:,:) = 0.0_wp
	1439
	1440	zutr_tlin (:,:) = 0.0_wp
	1441	zutr_tlout(:,:) = 0.0_wp
	1442	zvtr_tlin (:,:) = 0.0_wp
	1443	zvtr_tlout(:,:) = 0.0_wp
	1444	zutr_adout(:,:) = 0.0_wp
	1445	zutr_adin (:,:) = 0.0_wp
	1446	zvtr_adout(:,:) = 0.0_wp
	1447	zvtr_adin (:,:) = 0.0_wp
	1448
	1449	tsb_ad(:,:,:,:) = 0.0_wp
	1450	!--------------------------------------------------------------------
	1451	! Initialize the tangent input with random noise: dx
	1452	!--------------------------------------------------------------------
	1453
	1454	CALL grid_random( zrts(:,:,:,jp_tem), 'T', 0.0_wp, stdt )
	1455	CALL grid_random( zrts(:,:,:,jp_sal), 'T', 0.0_wp, stds )
	1456	DO jk = 1, jpk
	1457	DO jj = nldj, nlej
	1458	DO ji = nldi, nlei
	1459	ztsa_tlin(ji,jj,jk,:) = zrts(ji,jj,jk,:)
	1460	END DO
	1461	END DO
	1462	END DO
	1463
	1464	CALL grid_random( zrts(:,:,:,jp_tem), 'T', 0.0_wp, stdt )
	1465	CALL grid_random( zrts(:,:,:,jp_sal), 'T', 0.0_wp, stds )
	1466	DO jk = 1, jpk
	1467	DO jj = nldj, nlej
	1468	DO ji = nldi, nlei
	1469	ztsb_tlin(ji,jj,jk,:) = zrts(ji,jj,jk,:)
	1470	END DO
	1471	END DO
	1472	END DO
	1473
	1474	CALL grid_random( zr2(:,:), 'U', 0.0_wp, stdu )
	1475	DO jj = nldj, nlej
	1476	DO ji = nldi, nlei
	1477	zutr_tlin(ji,jj) = zr2(ji,jj)
	1478	END DO
	1479	END DO
	1480
	1481	CALL grid_random( zr2(:,:), 'V', 0.0_wp, stdv )
	1482	DO jj = nldj, nlej
	1483	DO ji = nldi, nlei
	1484	zvtr_tlin(ji,jj) = zr2(ji,jj)
	1485	END DO
	1486	END DO
	1487
	1488	tsa_tl(:,:,:,:) = ztsa_tlin(:,:,:,:)
	1489	tsb_tl(:,:,:,:) = ztsb_tlin(:,:,:,:)
	1490	utr_bbl_tl(:,:) = zutr_tlin(:,:)
	1491	vtr_bbl_tl(:,:) = zvtr_tlin(:,:)
	1492	CALL tra_bbl_tan ( nit000 )
	1493	ztsa_tlout(:,:,:,:) = tsa_tl(:,:,:,:)
	1494	zutr_tlout(:,:) = utr_bbl_tl(:,:)
	1495	zvtr_tlout(:,:) = vtr_bbl_tl(:,:)
	1496	!--------------------------------------------------------------------
	1497	! Initialize the adjoint variables: dy^* = W dy
	1498	!--------------------------------------------------------------------
	1499
	1500	DO jk = 1, jpk
	1501	DO jj = nldj, nlej
	1502	DO ji = nldi, nlei
	1503	ztsa_adin(ji,jj,jk,jp_tem) = ztsa_tlout(ji,jj,jk,jp_tem) &
	1504	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
	1505	& * tmask(ji,jj,jk)
	1506	ztsa_adin(ji,jj,jk,jp_sal) = ztsa_tlout(ji,jj,jk,jp_sal) &
	1507	& * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) &
	1508	& * tmask(ji,jj,jk)
	1509	END DO
	1510	END DO
	1511	END DO
	1512	DO jj = nldj, nlej
	1513	DO ji = nldi, nlei
	1514	zutr_adin(ji,jj) = zutr_tlout(ji,jj) &
	1515	& * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,1) &
	1516	& * vmask(ji,jj,jk)
	1517	zvtr_adin(ji,jj) = zvtr_tlout(ji,jj) &
	1518	& * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,1) &
	1519	& * vmask(ji,jj,jk)
	1520	END DO
	1521	END DO
	1522	!--------------------------------------------------------------------
	1523	! Compute the scalar product: ( L dx )^T W dy
	1524	!--------------------------------------------------------------------
	1525
	1526	zsp1 = DOT_PRODUCT( ztsa_tlout(:,:,:,jp_tem), ztsa_adin(:,:,:,jp_tem) ) &
	1527	& + DOT_PRODUCT( ztsa_tlout(:,:,:,jp_sal), ztsa_adin(:,:,:,jp_sal) ) &
	1528	& + DOT_PRODUCT( zutr_tlout, zutr_adin ) + DOT_PRODUCT( zvtr_tlout, zvtr_adin )
	1529
	1530	!--------------------------------------------------------------------
	1531	! Call the adjoint routine: dx^* = L^T dy^*
	1532	!--------------------------------------------------------------------
	1533
	1534	tsa_ad(:,:,:,:) = ztsa_adin(:,:,:,:)
	1535	utr_bbl_ad(:,:) = zutr_adin(:,:)
	1536	vtr_bbl_ad(:,:) = zvtr_adin(:,:)
	1537	CALL tra_bbl_adj ( nit000 )
	1538	ztsa_adout(:,:,:,:) = tsa_ad(:,:,:,:)
	1539	ztsb_adout(:,:,:,:) = tsb_ad(:,:,:,:)
	1540	zutr_adout(:,:) = utr_bbl_ad(:,:)
	1541	zvtr_adout(:,:) = vtr_bbl_ad(:,:)
	1542
	1543
	1544	zsp2 = DOT_PRODUCT( ztsa_tlin(:,:,:,jp_tem), ztsa_adout(:,:,:,jp_tem) ) &
	1545	& + DOT_PRODUCT( ztsa_tlin(:,:,:,jp_sal), ztsa_adout(:,:,:,jp_sal) ) &
	1546	& + DOT_PRODUCT( ztsb_tlin(:,:,:,jp_tem), ztsb_adout(:,:,:,jp_tem) ) &
	1547	& + DOT_PRODUCT( ztsb_tlin(:,:,:,jp_sal), ztsb_adout(:,:,:,jp_sal) ) &
	1548	& + DOT_PRODUCT( zutr_tlin, zutr_adout ) + DOT_PRODUCT( zvtr_tlin, zvtr_adout )
	1549
	1550	SELECT CASE ( jtst )
	1551	CASE ( 1 )
	1552	! 14 char:'12345678901234'
	1553	cl_name = 'trabbl_adj_dif'
	1554	CASE ( 2 )
	1555	! 14 char:'12345678901234'
	1556	cl_name = 'trabbl_ad_adv1'
	1557	CASE ( 3 )
	1558	! 14 char:'12345678901234'
	1559	cl_name = 'trabbl_ad_adv2'
	1560	END SELECT
	1561	CALL prntst_adj( cl_name, kumadt, zsp1, zsp2 )
	1562
	1563	END DO
	1564
	1565	CALL wrk_dealloc( jpi, jpj, jpk, jpts, ztsa_tlin , ztsa_tlout, ztsb_tlin , &
	1566	& ztsa_adout, ztsa_adin , ztsb_adout, &
	1567	& zrts )
	1568	CALL wrk_dealloc( jpi, jpj, zutr_tlin , zutr_tlout, zvtr_tlin , zvtr_tlout, &
	1569	& zutr_adout, zutr_adin , zvtr_adout, zvtr_adin , &
	1570	& zr2 )
	1571
	1572	nn_bbl_ldf = jsav1
	1573	nn_bbl_adv = jsav2
	1574
	1575
	1576	END SUBROUTINE tra_bbl_adj_tst
	1577
	1578	#else
	1579	!!----------------------------------------------------------------------
	1580	!! Dummy module : No bottom boundary layer scheme
	1581	!!----------------------------------------------------------------------
	1582	CONTAINS
	1583	SUBROUTINE tra_bbl_init_tam ! Dummy routine
	1584	END SUBROUTINE tra_bbl_init_tam
	1585	SUBROUTINE tra_bbl_tan( kt ) ! Dummy routine
	1586	WRITE(,) 'tra_bbl_tan: You should not have seen this print! error?', kt
	1587	END SUBROUTINE tra_bbl_tan
	1588	SUBROUTINE tra_bbl_adj( kt ) ! Dummy routine
	1589	WRITE(,) 'tra_bbl_adj: You should not have seen this print! error?', kt
	1590	END SUBROUTINE tra_bbl_adj
	1591	SUBROUTINE tra_bbl_adj_tst( kt ) ! Dummy routine
	1592	WRITE(,) 'tra_bbl_adj_tst: You should not have seen this print! error?', kt
	1593	END SUBROUTINE tra_bbl_adj_tst
	1594	#endif
	1595	!!======================================================================
	1596	END MODULE trabbl_tam

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

Download in other formats: