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.

agrif_opa_sponge.F90 in trunk/NEMOGCM/NEMO/NST_SRC – NEMO

Context Navigation

source: trunk/NEMOGCM/NEMO/NST_SRC/agrif_opa_sponge.F90 @ 3520

Last change on this file since 3520 was 3294, checked in by rblod, 12 years ago
Merge of 3.4beta into the trunk
Property svn:keywords set to `Id`
File size: 13.6 KB

Rev	Line
[390]	1	#define SPONGE
	2
[636]	3	Module agrif_opa_sponge
[2528]	4	#if defined key_agrif && ! defined key_offline
[636]	5	USE par_oce
	6	USE oce
	7	USE dom_oce
	8	USE in_out_manager
[782]	9	USE agrif_oce
[3294]	10	USE wrk_nemo
[390]	11
[636]	12	IMPLICIT NONE
	13	PRIVATE
[390]	14
[3294]	15	PUBLIC Agrif_Sponge_Tra, Agrif_Sponge_Dyn, interptsn, interpun, interpvn
[390]	16
[1156]	17	!!----------------------------------------------------------------------
[2528]	18	!! NEMO/NST 3.3 , NEMO Consortium (2010)
[1156]	19	!! $Id$
[2528]	20	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1156]	21	!!----------------------------------------------------------------------
[390]	22
[636]	23	CONTAINS
	24
[782]	25	SUBROUTINE Agrif_Sponge_Tra
[636]	26	!!---------------------------------------------
	27	!! * ROUTINE Agrif_Sponge_Tra *
	28	!!---------------------------------------------
	29	#include "domzgr_substitute.h90"
[2715]	30	!!
[3294]	31	INTEGER :: ji,jj,jk,jn
[390]	32	INTEGER :: spongearea
[636]	33	REAL(wp) :: timecoeff
[3294]	34	REAL(wp) :: ztsa, zabe1, zabe2, zbtr
	35	REAL(wp), POINTER, DIMENSION(:,: ) :: localviscsponge
	36	REAL(wp), POINTER, DIMENSION(:,: ) :: ztu, ztv
	37	REAL(wp), POINTER, DIMENSION(:,:,:,:) :: ztab
	38	REAL(wp), POINTER, DIMENSION(:,:,:,:) :: tsbdiff
[390]	39
	40	#if defined SPONGE
[3294]	41	CALL wrk_alloc( jpi, jpj, localviscsponge, ztu, ztv )
	42	CALL wrk_alloc( jpi, jpj, jpk, jpts, ztab, tsbdiff )
[390]	43
[636]	44	timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot()
	45
[390]	46	Agrif_SpecialValue=0.
	47	Agrif_UseSpecialValue = .TRUE.
[636]	48	ztab = 0.e0
[3294]	49	CALL Agrif_Bc_Variable(ztab, tsa_id,calledweight=timecoeff,procname=interptsn)
[390]	50	Agrif_UseSpecialValue = .FALSE.
	51
[3294]	52	tsbdiff(:,:,:,:) = tsb(:,:,:,:) - ztab(:,:,:,:)
[390]	53
	54	spongearea = 2 + 2 * Agrif_irhox()
	55
	56	localviscsponge = 0.
[782]	57
	58	IF (.NOT. spongedoneT) THEN
	59	spe1ur(:,:) = 0.
	60	spe2vr(:,:) = 0.
[390]	61
	62	IF ((nbondi == -1).OR.(nbondi == 2)) THEN
[636]	63	DO ji = 2, spongearea
	64	localviscsponge(ji,:) = visc_tra * (spongearea-ji)/real(spongearea-2)
	65	ENDDO
[782]	66
	67	spe1ur(2:spongearea-1,:)=0.5 * (localviscsponge(2:spongearea-1,:) + localviscsponge(3:spongearea,:)) &
	68	* e2u(2:spongearea-1,:) / e1u(2:spongearea-1,:)
	69
	70	spe2vr(2:spongearea,1:jpjm1) = 0.5 * (localviscsponge(2:spongearea,1:jpjm1) + &
	71	localviscsponge(2:spongearea,2:jpj)) &
	72	* e1v(2:spongearea,1:jpjm1) / e2v(2:spongearea,1:jpjm1)
[390]	73	ENDIF
	74
	75	IF ((nbondi == 1).OR.(nbondi == 2)) THEN
[636]	76	DO ji = nlci-spongearea + 1,nlci-1
	77	localviscsponge(ji,:) = visc_tra * (ji - (nlci-spongearea+1))/real(spongearea-2)
	78	ENDDO
[782]	79
	80	spe1ur(nlci-spongearea + 1:nlci-2,:)=0.5 * (localviscsponge(nlci-spongearea + 1:nlci-2,:) + &
	81	localviscsponge(nlci-spongearea + 2:nlci-1,:)) &
	82	* e2u(nlci-spongearea + 1:nlci-2,:) / e1u(nlci-spongearea + 1:nlci-2,:)
	83
	84	spe2vr(nlci-spongearea + 1:nlci-1,1:jpjm1) = 0.5 * (localviscsponge(nlci-spongearea + 1:nlci-1,1:jpjm1) &
	85	+ localviscsponge(nlci-spongearea + 1:nlci-1,2:jpj)) &
	86	* e1v(nlci-spongearea + 1:nlci-1,1:jpjm1) / e2v(nlci-spongearea + 1:nlci-1,1:jpjm1)
[390]	87	ENDIF
	88
	89
	90	IF ((nbondj == -1).OR.(nbondj == 2)) THEN
[636]	91	DO jj = 2, spongearea
	92	localviscsponge(:,jj) = visc_tra * (spongearea-jj)/real(spongearea-2)
	93	ENDDO
[782]	94
	95	spe1ur(1:jpim1,2:spongearea)=0.5 * (localviscsponge(1:jpim1,2:spongearea) + &
	96	localviscsponge(2:jpi,2:spongearea)) &
	97	* e2u(1:jpim1,2:spongearea) / e1u(1:jpim1,2:spongearea)
	98
	99	spe2vr(:,2:spongearea-1) = 0.5 * (localviscsponge(:,2:spongearea-1) + &
	100	localviscsponge(:,3:spongearea)) &
	101	* e1v(:,2:spongearea-1) / e2v(:,2:spongearea-1)
[390]	102	ENDIF
	103
	104	IF ((nbondj == 1).OR.(nbondj == 2)) THEN
[636]	105	DO jj = nlcj-spongearea + 1,nlcj-1
	106	localviscsponge(:,jj) = visc_tra * (jj - (nlcj-spongearea+1))/real(spongearea-2)
	107	ENDDO
[782]	108
	109	spe1ur(1:jpim1,nlcj-spongearea + 1:nlcj-1)=0.5 * (localviscsponge(1:jpim1,nlcj-spongearea + 1:nlcj-1) + &
	110	localviscsponge(2:jpi,nlcj-spongearea + 1:nlcj-1)) &
	111	* e2u(1:jpim1,nlcj-spongearea + 1:nlcj-1) / e1u(1:jpim1,nlcj-spongearea + 1:nlcj-1)
	112
	113	spe2vr(:,nlcj-spongearea + 1:nlcj-2) = 0.5 * (localviscsponge(:,nlcj-spongearea + 1:nlcj-2) + &
	114	localviscsponge(:,nlcj-spongearea + 2:nlcj-1)) &
	115	* e1v(:,nlcj-spongearea + 1:nlcj-2) / e2v(:,nlcj-spongearea + 1:nlcj-2)
[390]	116	ENDIF
[782]	117
	118	spbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:))
[390]	119
	120	spongedoneT = .TRUE.
	121	ENDIF
	122
[3294]	123	DO jn = 1, jpts
	124	DO jk = 1, jpkm1
	125	!
	126	DO jj = 1, jpjm1
	127	DO ji = 1, jpim1
	128	zabe1 = umask(ji,jj,jk) * spe1ur(ji,jj) * fse3u(ji,jj,jk)
	129	zabe2 = vmask(ji,jj,jk) * spe2vr(ji,jj) * fse3v(ji,jj,jk)
	130	ztu(ji,jj) = zabe1 * ( tsbdiff(ji+1,jj ,jk,jn) - tsbdiff(ji,jj,jk,jn) )
	131	ztv(ji,jj) = zabe2 * ( tsbdiff(ji ,jj+1,jk,jn) - tsbdiff(ji,jj,jk,jn) )
	132	ENDDO
[390]	133	ENDDO
	134
[3294]	135	DO jj = 2, jpjm1
	136	DO ji = 2, jpim1
	137	zbtr = spbtr2(ji,jj) / fse3t(ji,jj,jk)
	138	! horizontal diffusive trends
	139	ztsa = zbtr * ( ztu(ji,jj) - ztu(ji-1,jj ) &
	140	& + ztv(ji,jj) - ztv(ji ,jj-1) )
	141	! add it to the general tracer trends
	142	tsa(ji,jj,jk,jn) = tsa(ji,jj,jk,jn) + ztsa
	143	END DO
[390]	144	END DO
[3294]	145	!
	146	ENDDO
[636]	147	ENDDO
[390]	148
[3294]	149	CALL wrk_dealloc( jpi, jpj, localviscsponge, ztu, ztv )
	150	CALL wrk_dealloc( jpi, jpj, jpk, jpts, ztab, tsbdiff )
[390]	151	#endif
	152
[636]	153	END SUBROUTINE Agrif_Sponge_Tra
[390]	154
[782]	155	SUBROUTINE Agrif_Sponge_dyn
[636]	156	!!---------------------------------------------
	157	!! * ROUTINE Agrif_Sponge_dyn *
	158	!!---------------------------------------------
	159	#include "domzgr_substitute.h90"
[2715]	160	!!
[390]	161	INTEGER :: ji,jj,jk
	162	INTEGER :: spongearea
[636]	163	REAL(wp) :: timecoeff
[782]	164	REAL(wp) :: ze2u, ze1v, zua, zva, zbtr
[2715]	165	REAL(wp), POINTER, DIMENSION(:,:) :: localviscsponge
	166	REAL(wp), POINTER, DIMENSION(:,:,:) :: ubdiff, vbdiff
	167	REAL(wp), POINTER, DIMENSION(:,:,:) :: rotdiff, hdivdiff
	168	REAL(wp), POINTER, DIMENSION(:,:,:) :: ztab
[390]	169
	170	#if defined SPONGE
[3294]	171	CALL wrk_alloc( jpi, jpj, localviscsponge )
	172	CALL wrk_alloc( jpi, jpj, jpk, ztab, ubdiff, vbdiff, rotdiff, hdivdiff )
[390]	173
[636]	174	timecoeff = REAL(Agrif_NbStepint(),wp)/Agrif_rhot()
[390]	175
	176	Agrif_SpecialValue=0.
[782]	177	Agrif_UseSpecialValue = ln_spc_dyn
[636]	178	ztab = 0.e0
[2715]	179	CALL Agrif_Bc_Variable(ztab, ua_id,calledweight=timecoeff,procname=interpun)
[390]	180	Agrif_UseSpecialValue = .FALSE.
	181
[782]	182	ubdiff(:,:,:) = (ub(:,:,:) - ztab(:,:,:))*umask(:,:,:)
[390]	183
[636]	184	ztab = 0.e0
[390]	185	Agrif_SpecialValue=0.
[782]	186	Agrif_UseSpecialValue = ln_spc_dyn
[2715]	187	CALL Agrif_Bc_Variable(ztab, va_id,calledweight=timecoeff,procname=interpvn)
[390]	188	Agrif_UseSpecialValue = .FALSE.
	189
[782]	190	vbdiff(:,:,:) = (vb(:,:,:) - ztab(:,:,:))*vmask(:,:,:)
[390]	191
	192	spongearea = 2 + 2 * Agrif_irhox()
	193
	194	localviscsponge = 0.
	195
[782]	196	IF (.NOT. spongedoneU) THEN
	197	spe1ur2(:,:) = 0.
	198	spe2vr2(:,:) = 0.
	199
[390]	200	IF ((nbondi == -1).OR.(nbondi == 2)) THEN
[636]	201	DO ji = 2, spongearea
	202	localviscsponge(ji,:) = visc_dyn * (spongearea-ji)/real(spongearea-2)
	203	ENDDO
[782]	204
	205	spe1ur2(2:spongearea-1,:)=0.5 * (localviscsponge(2:spongearea-1,:) + localviscsponge(3:spongearea,:))
	206
	207	spe2vr2(2:spongearea,1:jpjm1) = 0.5 * (localviscsponge(2:spongearea,1:jpjm1) + &
	208	localviscsponge(2:spongearea,2:jpj))
[390]	209	ENDIF
	210
	211	IF ((nbondi == 1).OR.(nbondi == 2)) THEN
[636]	212	DO ji = nlci-spongearea + 1,nlci-1
	213	localviscsponge(ji,:) = visc_dyn * (ji - (nlci-spongearea+1))/real(spongearea-2)
	214	ENDDO
[782]	215
	216	spe1ur2(nlci-spongearea + 1:nlci-2,:)=0.5 * (localviscsponge(nlci-spongearea + 1:nlci-2,:) + &
	217	localviscsponge(nlci-spongearea + 2:nlci-1,:))
	218
	219	spe2vr2(nlci-spongearea + 1:nlci-1,1:jpjm1) = 0.5 * (localviscsponge(nlci-spongearea + 1:nlci-1,1:jpjm1) &
	220	+ localviscsponge(nlci-spongearea + 1:nlci-1,2:jpj))
[390]	221	ENDIF
	222
[782]	223
[390]	224	IF ((nbondj == -1).OR.(nbondj == 2)) THEN
[636]	225	DO jj = 2, spongearea
	226	localviscsponge(:,jj) = visc_dyn * (spongearea-jj)/real(spongearea-2)
	227	ENDDO
[782]	228
	229	spe1ur2(1:jpim1,2:spongearea)=0.5 * (localviscsponge(1:jpim1,2:spongearea) + &
	230	localviscsponge(2:jpi,2:spongearea))
	231
	232	spe2vr2(:,2:spongearea-1) = 0.5 * (localviscsponge(:,2:spongearea-1) + &
	233	localviscsponge(:,3:spongearea))
[390]	234	ENDIF
	235
	236	IF ((nbondj == 1).OR.(nbondj == 2)) THEN
[636]	237	DO jj = nlcj-spongearea + 1,nlcj-1
	238	localviscsponge(:,jj) = visc_dyn * (jj - (nlcj-spongearea+1))/real(spongearea-2)
	239	ENDDO
[782]	240
	241	spe1ur2(1:jpim1,nlcj-spongearea + 1:nlcj-1)=0.5 * (localviscsponge(1:jpim1,nlcj-spongearea + 1:nlcj-1) + &
	242	localviscsponge(2:jpi,nlcj-spongearea + 1:nlcj-1))
	243
	244	spe2vr2(:,nlcj-spongearea + 1:nlcj-2) = 0.5 * (localviscsponge(:,nlcj-spongearea + 1:nlcj-2) + &
	245	localviscsponge(:,nlcj-spongearea + 2:nlcj-1))
[636]	246	ENDIF
[390]	247
[782]	248	spongedoneU = .TRUE.
	249
[2715]	250	spbtr3(:,:) = 1./( e1f(:,:) * e2f(:,:))
[782]	251	ENDIF
	252
	253	IF (.NOT. spongedoneT) THEN
	254	spbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:))
	255	ENDIF
	256
	257	DO jk=1,jpkm1
	258	ubdiff(:,:,jk) = ubdiff(:,:,jk) * spe1ur2(:,:)
	259	vbdiff(:,:,jk) = vbdiff(:,:,jk) * spe2vr2(:,:)
	260	ENDDO
	261
[390]	262	hdivdiff = 0.
	263	rotdiff = 0.
	264
	265	DO jk = 1, jpkm1 ! Horizontal slab
	266	! ! ===============
	267
	268	! ! --------
	269	! Horizontal divergence ! div
	270	! ! --------
	271	DO jj = 2, jpjm1
	272	DO ji = 2, jpim1 ! vector opt.
[782]	273	zbtr = spbtr2(ji,jj) / fse3t(ji,jj,jk)
[390]	274	hdivdiff(ji,jj,jk) = &
	275	( e2u(ji,jj)fse3u(ji,jj,jk) &
	276	ubdiff(ji,jj,jk) - e2u(ji-1,jj )* &
	277	fse3u(ji-1,jj ,jk) * ubdiff(ji-1,jj ,jk) &
[636]	278	+ e1v(ji,jj)fse3v(ji,jj,jk) &
[390]	279	vbdiff(ji,jj,jk) - e1v(ji ,jj-1)* &
[782]	280	fse3v(ji ,jj-1,jk) * vbdiff(ji ,jj-1,jk) ) * zbtr
[390]	281	END DO
	282	END DO
[636]	283
[390]	284	DO jj = 1, jpjm1
	285	DO ji = 1, jpim1 ! vector opt.
[782]	286	zbtr = spbtr3(ji,jj) * fse3f(ji,jj,jk)
	287	rotdiff(ji,jj,jk) = ( e2v(ji+1,jj ) * vbdiff(ji+1,jj ,jk) - e2v(ji,jj) * vbdiff(ji,jj,jk) &
	288	& - e1u(ji ,jj+1) * ubdiff(ji ,jj+1,jk) + e1u(ji,jj) * ubdiff(ji,jj,jk) ) &
	289	& * fmask(ji,jj,jk) * zbtr
[390]	290	END DO
	291	END DO
[636]	292
	293	ENDDO
	294
[390]	295	! ! ===============
	296	DO jk = 1, jpkm1 ! Horizontal slab
	297	! ! ===============
	298	DO jj = 2, jpjm1
	299	DO ji = 2, jpim1 ! vector opt.
[469]	300	ze2u = rotdiff (ji,jj,jk)
	301	ze1v = hdivdiff(ji,jj,jk)
[390]	302	! horizontal diffusive trends
[782]	303	zua = - ( ze2u - rotdiff (ji,jj-1,jk)) / ( e2u(ji,jj) * fse3u(ji,jj,jk) ) &
[636]	304	+ ( hdivdiff(ji+1,jj,jk) - ze1v &
	305	) / e1u(ji,jj)
[390]	306
[782]	307	zva = + ( ze2u - rotdiff (ji-1,jj,jk)) / ( e1v(ji,jj) * fse3v(ji,jj,jk) ) &
[636]	308	+ ( hdivdiff(ji,jj+1,jk) - ze1v &
	309	) / e2v(ji,jj)
[390]	310
	311	! add it to the general momentum trends
	312	ua(ji,jj,jk) = ua(ji,jj,jk) + zua
	313	va(ji,jj,jk) = va(ji,jj,jk) + zva
	314	END DO
	315	END DO
	316	! ! ===============
	317	END DO ! End of slab
	318	! ! ===============
[3294]	319	CALL wrk_dealloc( jpi, jpj, localviscsponge )
	320	CALL wrk_dealloc( jpi, jpj, jpk, ztab, ubdiff, vbdiff, rotdiff, hdivdiff )
[390]	321
	322	#endif
	323
[636]	324	END SUBROUTINE Agrif_Sponge_dyn
[390]	325
[3294]	326	SUBROUTINE interptsn(tabres,i1,i2,j1,j2,k1,k2,n1,n2)
[636]	327	!!---------------------------------------------
[3294]	328	!! * ROUTINE interptsn *
[636]	329	!!---------------------------------------------
[390]	330	# include "domzgr_substitute.h90"
[636]	331
[3294]	332	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
	333	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
[390]	334
[3294]	335	tabres(i1:i2,j1:j2,k1:k2,n1:n2) = tsn(i1:i2,j1:j2,k1:k2,n1:n2)
[390]	336
[3294]	337	END SUBROUTINE interptsn
[636]	338
	339	SUBROUTINE interpun(tabres,i1,i2,j1,j2,k1,k2)
	340	!!---------------------------------------------
	341	!! * ROUTINE interpun *
	342	!!---------------------------------------------
[390]	343	# include "domzgr_substitute.h90"
[636]	344
	345	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2
	346	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
[390]	347
[636]	348	tabres(i1:i2,j1:j2,k1:k2) = un(i1:i2,j1:j2,k1:k2)
[390]	349
[636]	350	END SUBROUTINE interpun
	351
	352	SUBROUTINE interpvn(tabres,i1,i2,j1,j2,k1,k2)
	353	!!---------------------------------------------
	354	!! * ROUTINE interpvn *
	355	!!---------------------------------------------
[390]	356	# include "domzgr_substitute.h90"
[636]	357
	358	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2
	359	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres
[390]	360
[636]	361	tabres(i1:i2,j1:j2,k1:k2) = vn(i1:i2,j1:j2,k1:k2)
[390]	362
[636]	363	END SUBROUTINE interpvn
[390]	364
	365	#else
[636]	366	CONTAINS
	367
	368	SUBROUTINE agrif_opa_sponge_empty
	369	!!---------------------------------------------
	370	!! * ROUTINE agrif_OPA_sponge_empty *
	371	!!---------------------------------------------
	372	WRITE(,) 'agrif_opa_sponge : You should not have seen this print! error?'
	373	END SUBROUTINE agrif_opa_sponge_empty
[390]	374	#endif
[636]	375
	376	END MODULE agrif_opa_sponge

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

Download in other formats: