Context Navigation

cla.F90 @ 2392

Last change on this file since 2392 was 2392, checked in by gm, 13 years ago
v3.3beta: Cross Land Advection (ticket #127) full rewriting + MPP bug corrections
Property svn:keywords set to `Id`
File size: 41.6 KB

Rev	Line
[3]	1	MODULE cla
[2392]	2	!!======================================================================
	3	!! * MODULE cla *
	4	!! Cross Land Advection : specific update of the horizontal divergence,
	5	!! tracer trends and after velocity
	6	!!
	7	!! --- Specific to ORCA_R2 ---
	8	!!
	9	!!======================================================================
	10	!! History : 1.0 ! 2002-11 (A. Bozec) Original code
	11	!! 3.2 ! 2009-07 (G. Madec) merge cla, cla_div, tra_cla, cla_dynspg
	12	!! ! and correct a mpp bug reported by A.R. Porter
	13	!!----------------------------------------------------------------------
[3]	14	#if defined key_orca_r2
	15	!!----------------------------------------------------------------------
[2392]	16	!! 'key_orca_r2' global ocean model R2
[3]	17	!!----------------------------------------------------------------------
[2392]	18	!! cla_div : update of horizontal divergence at cla straits
	19	!! tra_cla : update of tracers at cla straits
	20	!! cla_dynspg : update of after horizontal velocities at cla straits
	21	!! cla_init : initialisation - control check
	22	!! cla_bab_el_mandeb : cross land advection for Bab-el-mandeb strait
	23	!! cla_gibraltar : cross land advection for Gibraltar strait
	24	!! cla_hormuz : cross land advection for Hormuz strait
[3]	25	!!----------------------------------------------------------------------
[2392]	26	USE oce ! ocean dynamics and tracers
	27	USE dom_oce ! ocean space and time domain
	28	USE sbc_oce ! surface boundary condition: ocean
	29	USE dynspg_oce ! ocean dynamics: surface pressure gradient variables
	30	USE in_out_manager ! I/O manager
	31	USE lib_mpp ! distributed memory computing library
	32	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
[3]	33
	34	IMPLICIT NONE
	35	PRIVATE
[2392]	36
	37	PUBLIC cla_init ! routine called by opa.F90
	38	PUBLIC cla_div ! routine called by divcur.F90
	39	PUBLIC cla_traadv ! routine called by traadv.F90
	40	PUBLIC cla_dynspg ! routine called by dynspg_flt.F90
[3]	41
[2392]	42	INTEGER :: nbab, ngib, nhor ! presence or not of required grid-point on local domain
	43	! ! for Bab-el-Mandeb, Gibraltar, and Hormuz straits
	44
	45	! !!! profile of hdiv for some straits
	46	REAL(wp), DIMENSION (jpk) :: hdiv_139_101, hdiv_139_101_kt ! Gibraltar strait, fixed & time evolving part (i,j)=(172,101)
	47	REAL(wp), DIMENSION (jpk) :: hdiv_139_102 ! Gibraltar strait, fixed part only (i,j)=(139,102)
	48	REAL(wp), DIMENSION (jpk) :: hdiv_141_102, hdiv_141_102_kt ! Gibraltar strait, fixed & time evolving part (i,j)=(141,102)
	49	REAL(wp), DIMENSION (jpk) :: hdiv_161_88 , hdiv_161_88_kt ! Bab-el-Mandeb strait, fixed & time evolving part (i,j)=(161,88)
	50	REAL(wp), DIMENSION (jpk) :: hdiv_161_87 ! Bab-el-Mandeb strait, fixed part only (i,j)=(161,87)
	51	REAL(wp), DIMENSION (jpk) :: hdiv_160_89 , hdiv_160_89_kt ! Bab-el-Mandeb strait, fixed & time evolving part (i,j)=(160,89)
	52	REAL(wp), DIMENSION (jpk) :: hdiv_172_94 ! Hormuz strait, fixed part only (i,j)=(172, 94)
[3]	53
[2392]	54	REAL(wp), DIMENSION (jpk) :: t_171_94_hor, s_171_94_hor ! Temperature, salinity in the Hormuz strait
[3]	55
	56	!! * Substitutions
	57	# include "domzgr_substitute.h90"
	58	!!----------------------------------------------------------------------
[2287]	59	!! NEMO/OPA 3.3 , NEMO Consortium (2010)
[888]	60	!! $Id$
[2392]	61	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[3]	62	!!----------------------------------------------------------------------
	63	CONTAINS
	64
[2392]	65	SUBROUTINE cla_div( kt )
[3]	66	!!----------------------------------------------------------------------
[2392]	67	!! * ROUTINE div_cla *
	68	!!
	69	!! ** Purpose : update the horizontal divergence of the velocity field
	70	!! at some straits ( Gibraltar, Bab el Mandeb and Hormuz ).
	71	!!
	72	!! ** Method : - first time-step: initialisation of cla
	73	!! - all time-step: using imposed transport at each strait,
	74	!! the now horizontal divergence is updated
	75	!!
	76	!! ** Action : phdivn updted now horizontal divergence at cla straits
	77	!!----------------------------------------------------------------------
	78	INTEGER, INTENT( in ) :: kt ! ocean time-step index
	79	!!----------------------------------------------------------------------
	80	!
	81	IF( kt == nit000 ) THEN
	82	!
	83	CALL cla_init ! control check
	84	!
	85	IF(lwp) WRITE(numout,*)
	86	IF(lwp) WRITE(numout,*) 'div_cla : cross land advection on hdiv '
	87	IF(lwp) WRITE(numout,*) '~~~~~~~~'
	88	!
	89	IF( nbab == 1 ) CALL cla_bab_el_mandeb('ini') ! Bab el Mandeb ( Red Sea - Indian ocean )
	90	IF( ngib == 1 ) CALL cla_gibraltar ('ini') ! Gibraltar strait (Med Sea - Atlantic ocean)
	91	IF( nhor == 1 ) CALL cla_hormuz ('ini') ! Hormuz Strait ( Persian Gulf - Indian ocean )
	92	!
	93	ENDIF
	94	!
	95	IF( nbab == 1 ) CALL cla_bab_el_mandeb('div') ! Bab el Mandeb ( Red Sea - Indian ocean )
	96	IF( ngib == 1 ) CALL cla_gibraltar ('div') ! Gibraltar strait (Med Sea - Atlantic ocean)
	97	IF( nhor == 1 ) CALL cla_hormuz ('div') ! Hormuz Strait ( Persian Gulf - Indian ocean )
	98	!
	99	!!gm lbc useless here, no?
	100	!!gm CALL lbc_lnk( hdivn, 'T', 1. ) ! Lateral boundary conditions on hdivn
	101	!
	102	END SUBROUTINE cla_div
	103
	104
	105	SUBROUTINE cla_traadv( kt )
	106	!!----------------------------------------------------------------------
[3]	107	!! * ROUTINE tra_cla *
	108	!!
	109	!! ** Purpose : Update the now trend due to the advection of tracers
	110	!! and add it to the general trend of passive tracer equations
	111	!! at some straits ( Bab el Mandeb, Gibraltar, Hormuz ).
	112	!!
[2392]	113	!! ** Method : using both imposed transport at each strait and T & S
	114	!! budget, the now tracer trends is updated
[3]	115	!!
[2392]	116	!! ** Action : (ta,sa) updated now tracer trends at cla straits
	117	!!----------------------------------------------------------------------
	118	INTEGER, INTENT( in ) :: kt ! ocean time-step index
	119	!!----------------------------------------------------------------------
	120	!
	121	IF( kt == nit000 ) THEN
	122	IF(lwp) WRITE(numout,*)
	123	IF(lwp) WRITE(numout,*) 'tra_cla : cross land advection on tracers '
	124	IF(lwp) WRITE(numout,*) '~~~~~~~~'
	125	ENDIF
	126	!
	127	IF( nbab == 1 ) CALL cla_bab_el_mandeb('tra') ! Bab el Mandeb strait
	128	IF( ngib == 1 ) CALL cla_gibraltar ('tra') ! Gibraltar strait
	129	IF( nhor == 1 ) CALL cla_hormuz ('tra') ! Hormuz Strait ( Persian Gulf)
	130	!
	131	END SUBROUTINE cla_traadv
	132
	133
	134	SUBROUTINE cla_dynspg( kt )
	135	!!----------------------------------------------------------------------
	136	!! * ROUTINE cla_dynspg *
	137	!!
	138	!! ** Purpose : Update the after velocity at some straits
	139	!! (Bab el Mandeb, Gibraltar, Hormuz).
[3]	140	!!
[2392]	141	!! ** Method : required to compute the filtered surface pressure gradient
	142	!!
	143	!! ** Action : (ua,va) after velocity at the cla straits
[3]	144	!!----------------------------------------------------------------------
	145	INTEGER, INTENT( in ) :: kt ! ocean time-step index
	146	!!----------------------------------------------------------------------
[2392]	147	!
	148	IF( kt == nit000 ) THEN
	149	IF(lwp) WRITE(numout,*)
	150	IF(lwp) WRITE(numout,*) 'cla_dynspg : cross land advection on (ua,va) '
	151	IF(lwp) WRITE(numout,*) '~~~~~~~~~~'
	152	ENDIF
	153	!
	154	IF( nbab == 1 ) CALL cla_bab_el_mandeb('spg') ! Bab el Mandeb strait
	155	IF( ngib == 1 ) CALL cla_gibraltar ('spg') ! Gibraltar strait
	156	IF( nhor == 1 ) CALL cla_hormuz ('spg') ! Hormuz Strait ( Persian Gulf)
	157	!
	158	!!gm lbc is needed here, not?
	159	!!gm CALL lbc_lnk( hdivn, 'U', -1. ) ; CALL lbc_lnk( hdivn, 'V', -1. ) ! Lateral boundary conditions
	160	!
	161	END SUBROUTINE cla_dynspg
[3]	162
	163
[2392]	164	SUBROUTINE cla_init
	165	!! -------------------------------------------------------------------
	166	!! * ROUTINE cla_init *
	167	!!
	168	!! ** Purpose : control check for mpp computation
	169	!!
	170	!! ** Method : - All the strait grid-points must be inside one of the
	171	!! local domain interior for the cla advection to work
	172	!! properly in mpp (i.e. inside (2:jpim1,2:jpjm1) ).
	173	!! Define the corresponding indicators (nbab, ngib, nhor)
	174	!! - The profiles of cross-land fluxes are currently hard
	175	!! coded for L31 levels. Stop if jpk/=31
	176	!!
	177	!! ** Action : nbab, ngib, nhor strait inside the local domain or not
	178	!!---------------------------------------------------------------------
	179	REAL(wp) :: ztemp
	180	!!---------------------------------------------------------------------
	181	!
	182	IF(lwp) WRITE(numout,*)
	183	IF(lwp) WRITE(numout,*) 'cla_init : cross land advection initialisation '
	184	IF(lwp) WRITE(numout,*) '~~~~~~~~~'
	185	!
	186	IF( .NOT.lk_dynspg_flt ) CALL ctl_stop( 'cla_init: Cross Land Advection works only with lk_dynspg_flt=T ' )
	187	!
	188	IF( lk_vvl ) CALL ctl_stop( 'cla_init: Cross Land Advection does not work with lk_vvl=T option' )
	189	!
	190	IF( jpk /= 31 ) CALL ctl_stop( 'cla_init: Cross Land Advection hard coded for ORCA_R2_L31' )
	191	!
	192	! _\|_______\|_______\|_
	193	! 89 \| \|///////\|
	194	! _\|_______\|_______\|_
	195	! ------------------------ ! 88 \|///////\| \|
	196	! Bab el Mandeb strait ! _\|_______\|_______\|_
	197	! ------------------------ ! 87 \|///////\| \|
	198	! _\|_______\|_______\|_
	199	! \| 160 \| 161 \|
	200	!
	201	! The 6 Bab el Mandeb grid-points must be inside one of the interior of the
	202	! local domain for the cla advection to work properly (i.e. (2:jpim1,2:jpjm1)
	203	nbab = 0
	204	IF( ( 1 <= mj0( 88) .AND. mj1( 89) <= jpj ) .AND. & !* (161,89), (161,88) and (161,88) on the local pocessor
	205	& ( 1 <= mi0(160) .AND. mi1(161) <= jpi ) ) nbab = 1
	206	!
	207	! test if there is no local domain that includes all required grid-points
	208	ztemp = REAL( nbab )
	209	IF( lk_mpp ) CALL mpp_sum( ztemp ) ! sum with other processors value
	210	IF( ztemp == 0 ) THEN ! Only 2 points in each direction, this should never be a problem
	211	CALL ctl_stop( ' cross land advection at Bab-el_Mandeb does not work with your processor cutting: change it' )
	212	ENDIF
	213	! ___________________________
	214	! ------------------------ ! 102 \| \|///////\| \|
	215	! Gibraltar strait ! _\|_______\|_______\|_______\|_
	216	! ------------------------ ! 101 \| \|///////\| \|
	217	! _\|_______\|_______\|_______\|_
	218	! \| 139 \| 140 \| 141 \|
	219	!
	220	! The 6 Gibraltar grid-points must be inside one of the interior of the
	221	! local domain for the cla advection to work properly (i.e. (2:jpim1,2:jpjm1)
	222	ngib = 0
	223	IF( ( 2 <= mj0(101) .AND. mj1(102) <= jpjm1 ) .AND. & !* (139:141,101:102) on the local pocessor
	224	& ( 2 <= mi0(139) .AND. mi1(141) <= jpim1 ) ) ngib = 1
	225	!
	226	! test if there is no local domain that includes all required grid-points
	227	ztemp = REAL( ngib )
	228	IF( lk_mpp ) CALL mpp_sum( ztemp ) ! sum with other processors value
	229	IF( ztemp == 0 ) THEN ! 3 points in i-direction, this may be a problem with some cutting
	230	CALL ctl_stop( ' cross land advection at Gibraltar does not work with your processor cutting: change it' )
	231	ENDIF
	232	! _______________
	233	! ------------------------ ! 94 \|/////\| \|
	234	! Hormuz strait ! _\|_____\|_____\|_
	235	! ------------------------ ! 171 172
	236	!
	237	! The 2 Hormuz grid-points must be inside one of the interior of the
	238	! local domain for the cla advection to work properly (i.e. (2:jpim1,2:jpjm1)
	239	nhor = 0
	240	IF( 2 <= mj0( 94) .AND. mj1( 94) <= jpjm1 .AND. &
	241	& 2 <= mi0(171) .AND. mi1(172) <= jpim1 ) nhor = 1
	242	!
	243	! test if there is no local domain that includes all required grid-points
	244	ztemp = REAL( nhor )
	245	IF( lk_mpp ) CALL mpp_sum( ztemp ) ! sum with other processors value
	246	IF( ztemp == 0 ) THEN ! 3 points in i-direction, this may be a problem with some cutting
	247	CALL ctl_stop( ' cross land advection at Hormuz does not work with your processor cutting: change it' )
	248	ENDIF
	249	!
	250	END SUBROUTINE cla_init
[3]	251
	252
[2392]	253	SUBROUTINE cla_bab_el_mandeb( cd_td )
	254	!!----------------------------------------------------------------------
	255	!! * ROUTINE cla_bab_el_mandeb *
	256	!!
	257	!! ** Purpose : update the now horizontal divergence, the tracer tendancy
	258	!! and the after velocity in vicinity of Bab el Mandeb ( Red Sea - Indian ocean).
[3]	259	!!
[2392]	260	!! ** Method : compute the exchanges at each side of the strait :
[3]	261	!!
[2392]	262	!! surf. zio_flow
	263	!! (+ balance of emp) /\ \|\\\\\\\\\\\\|
	264	!! \|\| \|\\\\\\\\\\\\|
	265	!! deep zio_flow \|\| \|\\\\\\\\\\\\|
	266	!! \| \|\| \|\| \|\\\\\\\\\\\\|
	267	!! 89 \| \|\| \|\| \|\\\\\\\\\\\\|
	268	!! \|__\/_v_\|\|__\|____________
	269	!! !\\\\\\\\\\\\| surf. zio_flow
	270	!! \|\\\\\\\\\\\\|<=== (+ balance of emp)
	271	!! \|\\\\\\\\\\\u
	272	!! 88 \|\\\\\\\\\\\\|<--- deep zrecirc (upper+deep at 2 different levels)
	273	!! \|___________\|__________
	274	!! !\\\\\\\\\\\\|
	275	!! \|\\\\\\\\\\\\| ---\ deep zrecirc (upper+deep)
	276	!! 87 !\\\\\\\\\\\u ===/ + deep zio_flow (all at the same level)
	277	!! !\\\\\\\\\\\\|
	278	!! !___________\|__________
	279	!! 160 161
[3]	280	!!
[2392]	281	!!----------------------------------------------------------------------
	282	CHARACTER(len=1), INTENT(in) :: cd_td ! ='div' update the divergence
	283	! ! ='tra' update the tracers
	284	! ! ='spg' update after velocity
	285	INTEGER :: ji, jj, jk ! dummy loop indices
	286	REAL(wp) :: zemp_red ! temporary scalar
	287	REAL(wp) :: zio_flow, zrecirc_upp, zrecirc_mid, zrecirc_bot
[3]	288	!!---------------------------------------------------------------------
[2392]	289	!
	290	SELECT CASE( cd_td )
	291	! ! ---------------- !
	292	CASE( 'ini' ) ! initialisation !
	293	! ! ---------------- !
	294	!
	295	zio_flow = 0.4e6 ! imposed in/out flow
	296	zrecirc_upp = 0.2e6 ! imposed upper recirculation water
	297	zrecirc_bot = 0.5e6 ! imposed bottom recirculation water
[3]	298
[2392]	299	hdiv_161_88(:) = 0.e0 ! (161,88) Gulf of Aden side, north point
	300	hdiv_161_87(:) = 0.e0 ! (161,87) Gulf of Aden side, south point
	301	hdiv_160_89(:) = 0.e0 ! (160,89) Red sea side
[3]	302
[2392]	303	DO jj = mj0(88), mj1(88) !** profile of hdiv at (161,88) (Gulf of Aden side, north point)
	304	DO ji = mi0(161), mi1(161) !------------------------------
	305	DO jk = 1, 8 ! surface in/out flow (Ind -> Red) (div >0)
	306	hdiv_161_88(jk) = + zio_flow / ( 8. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	307	END DO
	308	! ! recirculation water (Ind -> Red) (div >0)
	309	hdiv_161_88(20) = + zrecirc_upp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,20) )
	310	hdiv_161_88(21) = + ( zrecirc_bot - zrecirc_upp ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,21) )
	311	END DO
[3]	312	END DO
[2392]	313	!
	314	DO jj = mj0(87), mj1(87) !** profile of hdiv at (161,88) (Gulf of Aden side, south point)
	315	DO ji = mi0(161), mi1(161) !------------------------------
	316	! ! deep out flow + recirculation (Red -> Ind) (div <0)
	317	hdiv_161_87(21) = - ( zio_flow + zrecirc_bot ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,21) )
[3]	318	END DO
	319	END DO
[2392]	320	!
	321	DO jj = mj0(89), mj1(89) !** profile of hdiv at (161,88) (Red sea side)
	322	DO ji = mi0(160), mi1(160) !------------------------------
	323	DO jk = 1, 8 ! surface inflow (Ind -> Red) (div <0)
	324	hdiv_160_89(jk) = - zio_flow / ( 8. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	325	END DO
	326	! ! deep outflow (Red -> Ind) (div >0)
	327	hdiv_160_89(16) = + zio_flow / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,16) )
[3]	328	END DO
	329	END DO
[2392]	330	! ! ---------------- !
	331	CASE( 'div' ) ! update hdivn ! (call by divcur module)
	332	! ! ---------=====-- !
	333	! !** emp on the Red Sea (div >0)
	334	zemp_red = 0.e0 !---------------------
	335	DO jj = mj0(87), mj1(96) ! sum over the Red sea
	336	DO ji = mi0(148), mi1(160)
	337	zemp_red = zemp_red + emp(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
	338	END DO
	339	END DO
	340	IF( lk_mpp ) CALL mpp_sum( zemp_red ) ! sum with other processors value
	341	zemp_red = zemp_red * 1.e-3 ! convert in m3
	342	!
	343	! !** Correct hdivn (including emp adjustment)
	344	! !-------------------------------------------
	345	DO jj = mj0(88), mj1(88) !* profile of hdiv at (161,88) (Gulf of Aden side, north point)
[3]	346	DO ji = mi0(161), mi1(161)
[2392]	347	hdiv_161_88_kt(:) = hdiv_161_88(:)
	348	DO jk = 1, 8 ! increase the inflow from the Indian (div >0)
	349	hdiv_161_88_kt(jk) = hdiv_161_88(jk) + zemp_red / (8. * e2u(ji,jj) * fse3u(ji,jj,jk) )
	350	END DO
	351	hdivn(ji,jj,:) = hdivn(ji,jj,:) + hdiv_161_88_kt(:)
[3]	352	END DO
	353	END DO
[2392]	354	DO jj = mj0(87), mj1(87) !* profile of divergence at (161,87) (Gulf of Aden side, south point)
[3]	355	DO ji = mi0(161), mi1(161)
[2392]	356	hdivn(ji,jj,:) = hdivn(ji,jj,:) + hdiv_161_87(:)
[3]	357	END DO
	358	END DO
[2392]	359	DO jj = mj0(89), mj1(89) !* profile of divergence at (160,89) (Red sea side)
	360	DO ji = mi0(160), mi1(160)
	361	hdiv_160_89_kt(:) = hdiv_160_89(:)
	362	DO jk = 1, 18 ! increase the inflow from the Indian (div <0)
	363	hdiv_160_89_kt(jk) = hdiv_160_89(jk) - zemp_red / (10. * e1v(ji,jj) * fse3v(ji,jj,jk) )
	364	END DO
	365	hdivn(ji, jj,:) = hdivn(ji, jj,:) + hdiv_160_89_kt(:)
	366	END DO
[3]	367	END DO
[2392]	368	! ! ---------------- !
	369	CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module)
	370	! ! --------=======- !
	371	!
	372	DO jj = mj0(88), mj1(88) !** (161,88) (Gulf of Aden side, north point)
[3]	373	DO ji = mi0(161), mi1(161)
[2392]	374	DO jk = 1, jpkm1 ! surf inflow + reciculation (from Gulf of Aden)
	375	ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_161_88_kt(jk) * tn(ji,jj,jk)
	376	sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_161_88_kt(jk) * sn(ji,jj,jk)
	377	END DO
[3]	378	END DO
	379	END DO
[2392]	380	DO jj = mj0(87), mj1(87) !** (161,87) (Gulf of Aden side, south point)
[3]	381	DO ji = mi0(161), mi1(161)
[2392]	382	jk = 21 ! deep outflow + recirulation (combined flux)
	383	ta(ji,jj,jk) = ta(ji,jj,jk) + hdiv_161_88(20) * tn(ji ,jj+1,20) & ! upper recirculation from Gulf of Aden
	384	& + hdiv_161_88(21) * tn(ji ,jj+1,21) & ! deep recirculation from Gulf of Aden
	385	& + hdiv_160_89(16) * tn(ji-1,jj+2,16) ! deep inflow from Red sea
	386	sa(ji,jj,jk) = sa(ji,jj,jk) + hdiv_161_88(20) * sn(ji ,jj+1,20) &
	387	& + hdiv_161_88(21) * sn(ji ,jj+1,21) &
	388	& + hdiv_160_89(16) * sn(ji-1,jj+2,16)
[3]	389	END DO
	390	END DO
[2392]	391	DO jj = mj0(89), mj1(89) !** (161,88) (Red sea side)
	392	DO ji = mi0(160), mi1(160)
	393	DO jk = 1, 14 ! surface inflow (from Gulf of Aden)
	394	ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_160_89_kt(jk) * tn(ji+1,jj-1,jk)
	395	sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_160_89_kt(jk) * sn(ji+1,jj-1,jk)
	396	END DO
	397	! ! deep outflow (from Red sea)
	398	ta(ji,jj,16) = ta(ji,jj,16) - hdiv_160_89(jk) * tn(ji,jj,jk)
	399	sa(ji,jj,16) = sa(ji,jj,16) - hdiv_160_89(jk) * sn(ji,jj,jk)
	400	END DO
[3]	401	END DO
[2392]	402	!
	403	! ! ---------------- !
	404	CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module)
	405	! ! --------=======- !
	406	! at this stage, (ua,va) are the after velocity, not the tendancy
	407	! compute the velocity from the divergence at T-point
	408	!
	409	DO jj = mj0(88), mj1(88) !** (160,88) (Gulf of Aden side, north point)
	410	DO ji = mi0(160), mi1(160) ! 160, not 161 as it is a U-point)
	411	ua(ji,jj,:) = - hdiv_161_88_kt(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
	412	& * e2u(ji,jj) * fse3u(ji,jj,:)
[3]	413	END DO
	414	END DO
[2392]	415	DO jj = mj0(87), mj1(87) !** (160,87) (Gulf of Aden side, south point)
	416	DO ji = mi0(160), mi1(160) ! 160, not 161 as it is a U-point)
	417	ua(ji,jj,:) = - hdiv_161_87(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
	418	& * e2u(ji,jj) * fse3u(ji,jj,:)
	419	END DO
[3]	420	END DO
[2392]	421	DO jj = mj0(88), mj1(88) !** profile of divergence at (160,89) (Red sea side)
	422	DO ji = mi0(160), mi1(160) ! 88, not 89 as it is a V-point)
	423	va(ji,jj,:) = - hdiv_160_89_kt(:) / ( e1t(ji,jj+1) * e2t(ji,jj+1) * fse3t(ji,jj+1,:) ) &
	424	& * e1v(ji,jj) * fse3v(ji,jj,:)
	425	END DO
	426	END DO
	427	END SELECT
	428	!
	429	END SUBROUTINE cla_bab_el_mandeb
	430
[3]	431
[2392]	432	SUBROUTINE cla_gibraltar( cd_td )
	433	!! -------------------------------------------------------------------
	434	!! * ROUTINE cla_gibraltar *
	435	!!
	436	!! ** Purpose : update the now horizontal divergence, the tracer
	437	!! tendancyand the after velocity in vicinity of Gibraltar
	438	!! strait ( Persian Gulf - Indian ocean ).
[3]	439	!!
	440	!! ** Method :
[2392]	441	!! _______________________
	442	!! deep zio_flow /====\|///////\|====> surf. zio_flow
	443	!! + deep zrecirc \----\|///////\| (+balance of emp)
	444	!! 102 u///////u
	445	!! mid. recicul <--\|///////\|<==== deep zio_flow
	446	!! _____\|_______\|_____
	447	!! surf. zio_flow ====>\|///////\|
	448	!! (+balance of emp) \|///////\|
	449	!! 101 u///////\|
	450	!! mid. recicul -->\|///////\| Caution: zrecirc split into
	451	!! deep zrecirc ---->\|///////\| upper & bottom recirculation
	452	!! _______\|_______\|_______
	453	!! 139 140 141
[3]	454	!!
	455	!!---------------------------------------------------------------------
[2392]	456	CHARACTER(len=1), INTENT(in) :: cd_td ! ='div' update the divergence
	457	! ! ='tra' update the tracers
	458	! ! ='spg' update after velocity
	459	INTEGER :: ji, jj, jk ! dummy loop indices
	460	REAL(wp) :: zemp_med ! temporary scalar
	461	REAL(wp) :: zio_flow, zrecirc_upp, zrecirc_mid, zrecirc_bot
[3]	462	!!---------------------------------------------------------------------
[2392]	463	!
	464	SELECT CASE( cd_td )
	465	! ! ---------------- !
	466	CASE( 'ini' ) ! initialisation !
	467	! ! ---------------- !
	468	! !** initialization of the velocity
	469	hdiv_139_101(:) = 0.e0 ! 139,101 (Atlantic side, south point)
	470	hdiv_139_102(:) = 0.e0 ! 139,102 (Atlantic side, north point)
	471	hdiv_141_102(:) = 0.e0 ! 141,102 (Med sea side)
	472
	473	! !** imposed transport
	474	zio_flow = 0.8e6 ! inflow surface water
	475	zrecirc_mid = 0.7e6 ! middle recirculation water
	476	zrecirc_upp = 2.5e6 ! upper recirculation water
	477	zrecirc_bot = 3.5e6 ! bottom recirculation water
	478	!
	479	DO jj = mj0(101), mj1(101) !** profile of hdiv at 139,101 (Atlantic side, south point)
	480	DO ji = mi0(139), mi1(139) !-----------------------------
	481	DO jk = 1, 14 ! surface in/out flow (Atl -> Med) (div >0)
	482	hdiv_139_101(jk) = + zio_flow / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	483	END DO
	484	DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div >0)
	485	hdiv_139_101(jk) = + zrecirc_mid / ( 6. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	486	END DO
	487	! ! upper reciculation (Atl 101 -> Atl 101) (div >0)
	488	hdiv_139_101(21) = + zrecirc_upp / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	489	!
	490	! ! upper & bottom reciculation (Atl 101 -> Atl 101 & 102) (div >0)
	491	hdiv_139_101(22) = ( zrecirc_bot - zrecirc_upp ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	492	END DO
[3]	493	END DO
[2392]	494	DO jj = mj0(102), mj1(102) !** profile of hdiv at 139,102 (Atlantic side, north point)
	495	DO ji = mi0(139), mi1(139) !-----------------------------
	496	DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div <0)
	497	hdiv_139_102(jk) = - zrecirc_mid / ( 6. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	498	END DO
	499	! ! outflow of Mediterranean sea + deep recirculation (div <0)
	500	hdiv_139_102(22) = - ( zio_flow + zrecirc_bot ) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	501	END DO
[3]	502	END DO
[2392]	503	DO jj = mj0(102), mj1(102) !** velocity profile at 141,102 (Med sea side)
	504	DO ji = mi0(141), mi1(141) !------------------------------
	505	DO jk = 1, 14 ! surface inflow in the Med (div <0)
	506	hdiv_141_102(jk) = - zio_flow / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	507	END DO
	508	! ! deep outflow toward the Atlantic (div >0)
	509	hdiv_141_102(21) = + zio_flow / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
[3]	510	END DO
	511	END DO
[2392]	512	! ! ---------------- !
	513	CASE( 'div' ) ! update hdivn ! (call by divcur module)
	514	! ! ---------=====-- !
	515	! !** emp on the Mediterranean Sea (div >0)
	516	zemp_med = 0.e0 !-------------------------------
	517	DO jj = mj0(96), mj1(110) ! sum over the Med sea
	518	DO ji = mi0(141),mi1(181)
	519	zemp_med = zemp_med + emp(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
[3]	520	END DO
	521	END DO
[2392]	522	DO jj = mj0(96), mj1(96) ! minus 2 points in Red Sea
	523	DO ji = mi0(148),mi1(148)
	524	zemp_med = zemp_med - emp(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
[3]	525	END DO
[2392]	526	DO ji = mi0(149),mi1(149)
	527	zemp_med = zemp_med - emp(ji,jj) * e1t(ji,jj) * e2t(ji,jj) * tmask_i(ji,jj)
	528	END DO
[3]	529	END DO
[2392]	530	IF( lk_mpp ) CALL mpp_sum( zemp_med ) ! sum with other processors value
	531	zemp_med = zemp_med * 1.e-3 ! convert in m3
	532	!
	533	! !** Correct hdivn (including emp adjustment)
	534	! !-------------------------------------------
	535	DO jj = mj0(101), mj1(101) !* 139,101 (Atlantic side, south point)
[3]	536	DO ji = mi0(139), mi1(139)
[2392]	537	hdiv_139_101_kt(:) = hdiv_139_101(:)
	538	DO jk = 1, 14 ! increase the inflow from the Atlantic (div >0)
	539	hdiv_139_101_kt(jk) = hdiv_139_101(jk) + zemp_med / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	540	END DO
	541	hdivn(ji, jj,:) = hdivn(ji, jj,:) + hdiv_139_101_kt(:)
[3]	542	END DO
	543	END DO
[2392]	544	DO jj = mj0(102), mj1(102) !* 139,102 (Atlantic side, north point)
[3]	545	DO ji = mi0(139), mi1(139)
[2392]	546	hdivn(ji,jj,:) = hdivn(ji,jj,:) + hdiv_139_102(:)
[3]	547	END DO
	548	END DO
[2392]	549	DO jj = mj0(102), mj1(102) !* 141,102 (Med side)
	550	DO ji = mi0(141), mi1(141)
	551	hdiv_141_102(:) = hdiv_141_102(:)
	552	DO jk = 1, 14 ! increase the inflow from the Atlantic (div <0)
	553	hdiv_141_102_kt(jk) = hdiv_141_102(jk) - zemp_med / ( 14. * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	554	END DO
	555	hdivn(ji, jj,:) = hdivn(ji, jj,:) + hdiv_141_102_kt(:)
	556	END DO
[3]	557	END DO
[2392]	558	! ! ---------------- !
	559	CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module)
	560	! ! --------=======- !
	561	!
	562	DO jj = mj0(101), mj1(101) !** 139,101 (Atlantic side, south point) (div >0)
[3]	563	DO ji = mi0(139), mi1(139)
[2392]	564	DO jk = 1, jpkm1 ! surf inflow + mid. & bottom reciculation (from Atlantic)
	565	ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_139_101_kt(jk) * tn(ji,jj,jk)
	566	sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_139_101_kt(jk) * sn(ji,jj,jk)
	567	END DO
[3]	568	END DO
	569	END DO
[2392]	570	!
	571	DO jj = mj0(102), mj1(102) !** 139,102 (Atlantic side, north point) (div <0)
[3]	572	DO ji = mi0(139), mi1(139)
[2392]	573	DO jk = 15, 20 ! middle reciculation (Atl 101 -> Atl 102) (div <0)
	574	ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_139_102(jk) * tn(ji,jj-1,jk) ! middle Atlantic recirculation
	575	sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_139_102(jk) * sn(ji,jj-1,jk)
	576	END DO
	577	! ! upper & bottom Atl. reciculation (Atl 101 -> Atl 102) - (div <0)
	578	! ! deep Med flow (Med 102 -> Atl 102) - (div <0)
	579	ta(ji,jj,22) = ta(ji,jj,22) + hdiv_141_102(21) * tn(ji+2,jj ,21) & ! deep Med flow
	580	& + hdiv_139_101(21) * tn(ji ,jj-1,21) & ! upper Atlantic recirculation
	581	& + hdiv_139_101(22) * tn(ji ,jj-1,22) ! bottom Atlantic recirculation
	582	sa(ji,jj,22) = sa(ji,jj,22) + hdiv_141_102(21) * sn(ji+2,jj ,21) &
	583	& + hdiv_139_101(21) * sn(ji ,jj-1,21) &
	584	& + hdiv_139_101(22) * sn(ji ,jj-1,22)
[3]	585	END DO
	586	END DO
[2392]	587	DO jj = mj0(102), mj1(102) !* 141,102 (Med side) (div <0)
[3]	588	DO ji = mi0(141), mi1(141)
[2392]	589	DO jk = 1, 14 ! surface flow from Atlantic to Med sea
	590	ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_141_102_kt(jk) * tn(ji-2,jj-1,jk)
	591	sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_141_102_kt(jk) * sn(ji-2,jj-1,jk)
	592	END DO
	593	! ! deeper flow from Med sea to Atlantic
	594	ta(ji,jj,21) = ta(ji,jj,21) - hdiv_141_102(21) * tn(ji,jj,21)
	595	sa(ji,jj,21) = sa(ji,jj,21) - hdiv_141_102(21) * sn(ji,jj,21)
[3]	596	END DO
	597	END DO
[2392]	598	! ! ---------------- !
	599	CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module)
	600	! ! --------=======- !
	601	! at this stage, (ua,va) are the after velocity, not the tendancy
	602	! compute the velocity from the divergence at T-point
	603	!
	604	DO jj = mj0(101), mj1(101) !** 139,101 (Atlantic side, south point)
	605	DO ji = mi0(139), mi1(139) ! div >0 => ua >0, same sign
	606	ua(ji,jj,:) = hdiv_139_101_kt(:) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,:) ) &
	607	& * e2u(ji,jj) * fse3u(ji,jj,:)
	608	END DO
[3]	609	END DO
[2392]	610	DO jj = mj0(102), mj1(102) !** 139,102 (Atlantic side, north point)
	611	DO ji = mi0(139), mi1(139) ! div <0 => ua <0, same sign
	612	ua(ji,jj,:) = hdiv_139_102(:) / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,:) ) &
	613	& * e2u(ji,jj) * fse3u(ji,jj,:)
[3]	614	END DO
	615	END DO
[2392]	616	DO jj = mj0(102), mj1(102) !** 140,102 (Med side) (140 not 141 as it is a U-point)
	617	DO ji = mi0(140), mi1(140) ! div >0 => ua <0, opposite sign
	618	ua(ji,jj,:) = - hdiv_141_102(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
	619	& * e2u(ji,jj) * fse3u(ji,jj,:)
[3]	620	END DO
	621	END DO
[2392]	622	!
	623	END SELECT
	624	!
	625	END SUBROUTINE cla_gibraltar
[3]	626
	627
[2392]	628	SUBROUTINE cla_hormuz( cd_td )
	629	!! -------------------------------------------------------------------
	630	!! * ROUTINE div_hormuz *
	631	!!
	632	!! ** Purpose : update the now horizontal divergence, the tracer
	633	!! tendancyand the after velocity in vicinity of Hormuz
	634	!! strait ( Persian Gulf - Indian ocean ).
	635	!!
	636	!! ** Method : Hormuz strait
	637	!! ______________
	638	!! \|/////\|<== surface inflow
	639	!! 94 \|/////\|
	640	!! \|/////\|==> deep outflow
	641	!! \|_____\|_______
	642	!! 171 172
[3]	643	!!---------------------------------------------------------------------
[2392]	644	CHARACTER(len=1), INTENT(in) :: cd_td ! ='ini' initialisation
	645	!! ! ='div' update the divergence
	646	!! ! ='tra' update the tracers
	647	!! ! ='spg' update after velocity
[3]	648	!!
[2392]	649	INTEGER :: ji, jj, jk ! dummy loop indices
	650	REAL(wp) :: zio_flow ! temporary scalar
[3]	651	!!---------------------------------------------------------------------
[2392]	652	!
	653	SELECT CASE( cd_td )
	654	! ! ---------------- !
	655	CASE( 'ini' ) ! initialisation !
	656	! ! ---------------- !
	657	! !** profile of horizontal divergence due to cross-land advection
	658	zio_flow = 1.e6 ! imposed in/out flow
	659	!
	660	hdiv_172_94(:) = 0.e0
	661	!
	662	DO jj = mj0(94), mj1(94) ! in/out flow at (i,j) = (172,94)
	663	DO ji = mi0(172), mi1(172)
	664	DO jk = 1, 8 ! surface inflow (Indian ocean to Persian Gulf) (div<0)
	665	hdiv_172_94(jk) = - ( zio_flow / 8.e0 * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	666	END DO
	667	DO jk = 16, 18 ! deep outflow (Persian Gulf to Indian ocean) (div>0)
	668	hdiv_172_94(jk) = + ( zio_flow / 3.e0 * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
	669	END DO
[3]	670	END DO
	671	END DO
[2392]	672	! !** T & S profile in the Hormuz strait (use in deep outflow)
	673	! Temperature and Salinity
	674	t_171_94_hor(:) = 0.e0 ; s_171_94_hor(:) = 0.e0
	675	t_171_94_hor(16) = 18.4 ; s_171_94_hor(16) = 36.27
	676	t_171_94_hor(17) = 17.8 ; s_171_94_hor(17) = 36.4
	677	t_171_94_hor(18) = 16. ; s_171_94_hor(18) = 36.27
	678	!
	679	! ! ---------------- !
	680	CASE( 'div' ) ! update hdivn ! (call by divcur module)
	681	! ! ---------=====-- !
	682	!
	683	DO jj = mj0(94), mj1(94) !** 172,94 (Indian ocean side)
	684	DO ji = mi0(172), mi1(172)
	685	hdivn(ji,jj,:) = hdivn(ji,jj,:) + hdiv_172_94(:)
[3]	686	END DO
	687	END DO
[2392]	688	! ! ---------------- !
	689	CASE( 'tra' ) ! update (ta,sa) ! (call by traadv module)
	690	! ! --------=======- !
	691	!
	692	DO jj = mj0(94), mj1(94) !** 172,94 (Indian ocean side)
[3]	693	DO ji = mi0(172), mi1(172)
[2392]	694	DO jk = 1, 8 ! surface inflow (Indian ocean to Persian Gulf) (div<0)
	695	ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_172_94(jk) * tn(ji,jj,jk)
	696	sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_172_94(jk) * sn(ji,jj,jk)
	697	END DO
	698	DO jk = 16, 18 ! deep outflow (Persian Gulf to Indian ocean) (div>0)
	699	ta(ji,jj,jk) = ta(ji,jj,jk) - hdiv_172_94(jk) * t_171_94_hor(jk)
	700	sa(ji,jj,jk) = sa(ji,jj,jk) - hdiv_172_94(jk) * s_171_94_hor(jk)
	701	END DO
[3]	702	END DO
	703	END DO
[2392]	704	! ! ---------------- !
	705	CASE( 'spg' ) ! update (ua,va) ! (call by dynspg module)
	706	! ! --------=======- !
	707	! No barotropic flow through Hormuz strait
	708	! at this stage, (ua,va) are the after velocity, not the tendancy
	709	! compute the velocity from the divergence at T-point
	710	DO jj = mj0(94), mj1(94) !** 171,94 (Indian ocean side) (171 not 172 as it is the western U-point)
	711	DO ji = mi0(171), mi1(171) ! div >0 => ua >0, opposite sign
	712	ua(ji,jj,:) = - hdiv_172_94(:) / ( e1t(ji+1,jj) * e2t(ji+1,jj) * fse3t(ji+1,jj,:) ) &
	713	& * e2u(ji,jj) * fse3u(ji,jj,:)
[3]	714	END DO
	715	END DO
[2392]	716	!
	717	END SELECT
[3]	718	!
[2392]	719	END SUBROUTINE cla_hormuz
	720
[3]	721	#else
	722	!!----------------------------------------------------------------------
[2392]	723	!! Default key Dummy module
[3]	724	!!----------------------------------------------------------------------
[2392]	725	USE in_out_manager ! I/O manager
[3]	726	CONTAINS
[2392]	727	SUBROUTINE cla_init
	728	CALL ctl_stop( 'cla_init: Cross Land Advection hard coded for ORCA_R2 with 31 levels' )
	729	END SUBROUTINE cla_init
	730	SUBROUTINE cla_div( kt )
	731	WRITE(,) 'cla_div: You should have not see this print! error?', kt
	732	END SUBROUTINE cla_div
	733	SUBROUTINE cla_traadv( kt )
	734	WRITE(,) 'cla_traadv: You should have not see this print! error?', kt
	735	END SUBROUTINE cla_traadv
	736	SUBROUTINE cla_dynspg( kt )
	737	WRITE(,) 'dyn_spg_cla: You should have not see this print! error?', kt
	738	END SUBROUTINE cla_dynspg
[3]	739	#endif
[2392]	740
[3]	741	!!======================================================================
	742	END MODULE cla

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

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

Context Navigation

source: branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/cla.F90 @ 2392

Download in other formats: