Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

agrif_opa_sponge.F90 @ 7910

Last change on this file since 7910 was 7910, checked in by timgraham, 7 years ago
All wrk_alloc removed
Property svn:keywords set to `Id`
File size: 16.3 KB

Rev	Line
[3680]	1	#define SPONGE && define SPONGE_TOP
[390]	2
[5656]	3	MODULE agrif_opa_sponge
[6140]	4	!!======================================================================
	5	!! * MODULE agrif_opa_update *
	6	!! AGRIF :
	7	!!======================================================================
	8	!! History :
	9	!!----------------------------------------------------------------------
[7646]	10	#if defined key_agrif
[636]	11	USE par_oce
	12	USE oce
	13	USE dom_oce
	14	USE in_out_manager
[782]	15	USE agrif_oce
[5656]	16	USE lbclnk ! ocean lateral boundary conditions (or mpp link)
[390]	17
[636]	18	IMPLICIT NONE
	19	PRIVATE
[390]	20
[5656]	21	PUBLIC Agrif_Sponge, Agrif_Sponge_Tra, Agrif_Sponge_Dyn
	22	PUBLIC interptsn_sponge, interpun_sponge, interpvn_sponge
[390]	23
[1156]	24	!!----------------------------------------------------------------------
[6140]	25	!! NEMO/NST 3.7 , NEMO Consortium (2015)
[1156]	26	!! $Id$
[2528]	27	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1156]	28	!!----------------------------------------------------------------------
[5656]	29	CONTAINS
[636]	30
[782]	31	SUBROUTINE Agrif_Sponge_Tra
[636]	32	!!---------------------------------------------
	33	!! * ROUTINE Agrif_Sponge_Tra *
	34	!!---------------------------------------------
	35	REAL(wp) :: timecoeff
[6140]	36	!!---------------------------------------------
	37	!
[390]	38	#if defined SPONGE
[636]	39	timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot()
	40
[5656]	41	CALL Agrif_Sponge
[390]	42	Agrif_SpecialValue=0.
	43	Agrif_UseSpecialValue = .TRUE.
[5656]	44	tabspongedone_tsn = .FALSE.
[390]	45
[5656]	46	CALL Agrif_Bc_Variable(tsn_sponge_id,calledweight=timecoeff,procname=interptsn_sponge)
[390]	47
[5656]	48	Agrif_UseSpecialValue = .FALSE.
[390]	49	#endif
[6140]	50	!
[636]	51	END SUBROUTINE Agrif_Sponge_Tra
[390]	52
[6140]	53
[782]	54	SUBROUTINE Agrif_Sponge_dyn
[636]	55	!!---------------------------------------------
	56	!! * ROUTINE Agrif_Sponge_dyn *
	57	!!---------------------------------------------
	58	REAL(wp) :: timecoeff
[6140]	59	!!---------------------------------------------
[390]	60
	61	#if defined SPONGE
[636]	62	timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot()
[390]	63
	64	Agrif_SpecialValue=0.
[782]	65	Agrif_UseSpecialValue = ln_spc_dyn
[390]	66
[5656]	67	tabspongedone_u = .FALSE.
	68	tabspongedone_v = .FALSE.
	69	CALL Agrif_Bc_Variable(un_sponge_id,calledweight=timecoeff,procname=interpun_sponge)
[390]	70
[5656]	71	tabspongedone_u = .FALSE.
	72	tabspongedone_v = .FALSE.
	73	CALL Agrif_Bc_Variable(vn_sponge_id,calledweight=timecoeff,procname=interpvn_sponge)
	74
[390]	75	Agrif_UseSpecialValue = .FALSE.
	76	#endif
[6140]	77	!
[636]	78	END SUBROUTINE Agrif_Sponge_dyn
[390]	79
[6140]	80
[3680]	81	SUBROUTINE Agrif_Sponge
	82	!!---------------------------------------------
	83	!! * ROUTINE Agrif_Sponge *
	84	!!---------------------------------------------
	85	INTEGER :: ji,jj,jk
	86	INTEGER :: ispongearea, ilci, ilcj
[4153]	87	LOGICAL :: ll_spdone
	88	REAL(wp) :: z1spongearea, zramp
[7910]	89	REAL(wp), DIMENSION(jpi,jpj) :: ztabramp
[3680]	90
	91	#if defined SPONGE \|\| defined SPONGE_TOP
[4153]	92	ll_spdone=.TRUE.
	93	IF (( .NOT. spongedoneT ).OR.( .NOT. spongedoneU )) THEN
	94	! Define ramp from boundaries towards domain interior
	95	! at T-points
	96	! Store it in ztabramp
	97	ll_spdone=.FALSE.
[3680]	98
	99
[5656]	100	ispongearea = 2 + nn_sponge_len * Agrif_irhox()
[4153]	101	ilci = nlci - ispongearea
	102	ilcj = nlcj - ispongearea
	103	z1spongearea = 1._wp / REAL( ispongearea - 2 )
[3680]	104
[5656]	105	ztabramp(:,:) = 0._wp
[4153]	106
	107	IF( (nbondi == -1) .OR. (nbondi == 2) ) THEN
	108	DO jj = 1, jpj
	109	IF ( umask(2,jj,1) == 1._wp ) THEN
	110	DO ji = 2, ispongearea
	111	ztabramp(ji,jj) = ( ispongearea-ji ) * z1spongearea
	112	END DO
	113	ENDIF
	114	ENDDO
	115	ENDIF
	116
	117	IF( (nbondi == 1) .OR. (nbondi == 2) ) THEN
	118	DO jj = 1, jpj
	119	IF ( umask(nlci-2,jj,1) == 1._wp ) THEN
	120	DO ji = ilci+1,nlci-1
	121	zramp = (ji - (ilci+1) ) * z1spongearea
	122	ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp )
	123	ENDDO
	124	ENDIF
	125	ENDDO
	126	ENDIF
	127
	128	IF( (nbondj == -1) .OR. (nbondj == 2) ) THEN
	129	DO ji = 1, jpi
	130	IF ( vmask(ji,2,1) == 1._wp ) THEN
	131	DO jj = 2, ispongearea
	132	zramp = ( ispongearea-jj ) * z1spongearea
	133	ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp )
	134	END DO
	135	ENDIF
	136	ENDDO
	137	ENDIF
	138
	139	IF( (nbondj == 1) .OR. (nbondj == 2) ) THEN
	140	DO ji = 1, jpi
	141	IF ( vmask(ji,nlcj-2,1) == 1._wp ) THEN
	142	DO jj = ilcj+1,nlcj-1
	143	zramp = (jj - (ilcj+1) ) * z1spongearea
	144	ztabramp(ji,jj) = MAX( ztabramp(ji,jj), zramp )
	145	END DO
	146	ENDIF
	147	ENDDO
	148	ENDIF
	149
	150	ENDIF
	151
[3680]	152	! Tracers
	153	IF( .NOT. spongedoneT ) THEN
[5656]	154	fsaht_spu(:,:) = 0._wp
	155	fsaht_spv(:,:) = 0._wp
	156	DO jj = 2, jpjm1
	157	DO ji = 2, jpim1 ! vector opt.
	158	fsaht_spu(ji,jj) = 0.5_wp * visc_tra * (ztabramp(ji,jj) + ztabramp(ji+1,jj ))
	159	fsaht_spv(ji,jj) = 0.5_wp * visc_tra * (ztabramp(ji,jj) + ztabramp(ji ,jj+1))
	160	END DO
	161	END DO
[3680]	162
[5656]	163	CALL lbc_lnk( fsaht_spu, 'U', 1. ) ! Lateral boundary conditions
	164	CALL lbc_lnk( fsaht_spv, 'V', 1. )
[3680]	165	spongedoneT = .TRUE.
	166	ENDIF
	167
	168	! Dynamics
	169	IF( .NOT. spongedoneU ) THEN
[5656]	170	fsahm_spt(:,:) = 0._wp
	171	fsahm_spf(:,:) = 0._wp
	172	DO jj = 2, jpjm1
	173	DO ji = 2, jpim1 ! vector opt.
	174	fsahm_spt(ji,jj) = visc_dyn * ztabramp(ji,jj)
	175	fsahm_spf(ji,jj) = 0.25_wp * visc_dyn * ( ztabramp(ji,jj) + ztabramp(ji ,jj+1) &
	176	& +ztabramp(ji,jj) + ztabramp(ji+1,jj ) )
	177	END DO
	178	END DO
[3680]	179
[5656]	180	CALL lbc_lnk( fsahm_spt, 'T', 1. ) ! Lateral boundary conditions
	181	CALL lbc_lnk( fsahm_spf, 'F', 1. )
[3680]	182	spongedoneU = .TRUE.
	183	ENDIF
	184	!
	185	!
	186	#endif
[6140]	187	!
[3680]	188	END SUBROUTINE Agrif_Sponge
	189
[6140]	190
[5656]	191	SUBROUTINE interptsn_sponge(tabres,i1,i2,j1,j2,k1,k2,n1,n2,before)
[636]	192	!!---------------------------------------------
[5656]	193	!! * ROUTINE interptsn_sponge *
[636]	194	!!---------------------------------------------
[3294]	195	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
	196	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
[5656]	197	LOGICAL, INTENT(in) :: before
[6140]	198	!
[5656]	199	INTEGER :: ji, jj, jk, jn ! dummy loop indices
	200	INTEGER :: iku, ikv
	201	REAL(wp) :: ztsa, zabe1, zabe2, zbtr
	202	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ztu, ztv
	203	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2) ::tsbdiff
	204	!
[6140]	205	IF( before ) THEN
[5656]	206	tabres(i1:i2,j1:j2,k1:k2,n1:n2) = tsn(i1:i2,j1:j2,k1:k2,n1:n2)
	207	ELSE
[6140]	208	!
[5656]	209	tsbdiff(:,:,:,:) = tsb(i1:i2,j1:j2,:,:) - tabres(:,:,:,:)
	210	DO jn = 1, jpts
	211	DO jk = 1, jpkm1
	212	DO jj = j1,j2-1
	213	DO ji = i1,i2-1
[6140]	214	zabe1 = fsaht_spu(ji,jj) * umask(ji,jj,jk) * e2_e1u(ji,jj) * e3u_n(ji,jj,jk)
	215	zabe2 = fsaht_spv(ji,jj) * vmask(ji,jj,jk) * e1_e2v(ji,jj) * e3v_n(ji,jj,jk)
[5656]	216	ztu(ji,jj,jk) = zabe1 * ( tsbdiff(ji+1,jj ,jk,jn) - tsbdiff(ji,jj,jk,jn) )
	217	ztv(ji,jj,jk) = zabe2 * ( tsbdiff(ji ,jj+1,jk,jn) - tsbdiff(ji,jj,jk,jn) )
[6140]	218	END DO
	219	END DO
	220	!
[5656]	221	IF( ln_zps ) THEN ! set gradient at partial step level
	222	DO jj = j1,j2-1
	223	DO ji = i1,i2-1
	224	! last level
	225	iku = mbku(ji,jj)
	226	ikv = mbkv(ji,jj)
[6140]	227	IF( iku == jk ) ztu(ji,jj,jk) = 0._wp
	228	IF( ikv == jk ) ztv(ji,jj,jk) = 0._wp
[5656]	229	END DO
	230	END DO
	231	ENDIF
[6140]	232	END DO
	233	!
[5656]	234	DO jk = 1, jpkm1
	235	DO jj = j1+1,j2-1
	236	DO ji = i1+1,i2-1
	237	IF (.NOT. tabspongedone_tsn(ji,jj)) THEN
[6140]	238	zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk)
[5656]	239	! horizontal diffusive trends
	240	ztsa = zbtr * ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) + ztv(ji,jj,jk) - ztv(ji ,jj-1,jk) )
	241	! add it to the general tracer trends
	242	tsa(ji,jj,jk,jn) = tsa(ji,jj,jk,jn) + ztsa
	243	ENDIF
[6140]	244	END DO
	245	END DO
	246	END DO
	247	!
	248	END DO
	249	!
[5656]	250	tabspongedone_tsn(i1+1:i2-1,j1+1:j2-1) = .TRUE.
[6140]	251	!
[5656]	252	ENDIF
[6140]	253	!
[5656]	254	END SUBROUTINE interptsn_sponge
	255
[6140]	256
[5656]	257	SUBROUTINE interpun_sponge(tabres,i1,i2,j1,j2,k1,k2, before)
[636]	258	!!---------------------------------------------
[5656]	259	!! * ROUTINE interpun_sponge *
	260	!!---------------------------------------------
[636]	261	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2
	262	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
[5656]	263	LOGICAL, INTENT(in) :: before
[390]	264
[5656]	265	INTEGER :: ji,jj,jk
[390]	266
[5656]	267	! sponge parameters
	268	REAL(wp) :: ze2u, ze1v, zua, zva, zbtr
	269	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ubdiff
	270	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: rotdiff, hdivdiff
	271	INTEGER :: jmax
[6140]	272	!!---------------------------------------------
[5656]	273	!
[6140]	274	IF( before ) THEN
[5656]	275	tabres = un(i1:i2,j1:j2,:)
	276	ELSE
	277	ubdiff(i1:i2,j1:j2,:) = (ub(i1:i2,j1:j2,:) - tabres(:,:,:))*umask(i1:i2,j1:j2,:)
[6140]	278	!
[5656]	279	DO jk = 1, jpkm1 ! Horizontal slab
	280	! ! ===============
	281
	282	! ! --------
	283	! Horizontal divergence ! div
	284	! ! --------
	285	DO jj = j1,j2
	286	DO ji = i1+1,i2 ! vector opt.
[6140]	287	zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) * fsahm_spt(ji,jj)
	288	hdivdiff(ji,jj,jk) = ( e2u(ji ,jj)e3u_n(ji ,jj,jk) ubdiff(ji ,jj,jk) &
	289	& -e2u(ji-1,jj)e3u_n(ji-1,jj,jk) ubdiff(ji-1,jj,jk) ) * zbtr
[5656]	290	END DO
	291	END DO
	292
	293	DO jj = j1,j2-1
	294	DO ji = i1,i2 ! vector opt.
[6140]	295	zbtr = r1_e1e2f(ji,jj) * e3f_n(ji,jj,jk) * fsahm_spf(ji,jj)
[5656]	296	rotdiff(ji,jj,jk) = (-e1u(ji,jj+1) * ubdiff(ji,jj+1,jk) &
	297	+e1u(ji,jj ) * ubdiff(ji,jj ,jk) &
	298	& ) * fmask(ji,jj,jk) * zbtr
	299	END DO
	300	END DO
[6140]	301	END DO
[5656]	302	!
	303	DO jj = j1+1, j2-1
	304	DO ji = i1+1, i2-1 ! vector opt.
	305
	306	IF (.NOT. tabspongedone_u(ji,jj)) THEN
	307	DO jk = 1, jpkm1 ! Horizontal slab
	308	ze2u = rotdiff (ji,jj,jk)
	309	ze1v = hdivdiff(ji,jj,jk)
	310	! horizontal diffusive trends
[6140]	311	zua = - ( ze2u - rotdiff (ji,jj-1,jk)) / ( e2u(ji,jj) * e3u_n(ji,jj,jk) ) &
[5656]	312	+ ( hdivdiff(ji+1,jj,jk) - ze1v ) / e1u(ji,jj)
	313
	314	! add it to the general momentum trends
	315	ua(ji,jj,jk) = ua(ji,jj,jk) + zua
	316
	317	END DO
	318	ENDIF
	319
	320	END DO
	321	END DO
	322
	323	tabspongedone_u(i1+1:i2-1,j1+1:j2-1) = .TRUE.
	324
	325	jmax = j2-1
	326	IF ((nbondj == 1).OR.(nbondj == 2)) jmax = MIN(jmax,nlcj-3)
	327
	328	DO jj = j1+1, jmax
	329	DO ji = i1+1, i2 ! vector opt.
	330
	331	IF (.NOT. tabspongedone_v(ji,jj)) THEN
	332	DO jk = 1, jpkm1 ! Horizontal slab
	333	ze2u = rotdiff (ji,jj,jk)
	334	ze1v = hdivdiff(ji,jj,jk)
	335
	336	! horizontal diffusive trends
[6140]	337	zva = + ( ze2u - rotdiff (ji-1,jj,jk)) / ( e1v(ji,jj) * e3v_n(ji,jj,jk) ) &
[5656]	338	+ ( hdivdiff(ji,jj+1,jk) - ze1v ) / e2v(ji,jj)
	339
	340	! add it to the general momentum trends
	341	va(ji,jj,jk) = va(ji,jj,jk) + zva
	342	END DO
	343	ENDIF
[6140]	344	!
[5656]	345	END DO
	346	END DO
[6140]	347	!
[5656]	348	tabspongedone_v(i1+1:i2,j1+1:jmax) = .TRUE.
[6140]	349	!
[5656]	350	ENDIF
[6140]	351	!
[5656]	352	END SUBROUTINE interpun_sponge
	353
	354
	355	SUBROUTINE interpvn_sponge(tabres,i1,i2,j1,j2,k1,k2, before,nb,ndir)
[636]	356	!!---------------------------------------------
[5656]	357	!! * ROUTINE interpvn_sponge *
	358	!!---------------------------------------------
[636]	359	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2
	360	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
[5656]	361	LOGICAL, INTENT(in) :: before
	362	INTEGER, INTENT(in) :: nb , ndir
[6140]	363	!
	364	INTEGER :: ji, jj, jk
	365	REAL(wp) :: ze2u, ze1v, zua, zva, zbtr
[5656]	366	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: vbdiff
	367	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: rotdiff, hdivdiff
	368	INTEGER :: imax
[6140]	369	!!---------------------------------------------
[5656]	370
[6140]	371	IF( before ) THEN
[5656]	372	tabres = vn(i1:i2,j1:j2,:)
	373	ELSE
[6140]	374	!
[5656]	375	vbdiff(i1:i2,j1:j2,:) = (vb(i1:i2,j1:j2,:) - tabres(:,:,:))*vmask(i1:i2,j1:j2,:)
[6140]	376	!
[5656]	377	DO jk = 1, jpkm1 ! Horizontal slab
	378	! ! ===============
	379
	380	! ! --------
	381	! Horizontal divergence ! div
	382	! ! --------
	383	DO jj = j1+1,j2
	384	DO ji = i1,i2 ! vector opt.
[6140]	385	zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) * fsahm_spt(ji,jj)
	386	hdivdiff(ji,jj,jk) = ( e1v(ji,jj ) * e3v_n(ji,jj ,jk) * vbdiff(ji,jj ,jk) &
	387	& -e1v(ji,jj-1) * e3v_n(ji,jj-1,jk) * vbdiff(ji,jj-1,jk) ) * zbtr
[5656]	388	END DO
	389	END DO
	390	DO jj = j1,j2
	391	DO ji = i1,i2-1 ! vector opt.
[6140]	392	zbtr = r1_e1e2f(ji,jj) * e3f_n(ji,jj,jk) * fsahm_spf(ji,jj)
[5656]	393	rotdiff(ji,jj,jk) = ( e2v(ji+1,jj) * vbdiff(ji+1,jj,jk) &
[6140]	394	& -e2v(ji ,jj) * vbdiff(ji ,jj,jk) ) * fmask(ji,jj,jk) * zbtr
[5656]	395	END DO
	396	END DO
[6140]	397	END DO
[5656]	398
	399	! ! ===============
	400	!
	401
	402	imax = i2-1
[6140]	403	IF ((nbondi == 1).OR.(nbondi == 2)) imax = MIN(imax,nlci-3)
[5656]	404
	405	DO jj = j1+1, j2
	406	DO ji = i1+1, imax ! vector opt.
[6140]	407	IF( .NOT. tabspongedone_u(ji,jj) ) THEN
	408	DO jk = 1, jpkm1
	409	ua(ji,jj,jk) = ua(ji,jj,jk) &
	410	& - ( rotdiff (ji ,jj,jk) - rotdiff (ji,jj-1,jk)) / ( e2u(ji,jj) * e3u_n(ji,jj,jk) ) &
	411	& + ( hdivdiff(ji+1,jj,jk) - hdivdiff(ji,jj ,jk)) * r1_e1u(ji,jj)
[5656]	412	END DO
	413	ENDIF
	414	END DO
	415	END DO
[6140]	416	!
[5656]	417	tabspongedone_u(i1+1:imax,j1+1:j2) = .TRUE.
[6140]	418	!
[5656]	419	DO jj = j1+1, j2-1
	420	DO ji = i1+1, i2-1 ! vector opt.
[6140]	421	IF( .NOT. tabspongedone_v(ji,jj) ) THEN
	422	DO jk = 1, jpkm1
	423	va(ji,jj,jk) = va(ji,jj,jk) &
	424	& + ( rotdiff (ji,jj ,jk) - rotdiff (ji-1,jj,jk) ) / ( e1v(ji,jj) * e3v_n(ji,jj,jk) ) &
	425	& + ( hdivdiff(ji,jj+1,jk) - hdivdiff(ji ,jj,jk) ) * r1_e2v(ji,jj)
[5656]	426	END DO
	427	ENDIF
	428	END DO
	429	END DO
	430	tabspongedone_v(i1+1:i2-1,j1+1:j2-1) = .TRUE.
	431	ENDIF
[6140]	432	!
[5656]	433	END SUBROUTINE interpvn_sponge
	434
[390]	435	#else
[636]	436	CONTAINS
	437	SUBROUTINE agrif_opa_sponge_empty
	438	!!---------------------------------------------
	439	!! * ROUTINE agrif_OPA_sponge_empty *
	440	!!---------------------------------------------
	441	WRITE(,) 'agrif_opa_sponge : You should not have seen this print! error?'
	442	END SUBROUTINE agrif_opa_sponge_empty
[390]	443	#endif
[636]	444
[6140]	445	!!======================================================================
[636]	446	END MODULE agrif_opa_sponge

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

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

Context Navigation

source: branches/2017/dev_r7881_no_wrk_alloc/NEMOGCM/NEMO/NST_SRC/agrif_opa_sponge.F90 @ 7910

Download in other formats: