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_update.F90 in branches/2015/dev_r5803_UKMO_AGRIF_Vert_interp/NEMOGCM/NEMO/NST_SRC – NEMO

Context Navigation

source: branches/2015/dev_r5803_UKMO_AGRIF_Vert_interp/NEMOGCM/NEMO/NST_SRC/agrif_opa_update.F90 @ 6401

Last change on this file since 6401 was 6401, checked in by timgraham, 8 years ago
Reinstate svn keywords before update to head of trunk
Property svn:keywords set to `Id`
File size: 29.8 KB

Rev	Line
[6258]	1	#undef TWO_WAY /* TWO WAY NESTING */
[5656]	2	#undef DECAL_FEEDBACK /* SEPARATION of INTERFACES*/
	3
[636]	4	MODULE agrif_opa_update
[2528]	5	#if defined key_agrif && ! defined key_offline
[636]	6	USE par_oce
	7	USE oce
	8	USE dom_oce
[782]	9	USE agrif_oce
[2715]	10	USE in_out_manager ! I/O manager
	11	USE lib_mpp
[3294]	12	USE wrk_nemo
[4486]	13	USE dynspg_oce
[5656]	14	USE zdf_oce ! vertical physics: ocean variables
[3294]	15
[636]	16	IMPLICIT NONE
	17	PRIVATE
[390]	18
[5656]	19	PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn,Update_Scales
[6258]	20
	21	! VERTICAL REFINEMENT BEGIN
	22	PUBLIC Agrif_Init_UpdateScales
	23	REAL, DIMENSION(:,:,:), ALLOCATABLE :: update_scales_t, update_scales_u, update_scales_v
	24	!$AGRIF_DO_NOT_TREAT
	25	LOGICAL :: scaleT, scaleU, scaleV = .FALSE.
	26	!$AGRIF_END_DO_NOT_TREAT
	27	! VERTICAL REFINEMENT END
	28
[5656]	29	# if defined key_zdftke
	30	PUBLIC Agrif_Update_Tke
	31	# endif
[1156]	32	!!----------------------------------------------------------------------
[5656]	33	!! NEMO/NST 3.6 , NEMO Consortium (2010)
[6401]	34	!! $Id$
[2528]	35	!! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
[1156]	36	!!----------------------------------------------------------------------
	37
[6258]	38	# include "domzgr_substitute.h90"
	39
[636]	40	CONTAINS
	41
[6258]	42	! VERTICAL REFINEMENT BEGIN
	43	SUBROUTINE Agrif_Init_UpdateScales
	44
	45	scaleT = .TRUE.
	46	Agrif_UseSpecialValueInUpdate = .FALSE.
	47	Agrif_SpecialValueFineGrid = 0.
	48	Call Agrif_Update_Variable(scales_t_id,procname=updatescales)
	49	Agrif_UseSpecialValueInUpdate = .FALSE.
	50	scaleT = .FALSE.
	51
	52	scaleU = .TRUE.
	53	Call Agrif_Update_Variable(scales_u_id,procname=updatescales)
	54	scaleU = .FALSE.
	55
	56	scaleV = .TRUE.
	57	Call Agrif_Update_Variable(scales_v_id,procname=updatescales)
	58	scaleV = .FALSE.
	59
	60	END SUBROUTINE Agrif_Init_UpdateScales
	61
	62	SUBROUTINE updatescales( ptab, i1, i2, j1, j2, k1, k2, before )
	63
	64	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2
	65	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
	66	LOGICAL, iNTENT(in) :: before
	67
	68	INTEGER :: ji,jj,jk
	69
	70	IF (before) THEN
	71	IF (scaleT) THEN
	72	DO jk=k1,k2
	73	DO jj=j1,j2
	74	DO ji=i1,i2
	75	ptab(ji,jj,jk) = fse3t(ji,jj,jk)*tmask(ji,jj,jk)
	76	END DO
	77	END DO
	78	END DO
	79	ELSEIF (scaleU) THEN
	80	DO jk=k1,k2
	81	DO jj=j1,j2
	82	DO ji=i1,i2
	83	ptab(ji,jj,jk) = fse3u(ji,jj,jk)*umask(ji,jj,jk)
	84	END DO
	85	END DO
	86	END DO
	87	ELSEIF (scaleV) THEN
	88	DO jk=k1,k2
	89	DO jj=j1,j2
	90	DO ji=i1,i2
	91	ptab(ji,jj,jk) = fse3v(ji,jj,jk)*vmask(ji,jj,jk)
	92	END DO
	93	END DO
	94	END DO
	95	ENDIF
	96	ELSE
	97	IF (scaleT) THEN
	98	Allocate(update_scales_t(i1:i2,j1:j2,k1:k2))
	99	update_scales_t(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2)
	100	ELSEIF (scaleU) THEN
	101	Allocate(update_scales_u(i1:i2,j1:j2,k1:k2))
	102	update_scales_u(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2)
	103	ELSEIF (scaleV) THEN
	104	Allocate(update_scales_v(i1:i2,j1:j2,k1:k2))
	105	update_scales_v(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2)
	106	ENDIF
	107	ENDIF
	108
	109	END SUBROUTINE updatescales
	110	! VERTICAL REFINEMENT END
	111
[5656]	112	RECURSIVE SUBROUTINE Agrif_Update_Tra( )
[636]	113	!!---------------------------------------------
	114	!! * ROUTINE Agrif_Update_Tra *
	115	!!---------------------------------------------
[5656]	116	!
	117	IF (Agrif_Root()) RETURN
	118	!
	119	#if defined TWO_WAY
	120	IF (lwp.AND.lk_agrif_debug) Write(,) 'Update tracers from grid Number',Agrif_Fixed(), 'nbcline', nbcline
[636]	121
[390]	122	Agrif_UseSpecialValueInUpdate = .TRUE.
	123	Agrif_SpecialValueFineGrid = 0.
[5656]	124	!
[636]	125	IF (MOD(nbcline,nbclineupdate) == 0) THEN
[5656]	126	# if ! defined DECAL_FEEDBACK
	127	CALL Agrif_Update_Variable(tsn_id, procname=updateTS)
	128	# else
	129	CALL Agrif_Update_Variable(tsn_id, locupdate=(/1,0/),procname=updateTS)
	130	# endif
[390]	131	ELSE
[5656]	132	# if ! defined DECAL_FEEDBACK
	133	CALL Agrif_Update_Variable(tsn_id,locupdate=(/0,2/), procname=updateTS)
	134	# else
	135	CALL Agrif_Update_Variable(tsn_id,locupdate=(/1,2/), procname=updateTS)
	136	# endif
[390]	137	ENDIF
[5656]	138	!
[390]	139	Agrif_UseSpecialValueInUpdate = .FALSE.
[5656]	140	!
	141	IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update:
	142	CALL Agrif_ChildGrid_To_ParentGrid()
	143	CALL Agrif_Update_Tra()
	144	CALL Agrif_ParentGrid_To_ChildGrid()
	145	ENDIF
	146	!
[390]	147	#endif
[5656]	148	!
[636]	149	END SUBROUTINE Agrif_Update_Tra
[390]	150
[5656]	151	RECURSIVE SUBROUTINE Agrif_Update_Dyn( )
[636]	152	!!---------------------------------------------
	153	!! * ROUTINE Agrif_Update_Dyn *
	154	!!---------------------------------------------
[5656]	155	!
	156	IF (Agrif_Root()) RETURN
	157	!
[390]	158	#if defined TWO_WAY
[5656]	159	IF (lwp.AND.lk_agrif_debug) Write(,) 'Update momentum from grid Number',Agrif_Fixed(), 'nbcline', nbcline
[390]	160
[5656]	161	Agrif_UseSpecialValueInUpdate = .FALSE.
	162	Agrif_SpecialValueFineGrid = 0.
	163	!
[390]	164	IF (mod(nbcline,nbclineupdate) == 0) THEN
[5656]	165	# if ! defined DECAL_FEEDBACK
	166	CALL Agrif_Update_Variable(un_update_id,procname = updateU)
	167	CALL Agrif_Update_Variable(vn_update_id,procname = updateV)
	168	# else
	169	CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU)
	170	CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV)
	171	# endif
[390]	172	ELSE
[5656]	173	# if ! defined DECAL_FEEDBACK
	174	CALL Agrif_Update_Variable(un_update_id,locupdate=(/0,1/),procname = updateU)
	175	CALL Agrif_Update_Variable(vn_update_id,locupdate=(/0,1/),procname = updateV)
	176	# else
	177	CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateU)
	178	CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updateV)
	179	# endif
[390]	180	ENDIF
	181
[5656]	182	# if ! defined DECAL_FEEDBACK
	183	CALL Agrif_Update_Variable(e1u_id,procname = updateU2d)
	184	CALL Agrif_Update_Variable(e2v_id,procname = updateV2d)
	185	# else
	186	CALL Agrif_Update_Variable(e1u_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU2d)
	187	CALL Agrif_Update_Variable(e2v_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV2d)
	188	# endif
[636]	189
[5656]	190	# if defined key_dynspg_ts
[4486]	191	IF (ln_bt_fw) THEN
	192	! Update time integrated transports
	193	IF (mod(nbcline,nbclineupdate) == 0) THEN
[5656]	194	# if ! defined DECAL_FEEDBACK
	195	CALL Agrif_Update_Variable(ub2b_update_id,procname = updateub2b)
	196	CALL Agrif_Update_Variable(vb2b_update_id,procname = updatevb2b)
	197	# else
	198	CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateub2b)
	199	CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updatevb2b)
	200	# endif
[4486]	201	ELSE
[5656]	202	# if ! defined DECAL_FEEDBACK
	203	CALL Agrif_Update_Variable(ub2b_update_id,locupdate=(/0,1/),procname = updateub2b)
	204	CALL Agrif_Update_Variable(vb2b_update_id,locupdate=(/0,1/),procname = updatevb2b)
	205	# else
	206	CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateub2b)
	207	CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updatevb2b)
	208	# endif
[4486]	209	ENDIF
[5656]	210	END IF
	211	# endif
	212	!
[390]	213	nbcline = nbcline + 1
[5656]	214	!
	215	Agrif_UseSpecialValueInUpdate = .TRUE.
[636]	216	Agrif_SpecialValueFineGrid = 0.
[5656]	217	# if ! defined DECAL_FEEDBACK
	218	CALL Agrif_Update_Variable(sshn_id,procname = updateSSH)
	219	# else
	220	CALL Agrif_Update_Variable(sshn_id,locupdate=(/1,0/),procname = updateSSH)
	221	# endif
[636]	222	Agrif_UseSpecialValueInUpdate = .FALSE.
[5656]	223	!
[390]	224	#endif
[5656]	225	!
	226	! Do recursive update:
	227	IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update:
	228	CALL Agrif_ChildGrid_To_ParentGrid()
	229	CALL Agrif_Update_Dyn()
	230	CALL Agrif_ParentGrid_To_ChildGrid()
	231	ENDIF
	232	!
[636]	233	END SUBROUTINE Agrif_Update_Dyn
	234
[5656]	235	# if defined key_zdftke
	236	SUBROUTINE Agrif_Update_Tke( kt )
[636]	237	!!---------------------------------------------
[5656]	238	!! * ROUTINE Agrif_Update_Tke *
[636]	239	!!---------------------------------------------
[5656]	240	!!
[636]	241	INTEGER, INTENT(in) :: kt
[5656]	242	!
[6258]	243	return
	244
[5656]	245	IF( ( Agrif_NbStepint() .NE. 0 ) .AND. (kt /= 0) ) RETURN
	246	# if defined TWO_WAY
[636]	247
[5656]	248	Agrif_UseSpecialValueInUpdate = .TRUE.
	249	Agrif_SpecialValueFineGrid = 0.
[636]	250
[5656]	251	CALL Agrif_Update_Variable( en_id, locupdate=(/0,0/), procname=updateEN )
	252	CALL Agrif_Update_Variable(avt_id, locupdate=(/0,0/), procname=updateAVT )
	253	CALL Agrif_Update_Variable(avm_id, locupdate=(/0,0/), procname=updateAVM )
	254
	255	Agrif_UseSpecialValueInUpdate = .FALSE.
	256
	257	# endif
	258
	259	END SUBROUTINE Agrif_Update_Tke
	260	# endif /* key_zdftke */
	261
[6258]	262	SUBROUTINE updateTS( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before )
[636]	263	!!---------------------------------------------
	264	!! * ROUTINE updateT *
	265	!!---------------------------------------------
[6258]	266
[3294]	267	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
[6258]	268	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: ptab
[5656]	269	LOGICAL, INTENT(in) :: before
	270	!!
[3294]	271	INTEGER :: ji,jj,jk,jn
[6258]	272	! VERTICAL REFINEMENT BEGIN
	273	REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk,n1:n2) :: ptab_child
	274	REAL(wp) :: h_in(k1:k2)
	275	REAL(wp) :: h_out(1:jpk)
	276	INTEGER :: N_in, N_out
	277	REAL(wp) :: h_diff
	278	REAL(wp) :: tabin(k1:k2,n1:n2)
	279	! VERTICAL REFINEMENT END
[5656]	280	!!---------------------------------------------
	281	!
[636]	282	IF (before) THEN
[3294]	283	DO jn = n1,n2
	284	DO jk=k1,k2
	285	DO jj=j1,j2
	286	DO ji=i1,i2
[6258]	287	ptab(ji,jj,jk,jn) = tsn(ji,jj,jk,jn)
[3294]	288	END DO
[636]	289	END DO
	290	END DO
	291	END DO
	292	ELSE
[6258]	293	! VERTICAL REFINEMENT BEGIN
	294	ptab_child(:,:,:,:) = 0.
	295
	296	DO jj=j1,j2
	297	DO ji=i1,i2
	298	N_in = 0
	299	DO jk=k1,k2
	300	IF (update_scales_t(ji,jj,jk) == 0) EXIT
	301	N_in = N_in + 1
	302	tabin(jk,:) = ptab(ji,jj,jk,:)
	303	h_in(N_in) = update_scales_t(ji,jj,jk)
	304	ENDDO
	305	N_out = 0
	306	DO jk=1,jpk
	307	IF (tmask(ji,jj,jk) == 0) EXIT
	308	N_out = N_out + 1
	309	h_out(N_out) = fse3t(ji,jj,jk)
	310	ENDDO
	311	IF (N_in > 0) THEN
	312	h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
	313	IF (h_diff > 0) THEN
	314	N_in = N_in + 1
	315	h_in(N_in) = h_diff
	316	tabin(N_in,:) = tabin(N_in-1,:)
	317	ELSEIF (h_diff < 0) THEN
	318	print *,'CHECK YOUR bathy T points ...',ji,jj,h_diff,sum(h_in(1:N_in)),sum(h_out(1:N_out))
	319	print *,'Nval = ',N_out,mbathy(ji,jj)
	320	print *,'BATHY = ',gdepw_0(ji,jj,mbathy(ji,jj)+1),sum(e3t_0(ji,jj,1:mbathy(ji,jj)))
	321	! STOP
	322	N_out = N_out + 1
	323	h_out(N_out) = - h_diff
	324	ENDIF
	325	DO jn=n1,n2
	326	CALL reconstructandremap(tabin(1:N_in,jn),h_in(1:N_in),ptab_child(ji,jj,1:N_out,jn),h_out(1:N_out),N_in,N_out)
	327	ENDDO
	328	ENDIF
	329	ENDDO
	330	ENDDO
	331
	332	! WARNING :
	333	! ptab replaced by ptab_child in the following
	334	! k1 k2 replaced by 1 jpk
	335	! VERTICAL REFINEMENT END
	336
[4491]	337	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
[5656]	338	! Add asselin part
[4491]	339	DO jn = n1,n2
[6258]	340	DO jk=1,jpk
[4491]	341	DO jj=j1,j2
	342	DO ji=i1,i2
[6258]	343	IF( ptab_child(ji,jj,jk,jn) .NE. 0. ) THEN
[4491]	344	tsb(ji,jj,jk,jn) = tsb(ji,jj,jk,jn) &
[6258]	345	& + atfp * ( ptab_child(ji,jj,jk,jn) &
[5656]	346	& - tsn(ji,jj,jk,jn) ) * tmask(ji,jj,jk)
[4491]	347	ENDIF
	348	ENDDO
	349	ENDDO
	350	ENDDO
	351	ENDDO
	352	ENDIF
[3294]	353	DO jn = n1,n2
[6258]	354	DO jk=1,jpk
[3294]	355	DO jj=j1,j2
	356	DO ji=i1,i2
[6258]	357	IF( ptab_child(ji,jj,jk,jn) .NE. 0. ) THEN
	358	tsn(ji,jj,jk,jn) = ptab_child(ji,jj,jk,jn) * tmask(ji,jj,jk)
[3294]	359	END IF
	360	END DO
[636]	361	END DO
	362	END DO
	363	END DO
	364	ENDIF
[5656]	365	!
[3294]	366	END SUBROUTINE updateTS
[390]	367
[6258]	368	SUBROUTINE updateu( ptab, i1, i2, j1, j2, k1, k2, before )
[636]	369	!!---------------------------------------------
	370	!! * ROUTINE updateu *
	371	!!---------------------------------------------
[5656]	372	!!
[636]	373	INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
[6258]	374	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
[636]	375	LOGICAL, INTENT(in) :: before
[5656]	376	!!
[636]	377	INTEGER :: ji, jj, jk
	378	REAL(wp) :: zrhoy
[6258]	379	! VERTICAL REFINEMENT BEGIN
	380	REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: ptab_child
	381	REAL(wp) :: h_in(k1:k2)
	382	REAL(wp) :: h_out(1:jpk)
	383	INTEGER :: N_in, N_out
	384	REAL(wp) :: h_diff
	385	REAL(wp) :: tabin(k1:k2)
	386	! VERTICAL REFINEMENT END
[5656]	387	!!---------------------------------------------
	388	!
[636]	389	IF (before) THEN
	390	zrhoy = Agrif_Rhoy()
[390]	391	DO jk=k1,k2
[636]	392	DO jj=j1,j2
	393	DO ji=i1,i2
[6258]	394	ptab(ji,jj,jk) = e2u(ji,jj) * un(ji,jj,jk)
	395	ptab(ji,jj,jk) = ptab(ji,jj,jk) * fse3u_n(ji,jj,jk)
[636]	396	END DO
	397	END DO
	398	END DO
[6258]	399	ptab = zrhoy * ptab
[636]	400	ELSE
[6258]	401	! VERTICAL REFINEMENT BEGIN
	402	ptab_child(:,:,:) = 0.
	403
	404	DO jj=j1,j2
	405	DO ji=i1,i2
	406	N_in = 0
	407	DO jk=k1,k2
	408	IF (update_scales_u(ji,jj,jk) == 0) EXIT
	409	N_in = N_in + 1
	410	tabin(jk) = ptab(ji,jj,jk)/update_scales_u(ji,jj,jk)
	411	h_in(N_in) = update_scales_u(ji,jj,jk)
	412	ENDDO
	413	N_out = 0
	414	DO jk=1,jpk
	415	IF (umask(ji,jj,jk) == 0) EXIT
	416	N_out = N_out + 1
	417	h_out(N_out) = fse3u(ji,jj,jk)
	418	ENDDO
	419	IF (N_in * N_out > 0) THEN
	420	h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
	421	if (h_diff < 0.) then
	422	print *,'CHECK YOUR BATHY ...'
	423	stop
	424	else ! Extends with 0
	425	N_in = N_in + 1
	426	tabin(N_in) = 0.
	427	h_in(N_in) = h_diff
	428	endif
	429	CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),ptab_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out)
	430	ENDIF
	431	ENDDO
	432	ENDDO
	433
	434	! WARNING :
	435	! ptab replaced by ptab_child in the following
	436	! k1 k2 replaced by 1 jpk
	437	! remove division by fs e3u (already done)
	438	! VERTICAL REFINEMENT END
	439
	440	DO jk=1,jpk
[636]	441	DO jj=j1,j2
	442	DO ji=i1,i2
[6258]	443	ptab_child(ji,jj,jk) = ptab_child(ji,jj,jk) / e2u(ji,jj)
[4491]	444	!
	445	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
	446	ub(ji,jj,jk) = ub(ji,jj,jk) &
[6258]	447	& + atfp * ( ptab_child(ji,jj,jk) - un(ji,jj,jk) ) * umask(ji,jj,jk)
[4491]	448	ENDIF
	449	!
[6258]	450	un(ji,jj,jk) = ptab_child(ji,jj,jk) * umask(ji,jj,jk)
[636]	451	END DO
	452	END DO
	453	END DO
	454	ENDIF
[5656]	455	!
[636]	456	END SUBROUTINE updateu
[390]	457
[6258]	458	SUBROUTINE updatev( ptab, i1, i2, j1, j2, k1, k2, before )
[636]	459	!!---------------------------------------------
	460	!! * ROUTINE updatev *
	461	!!---------------------------------------------
[5656]	462	!!
[636]	463	INTEGER :: i1,i2,j1,j2,k1,k2
	464	INTEGER :: ji,jj,jk
[6258]	465	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: ptab
[636]	466	LOGICAL :: before
[5656]	467	!!
[636]	468	REAL(wp) :: zrhox
[6258]	469	! VERTICAL REFINEMENT BEGIN
	470	REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: ptab_child
	471	REAL(wp) :: h_in(k1:k2)
	472	REAL(wp) :: h_out(1:jpk)
	473	INTEGER :: N_in, N_out
	474	REAL(wp) :: h_diff
	475	REAL(wp) :: tabin(k1:k2)
	476	! VERTICAL REFINEMENT END
[5656]	477	!!---------------------------------------------
	478	!
[636]	479	IF (before) THEN
	480	zrhox = Agrif_Rhox()
[390]	481	DO jk=k1,k2
[636]	482	DO jj=j1,j2
	483	DO ji=i1,i2
[6258]	484	ptab(ji,jj,jk) = e1v(ji,jj) * vn(ji,jj,jk)
	485	ptab(ji,jj,jk) = ptab(ji,jj,jk) * fse3v_n(ji,jj,jk)
[636]	486	END DO
	487	END DO
	488	END DO
[6258]	489	ptab = zrhox * ptab
[636]	490	ELSE
[6258]	491	! VERTICAL REFINEMENT BEGIN
	492	ptab_child(:,:,:) = 0.
	493
	494	DO jj=j1,j2
	495	DO ji=i1,i2
	496	N_in = 0
	497	DO jk=k1,k2
	498	IF (update_scales_v(ji,jj,jk) == 0) EXIT
	499	N_in = N_in + 1
	500	tabin(jk) = ptab(ji,jj,jk)/update_scales_v(ji,jj,jk)
	501	h_in(N_in) = update_scales_v(ji,jj,jk)
	502	ENDDO
	503	N_out = 0
	504	DO jk=1,jpk
	505	IF (vmask(ji,jj,jk) == 0) EXIT
	506	N_out = N_out + 1
	507	h_out(N_out) = fse3v(ji,jj,jk)
	508	ENDDO
	509	IF (N_in * N_out > 0) THEN
	510	h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))
	511	if (h_diff < 0.) then
	512	print *,'CHECK YOUR BATHY ...'
	513	stop
	514	else ! Extends with 0
	515	N_in = N_in + 1
	516	tabin(N_in) = 0.
	517	h_in(N_in) = h_diff
	518	endif
	519	CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),ptab_child(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out)
	520	ENDIF
	521	ENDDO
	522	ENDDO
	523
	524	! WARNING :
	525	! ptab replaced by ptab_child in the following
	526	! k1 k2 replaced by 1 jpk
	527	! remove division by fs e3v (already done)
	528	! VERTICAL REFINEMENT END
	529
[390]	530	DO jk=k1,k2
[636]	531	DO jj=j1,j2
	532	DO ji=i1,i2
[6258]	533	ptab_child(ji,jj,jk) = ptab_child(ji,jj,jk) / e1v(ji,jj)
[4491]	534	!
	535	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
	536	vb(ji,jj,jk) = vb(ji,jj,jk) &
[6258]	537	& + atfp * ( ptab_child(ji,jj,jk) - vn(ji,jj,jk) ) * vmask(ji,jj,jk)
[4491]	538	ENDIF
	539	!
[6258]	540	vn(ji,jj,jk) = ptab_child(ji,jj,jk) * vmask(ji,jj,jk)
[636]	541	END DO
	542	END DO
	543	END DO
	544	ENDIF
[5656]	545	!
[636]	546	END SUBROUTINE updatev
[390]	547
[636]	548	SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before )
	549	!!---------------------------------------------
	550	!! * ROUTINE updateu2d *
	551	!!---------------------------------------------
[390]	552	# include "domzgr_substitute.h90"
[5656]	553	!!
[636]	554	INTEGER, INTENT(in) :: i1, i2, j1, j2
	555	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	556	LOGICAL, INTENT(in) :: before
[5656]	557	!!
[636]	558	INTEGER :: ji, jj, jk
	559	REAL(wp) :: zrhoy
[4491]	560	REAL(wp) :: zcorr
[5656]	561	!!---------------------------------------------
	562	!
[636]	563	IF (before) THEN
	564	zrhoy = Agrif_Rhoy()
	565	DO jj=j1,j2
	566	DO ji=i1,i2
[4486]	567	tabres(ji,jj) = un_b(ji,jj) * hu(ji,jj) * e2u(ji,jj)
[636]	568	END DO
	569	END DO
	570	tabres = zrhoy * tabres
	571	ELSE
	572	DO jj=j1,j2
	573	DO ji=i1,i2
[4491]	574	tabres(ji,jj) = tabres(ji,jj) * hur(ji,jj) / e2u(ji,jj)
	575	!
	576	! Update "now" 3d velocities:
	577	spgu(ji,jj) = 0.e0
	578	DO jk=1,jpkm1
	579	spgu(ji,jj) = spgu(ji,jj) + fse3u_n(ji,jj,jk) * un(ji,jj,jk)
	580	END DO
	581	spgu(ji,jj) = spgu(ji,jj) * hur(ji,jj)
	582	!
	583	zcorr = tabres(ji,jj) - spgu(ji,jj)
	584	DO jk=1,jpkm1
	585	un(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk)
	586	END DO
	587	!
[4486]	588	! Update barotropic velocities:
[4491]	589	#if defined key_dynspg_ts
	590	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
	591	zcorr = tabres(ji,jj) - un_b(ji,jj)
	592	ub_b(ji,jj) = ub_b(ji,jj) + atfp * zcorr * umask(ji,jj,1)
	593	END IF
	594	#endif
	595	un_b(ji,jj) = tabres(ji,jj) * umask(ji,jj,1)
	596	!
	597	! Correct "before" velocities to hold correct bt component:
	598	spgu(ji,jj) = 0.e0
	599	DO jk=1,jpkm1
	600	spgu(ji,jj) = spgu(ji,jj) + fse3u_b(ji,jj,jk) * ub(ji,jj,jk)
	601	END DO
	602	spgu(ji,jj) = spgu(ji,jj) * hur_b(ji,jj)
	603	!
	604	zcorr = ub_b(ji,jj) - spgu(ji,jj)
	605	DO jk=1,jpkm1
	606	ub(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk)
	607	END DO
	608	!
[636]	609	END DO
	610	END DO
	611	ENDIF
[5656]	612	!
[636]	613	END SUBROUTINE updateu2d
[390]	614
[636]	615	SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before )
	616	!!---------------------------------------------
	617	!! * ROUTINE updatev2d *
	618	!!---------------------------------------------
	619	INTEGER, INTENT(in) :: i1, i2, j1, j2
	620	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	621	LOGICAL, INTENT(in) :: before
[5656]	622	!!
[636]	623	INTEGER :: ji, jj, jk
	624	REAL(wp) :: zrhox
[4491]	625	REAL(wp) :: zcorr
[5656]	626	!!---------------------------------------------
	627	!
[636]	628	IF (before) THEN
	629	zrhox = Agrif_Rhox()
	630	DO jj=j1,j2
	631	DO ji=i1,i2
[4486]	632	tabres(ji,jj) = vn_b(ji,jj) * hv(ji,jj) * e1v(ji,jj)
[636]	633	END DO
	634	END DO
	635	tabres = zrhox * tabres
	636	ELSE
	637	DO jj=j1,j2
	638	DO ji=i1,i2
[4491]	639	tabres(ji,jj) = tabres(ji,jj) * hvr(ji,jj) / e1v(ji,jj)
	640	!
	641	! Update "now" 3d velocities:
	642	spgv(ji,jj) = 0.e0
	643	DO jk=1,jpkm1
	644	spgv(ji,jj) = spgv(ji,jj) + fse3v_n(ji,jj,jk) * vn(ji,jj,jk)
	645	END DO
	646	spgv(ji,jj) = spgv(ji,jj) * hvr(ji,jj)
	647	!
	648	zcorr = tabres(ji,jj) - spgv(ji,jj)
	649	DO jk=1,jpkm1
	650	vn(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk)
	651	END DO
	652	!
[4486]	653	! Update barotropic velocities:
[4491]	654	#if defined key_dynspg_ts
	655	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part
	656	zcorr = tabres(ji,jj) - vn_b(ji,jj)
	657	vb_b(ji,jj) = vb_b(ji,jj) + atfp * zcorr * vmask(ji,jj,1)
	658	END IF
	659	#endif
	660	vn_b(ji,jj) = tabres(ji,jj) * vmask(ji,jj,1)
	661	!
	662	! Correct "before" velocities to hold correct bt component:
	663	spgv(ji,jj) = 0.e0
	664	DO jk=1,jpkm1
	665	spgv(ji,jj) = spgv(ji,jj) + fse3v_b(ji,jj,jk) * vb(ji,jj,jk)
	666	END DO
	667	spgv(ji,jj) = spgv(ji,jj) * hvr_b(ji,jj)
	668	!
	669	zcorr = vb_b(ji,jj) - spgv(ji,jj)
	670	DO jk=1,jpkm1
	671	vb(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk)
	672	END DO
	673	!
[636]	674	END DO
	675	END DO
	676	ENDIF
[5656]	677	!
[636]	678	END SUBROUTINE updatev2d
[390]	679
[5656]	680
[636]	681	SUBROUTINE updateSSH( tabres, i1, i2, j1, j2, before )
	682	!!---------------------------------------------
	683	!! * ROUTINE updateSSH *
	684	!!---------------------------------------------
	685	INTEGER, INTENT(in) :: i1, i2, j1, j2
	686	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	687	LOGICAL, INTENT(in) :: before
[5656]	688	!!
[636]	689	INTEGER :: ji, jj
[5656]	690	!!---------------------------------------------
	691	!
[636]	692	IF (before) THEN
	693	DO jj=j1,j2
	694	DO ji=i1,i2
[4486]	695	tabres(ji,jj) = sshn(ji,jj)
[636]	696	END DO
	697	END DO
	698	ELSE
[4491]	699	#if ! defined key_dynspg_ts
	700	IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN
	701	DO jj=j1,j2
	702	DO ji=i1,i2
[5656]	703	sshb(ji,jj) = sshb(ji,jj) &
	704	& + atfp * ( tabres(ji,jj) - sshn(ji,jj) ) * tmask(ji,jj,1)
[4491]	705	END DO
	706	END DO
	707	ENDIF
	708	#endif
[636]	709	DO jj=j1,j2
	710	DO ji=i1,i2
[4486]	711	sshn(ji,jj) = tabres(ji,jj) * tmask(ji,jj,1)
[636]	712	END DO
	713	END DO
	714	ENDIF
[5656]	715	!
[636]	716	END SUBROUTINE updateSSH
	717
[4486]	718	SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before )
	719	!!---------------------------------------------
	720	!! * ROUTINE updateub2b *
	721	!!---------------------------------------------
	722	INTEGER, INTENT(in) :: i1, i2, j1, j2
	723	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	724	LOGICAL, INTENT(in) :: before
[5656]	725	!!
[4486]	726	INTEGER :: ji, jj
	727	REAL(wp) :: zrhoy
[5656]	728	!!---------------------------------------------
	729	!
[4486]	730	IF (before) THEN
	731	zrhoy = Agrif_Rhoy()
	732	DO jj=j1,j2
	733	DO ji=i1,i2
	734	tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj)
	735	END DO
	736	END DO
	737	tabres = zrhoy * tabres
	738	ELSE
	739	DO jj=j1,j2
	740	DO ji=i1,i2
	741	ub2_b(ji,jj) = tabres(ji,jj) / e2u(ji,jj)
	742	END DO
	743	END DO
	744	ENDIF
[5656]	745	!
[4486]	746	END SUBROUTINE updateub2b
	747
	748	SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before )
	749	!!---------------------------------------------
	750	!! * ROUTINE updatevb2b *
	751	!!---------------------------------------------
	752	INTEGER, INTENT(in) :: i1, i2, j1, j2
	753	REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres
	754	LOGICAL, INTENT(in) :: before
[5656]	755	!!
[4486]	756	INTEGER :: ji, jj
	757	REAL(wp) :: zrhox
[5656]	758	!!---------------------------------------------
	759	!
[4486]	760	IF (before) THEN
	761	zrhox = Agrif_Rhox()
	762	DO jj=j1,j2
	763	DO ji=i1,i2
	764	tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj)
	765	END DO
	766	END DO
	767	tabres = zrhox * tabres
	768	ELSE
	769	DO jj=j1,j2
	770	DO ji=i1,i2
	771	vb2_b(ji,jj) = tabres(ji,jj) / e1v(ji,jj)
	772	END DO
	773	END DO
	774	ENDIF
[5656]	775	!
[4486]	776	END SUBROUTINE updatevb2b
	777
[5656]	778
	779	SUBROUTINE update_scales( tabres, i1, i2, j1, j2, k1, k2, n1,n2, before )
	780	! currently not used
	781	!!---------------------------------------------
	782	!! * ROUTINE updateT *
	783	!!---------------------------------------------
	784	# include "domzgr_substitute.h90"
	785
	786	INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2
	787	REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres
	788	LOGICAL, iNTENT(in) :: before
	789
	790	INTEGER :: ji,jj,jk
	791	REAL(wp) :: ztemp
	792
	793	IF (before) THEN
	794	DO jk=k1,k2
	795	DO jj=j1,j2
	796	DO ji=i1,i2
	797	tabres(ji,jj,jk,1) = e1t(ji,jj)e2t(ji,jj)tmask(ji,jj,jk)
	798	tabres(ji,jj,jk,2) = e1t(ji,jj)*tmask(ji,jj,jk)
	799	tabres(ji,jj,jk,3) = e2t(ji,jj)*tmask(ji,jj,jk)
	800	END DO
	801	END DO
	802	END DO
	803	tabres(:,:,:,1)=tabres(:,:,:,1)Agrif_Rhox()Agrif_Rhoy()
	804	tabres(:,:,:,2)=tabres(:,:,:,2)*Agrif_Rhox()
	805	tabres(:,:,:,3)=tabres(:,:,:,3)*Agrif_Rhoy()
	806	ELSE
	807	DO jk=k1,k2
	808	DO jj=j1,j2
	809	DO ji=i1,i2
	810	IF( tabres(ji,jj,jk,1) .NE. 0. ) THEN
	811	print ,'VAL = ',ji,jj,jk,tabres(ji,jj,jk,1),e1t(ji,jj)e2t(ji,jj)*tmask(ji,jj,jk)
	812	print ,'VAL2 = ',ji,jj,jk,tabres(ji,jj,jk,2),e1t(ji,jj)tmask(ji,jj,jk)
	813	print ,'VAL3 = ',ji,jj,jk,tabres(ji,jj,jk,3),e2t(ji,jj)tmask(ji,jj,jk)
	814	ztemp = sqrt(tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3)))
	815	print *,'CORR = ',ztemp-1.
	816	print ,'NEW VALS = ',tabres(ji,jj,jk,2)ztemp,tabres(ji,jj,jk,3)*ztemp, &
	817	tabres(ji,jj,jk,2)ztemptabres(ji,jj,jk,3)*ztemp
	818	e1t(ji,jj) = tabres(ji,jj,jk,2)*ztemp
	819	e2t(ji,jj) = tabres(ji,jj,jk,3)*ztemp
	820	END IF
	821	END DO
	822	END DO
	823	END DO
	824	ENDIF
	825	!
	826	END SUBROUTINE update_scales
	827
	828	# if defined key_zdftke
	829	SUBROUTINE updateEN( ptab, i1, i2, j1, j2, k1, k2, before )
	830	!!---------------------------------------------
	831	!! * ROUTINE updateen *
	832	!!---------------------------------------------
	833	INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
	834	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
	835	LOGICAL, INTENT(in) :: before
	836	!!---------------------------------------------
	837	!
	838	IF (before) THEN
	839	ptab (i1:i2,j1:j2,k1:k2) = en(i1:i2,j1:j2,k1:k2)
	840	ELSE
	841	en(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2)
	842	ENDIF
	843	!
	844	END SUBROUTINE updateEN
	845
	846
	847	SUBROUTINE updateAVT( ptab, i1, i2, j1, j2, k1, k2, before )
	848	!!---------------------------------------------
	849	!! * ROUTINE updateavt *
	850	!!---------------------------------------------
	851	INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
	852	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
	853	LOGICAL, INTENT(in) :: before
	854	!!---------------------------------------------
	855	!
	856	IF (before) THEN
	857	ptab (i1:i2,j1:j2,k1:k2) = avt_k(i1:i2,j1:j2,k1:k2)
	858	ELSE
	859	avt_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2)
	860	ENDIF
	861	!
	862	END SUBROUTINE updateAVT
	863
	864
	865	SUBROUTINE updateAVM( ptab, i1, i2, j1, j2, k1, k2, before )
	866	!!---------------------------------------------
	867	!! * ROUTINE updateavm *
	868	!!---------------------------------------------
	869	INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2
	870	REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab
	871	LOGICAL, INTENT(in) :: before
	872	!!---------------------------------------------
	873	!
	874	IF (before) THEN
	875	ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2)
	876	ELSE
	877	avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2)
	878	ENDIF
	879	!
	880	END SUBROUTINE updateAVM
	881
	882	# endif /* key_zdftke */
	883
[390]	884	#else
[636]	885	CONTAINS
	886	SUBROUTINE agrif_opa_update_empty
	887	!!---------------------------------------------
	888	!! * ROUTINE agrif_opa_update_empty *
	889	!!---------------------------------------------
	890	WRITE(,) 'agrif_opa_update : You should not have seen this print! error?'
	891	END SUBROUTINE agrif_opa_update_empty
[390]	892	#endif
[636]	893	END MODULE agrif_opa_update
[4486]	894

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

Download in other formats: