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

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

trabbl.F90 in NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA – NEMO

Context Navigation

source: NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trabbl.F90 @ 10946

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

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

Download in other formats: