Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

trabbl.F90 @ 12286

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

Download in other formats: