Home

Context Navigation

diffusion_v_mod.f90 @ 1867

Last change on this file since 1867 was 1851, checked in by maronga, 9 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 26.6 KB

Rev	Line
[1850]	1	!> @file diffusion_v_mod.f90
[1036]	2	!--------------------------------------------------------------------------------!
	3	! This file is part of PALM.
	4	!
	5	! PALM is free software: you can redistribute it and/or modify it under the terms
	6	! of the GNU General Public License as published by the Free Software Foundation,
	7	! either version 3 of the License, or (at your option) any later version.
	8	!
	9	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
	10	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
	11	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
	12	!
	13	! You should have received a copy of the GNU General Public License along with
	14	! PALM. If not, see <http://www.gnu.org/licenses/>.
	15	!
[1818]	16	! Copyright 1997-2016 Leibniz Universitaet Hannover
[1036]	17	!--------------------------------------------------------------------------------!
	18	!
[484]	19	! Current revisions:
[1]	20	! -----------------
[1341]	21	!
[1851]	22	!
[1321]	23	! Former revisions:
	24	! -----------------
	25	! $Id: diffusion_v_mod.f90 1851 2016-04-08 13:32:50Z raasch $
	26	!
[1851]	27	! 1850 2016-04-08 13:29:27Z maronga
	28	! Module renamed
	29	!
	30	!
[1741]	31	! 1740 2016-01-13 08:19:40Z raasch
	32	! unnecessary calculations of kmzm and kmzp in wall bounded parts removed
	33	!
[1683]	34	! 1682 2015-10-07 23:56:08Z knoop
	35	! Code annotations made doxygen readable
	36	!
[1341]	37	! 1340 2014-03-25 19:45:13Z kanani
	38	! REAL constants defined as wp-kind
	39	!
[1321]	40	! 1320 2014-03-20 08:40:49Z raasch
[1320]	41	! ONLY-attribute added to USE-statements,
	42	! kind-parameters added to all INTEGER and REAL declaration statements,
	43	! kinds are defined in new module kinds,
	44	! revision history before 2012 removed,
	45	! comment fields (!:) to be used for variable explanations added to
	46	! all variable declaration statements
[1321]	47	!
[1258]	48	! 1257 2013-11-08 15:18:40Z raasch
	49	! openacc loop and loop vector clauses removed, declare create moved after
	50	! the FORTRAN declaration statement
	51	!
[1132]	52	! 1128 2013-04-12 06:19:32Z raasch
	53	! loop index bounds in accelerator version replaced by i_left, i_right, j_south,
	54	! j_north
	55	!
[1037]	56	! 1036 2012-10-22 13:43:42Z raasch
	57	! code put under GPL (PALM 3.9)
	58	!
[1017]	59	! 1015 2012-09-27 09:23:24Z raasch
	60	! accelerator version (*_acc) added
	61	!
[1002]	62	! 1001 2012-09-13 14:08:46Z raasch
	63	! arrays comunicated by module instead of parameter list
	64	!
[979]	65	! 978 2012-08-09 08:28:32Z fricke
	66	! outflow damping layer removed
	67	! kmxm_x/_y and kmxp_x/_y change to kmxm and kmxp
	68	!
[1]	69	! Revision 1.1 1997/09/12 06:24:01 raasch
	70	! Initial revision
	71	!
	72	!
	73	! Description:
	74	! ------------
[1682]	75	!> Diffusion term of the v-component
[1]	76	!------------------------------------------------------------------------------!
[1682]	77	MODULE diffusion_v_mod
	78
[1]	79
[56]	80	USE wall_fluxes_mod
	81
[1]	82	PRIVATE
[1015]	83	PUBLIC diffusion_v, diffusion_v_acc
[1]	84
	85	INTERFACE diffusion_v
	86	MODULE PROCEDURE diffusion_v
	87	MODULE PROCEDURE diffusion_v_ij
	88	END INTERFACE diffusion_v
	89
[1015]	90	INTERFACE diffusion_v_acc
	91	MODULE PROCEDURE diffusion_v_acc
	92	END INTERFACE diffusion_v_acc
	93
[1]	94	CONTAINS
	95
	96
	97	!------------------------------------------------------------------------------!
[1682]	98	! Description:
	99	! ------------
	100	!> Call for all grid points
[1]	101	!------------------------------------------------------------------------------!
[1001]	102	SUBROUTINE diffusion_v
[1]	103
[1320]	104	USE arrays_3d, &
	105	ONLY: ddzu, ddzw, km, tend, u, v, vsws, vswst, w
	106
	107	USE control_parameters, &
	108	ONLY: constant_top_momentumflux, topography, use_surface_fluxes, &
	109	use_top_fluxes
	110
	111	USE grid_variables, &
	112	ONLY: ddx, ddy, ddy2, fxm, fxp, wall_v
	113
	114	USE indices, &
	115	ONLY: nxl, nxr, nyn, nys, nysv, nzb, nzb_diff_v, nzb_v_inner, &
	116	nzb_v_outer, nzt, nzt_diff
	117
	118	USE kinds
[1]	119
	120	IMPLICIT NONE
	121
[1682]	122	INTEGER(iwp) :: i !<
	123	INTEGER(iwp) :: j !<
	124	INTEGER(iwp) :: k !<
	125	REAL(wp) :: kmxm !<
	126	REAL(wp) :: kmxp !<
	127	REAL(wp) :: kmzm !<
	128	REAL(wp) :: kmzp !<
[1001]	129
[1682]	130	REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: vsus !<
[1]	131
[56]	132	!
	133	!-- First calculate horizontal momentum flux v'u' at vertical walls,
	134	!-- if neccessary
	135	IF ( topography /= 'flat' ) THEN
[1320]	136	CALL wall_fluxes( vsus, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, nzb_v_inner, &
[56]	137	nzb_v_outer, wall_v )
	138	ENDIF
	139
[1]	140	DO i = nxl, nxr
[106]	141	DO j = nysv, nyn
[1]	142	!
	143	!-- Compute horizontal diffusion
	144	DO k = nzb_v_outer(j,i)+1, nzt
	145	!
	146	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	147	kmxp = 0.25_wp * &
[978]	148	( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
[1340]	149	kmxm = 0.25_wp * &
[978]	150	( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
[1]	151
[1320]	152	tend(k,j,i) = tend(k,j,i) &
	153	& + ( kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
	154	& + kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
	155	& - kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
	156	& - kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
	157	& ) * ddx &
[1340]	158	& + 2.0_wp * ( &
[1320]	159	& km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
	160	& - km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
[1340]	161	& ) * ddy2
[1]	162	ENDDO
	163
	164	!
	165	!-- Wall functions at the left and right walls, respectively
[1340]	166	IF ( wall_v(j,i) /= 0.0_wp ) THEN
[51]	167
[1]	168	DO k = nzb_v_inner(j,i)+1, nzb_v_outer(j,i)
[1340]	169	kmxp = 0.25_wp * &
[978]	170	( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
[1340]	171	kmxm = 0.25_wp * &
[978]	172	( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
	173
[1]	174	tend(k,j,i) = tend(k,j,i) &
[1340]	175	+ 2.0_wp * ( &
[1]	176	km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
	177	- km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
[1340]	178	) * ddy2 &
[1]	179	+ ( fxp(j,i) * ( &
[978]	180	kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
	181	+ kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
[1]	182	) &
	183	- fxm(j,i) * ( &
[978]	184	kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
	185	+ kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
[1]	186	) &
[56]	187	+ wall_v(j,i) * vsus(k,j,i) &
[1]	188	) * ddx
	189	ENDDO
	190	ENDIF
	191
	192	!
	193	!-- Compute vertical diffusion. In case of simulating a Prandtl
	194	!-- layer, index k starts at nzb_v_inner+2.
[102]	195	DO k = nzb_diff_v(j,i), nzt_diff
[1]	196	!
	197	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	198	kmzp = 0.25_wp * &
[1]	199	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
[1340]	200	kmzm = 0.25_wp * &
[1]	201	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
	202
[1320]	203	tend(k,j,i) = tend(k,j,i) &
	204	& + ( kmzp * ( ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
	205	& + ( w(k,j,i) - w(k,j-1,i) ) * ddy &
	206	& ) &
	207	& - kmzm * ( ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
	208	& + ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
	209	& ) &
[1]	210	& ) * ddzw(k)
	211	ENDDO
	212
	213	!
	214	!-- Vertical diffusion at the first grid point above the surface,
	215	!-- if the momentum flux at the bottom is given by the Prandtl law
	216	!-- or if it is prescribed by the user.
	217	!-- Difference quotient of the momentum flux is not formed over
	218	!-- half of the grid spacing (2.0*ddzw(k)) any more, since the
[1320]	219	!-- comparison with other (LES) models showed that the values of
[1]	220	!-- the momentum flux becomes too large in this case.
	221	!-- The term containing w(k-1,..) (see above equation) is removed here
	222	!-- because the vertical velocity is assumed to be zero at the surface.
	223	IF ( use_surface_fluxes ) THEN
	224	k = nzb_v_inner(j,i)+1
	225	!
	226	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	227	kmzp = 0.25_wp * &
[1]	228	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
	229
[1320]	230	tend(k,j,i) = tend(k,j,i) &
	231	& + ( kmzp * ( w(k,j,i) - w(k,j-1,i) ) * ddy &
	232	& ) * ddzw(k) &
	233	& + ( kmzp * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
	234	& + vsws(j,i) &
[1]	235	& ) * ddzw(k)
	236	ENDIF
	237
[102]	238	!
	239	!-- Vertical diffusion at the first gridpoint below the top boundary,
	240	!-- if the momentum flux at the top is prescribed by the user
[103]	241	IF ( use_top_fluxes .AND. constant_top_momentumflux ) THEN
[102]	242	k = nzt
	243	!
	244	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	245	kmzm = 0.25_wp * &
[102]	246	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
	247
[1320]	248	tend(k,j,i) = tend(k,j,i) &
	249	& - ( kmzm * ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
	250	& ) * ddzw(k) &
	251	& + ( -vswst(j,i) &
	252	& - kmzm * ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
[102]	253	& ) * ddzw(k)
	254	ENDIF
	255
[1]	256	ENDDO
	257	ENDDO
	258
	259	END SUBROUTINE diffusion_v
	260
	261
	262	!------------------------------------------------------------------------------!
[1682]	263	! Description:
	264	! ------------
	265	!> Call for all grid points - accelerator version
[1015]	266	!------------------------------------------------------------------------------!
	267	SUBROUTINE diffusion_v_acc
	268
[1320]	269	USE arrays_3d, &
	270	ONLY: ddzu, ddzw, km, tend, u, v, vsws, vswst, w
	271
	272	USE control_parameters, &
	273	ONLY: constant_top_momentumflux, topography, use_surface_fluxes, &
	274	use_top_fluxes
	275
	276	USE grid_variables, &
	277	ONLY: ddx, ddy, ddy2, fxm, fxp, wall_v
	278
	279	USE indices, &
	280	ONLY: i_left, i_right, j_north, j_south, nxl, nxr, nyn, nys, nzb, &
	281	nzb_diff_v, nzb_v_inner, nzb_v_outer, nzt, nzt_diff
	282
	283	USE kinds
[1015]	284
	285	IMPLICIT NONE
	286
[1682]	287	INTEGER(iwp) :: i !<
	288	INTEGER(iwp) :: j !<
	289	INTEGER(iwp) :: k !<
	290	REAL(wp) :: kmxm !<
	291	REAL(wp) :: kmxp !<
	292	REAL(wp) :: kmzm !<
	293	REAL(wp) :: kmzp !<
[1015]	294
[1682]	295	REAL(wp), DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: vsus !<
[1015]	296	!$acc declare create ( vsus )
	297
	298	!
	299	!-- First calculate horizontal momentum flux v'u' at vertical walls,
	300	!-- if neccessary
	301	IF ( topography /= 'flat' ) THEN
[1320]	302	CALL wall_fluxes_acc( vsus, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp, &
	303	nzb_v_inner, nzb_v_outer, wall_v )
[1015]	304	ENDIF
	305
[1320]	306	!$acc kernels present ( u, v, w, km, tend, vsws, vswst ) &
	307	!$acc present ( ddzu, ddzw, fxm, fxp, wall_v ) &
[1015]	308	!$acc present ( nzb_v_inner, nzb_v_outer, nzb_diff_v )
[1128]	309	DO i = i_left, i_right
	310	DO j = j_south, j_north
[1015]	311	!
	312	!-- Compute horizontal diffusion
	313	DO k = 1, nzt
	314	IF ( k > nzb_v_outer(j,i) ) THEN
	315	!
	316	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	317	kmxp = 0.25_wp * &
[1015]	318	( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
[1340]	319	kmxm = 0.25_wp * &
[1015]	320	( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
	321
	322	tend(k,j,i) = tend(k,j,i) &
	323	& + ( kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
	324	& + kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
	325	& - kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
	326	& - kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
	327	& ) * ddx &
[1340]	328	& + 2.0_wp * ( &
[1015]	329	& km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
	330	& - km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
[1340]	331	& ) * ddy2
[1015]	332	ENDIF
	333	ENDDO
	334
	335	!
	336	!-- Wall functions at the left and right walls, respectively
	337	DO k = 1, nzt
[1320]	338	IF( k > nzb_v_inner(j,i) .AND. k <= nzb_v_outer(j,i) .AND. &
[1340]	339	wall_v(j,i) /= 0.0_wp ) THEN
[1015]	340
[1340]	341	kmxp = 0.25_wp * &
[1015]	342	( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
[1340]	343	kmxm = 0.25_wp * &
[1015]	344	( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
	345
	346	tend(k,j,i) = tend(k,j,i) &
[1340]	347	+ 2.0_wp * ( &
[1015]	348	km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
	349	- km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
[1340]	350	) * ddy2 &
[1015]	351	+ ( fxp(j,i) * ( &
	352	kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
	353	+ kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
	354	) &
	355	- fxm(j,i) * ( &
	356	kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
	357	+ kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
	358	) &
	359	+ wall_v(j,i) * vsus(k,j,i) &
	360	) * ddx
	361	ENDIF
	362	ENDDO
	363
	364	!
	365	!-- Compute vertical diffusion. In case of simulating a Prandtl
	366	!-- layer, index k starts at nzb_v_inner+2.
	367	DO k = 1, nzt_diff
	368	IF ( k >= nzb_diff_v(j,i) ) THEN
	369	!
	370	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	371	kmzp = 0.25_wp * &
[1015]	372	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
[1340]	373	kmzm = 0.25_wp * &
[1015]	374	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
	375
	376	tend(k,j,i) = tend(k,j,i) &
	377	& + ( kmzp * ( ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1)&
	378	& + ( w(k,j,i) - w(k,j-1,i) ) * ddy &
	379	& ) &
	380	& - kmzm * ( ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k)&
	381	& + ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
	382	& ) &
	383	& ) * ddzw(k)
	384	ENDIF
	385	ENDDO
	386
	387	ENDDO
	388	ENDDO
	389
	390	!
	391	!-- Vertical diffusion at the first grid point above the surface,
	392	!-- if the momentum flux at the bottom is given by the Prandtl law
	393	!-- or if it is prescribed by the user.
	394	!-- Difference quotient of the momentum flux is not formed over
	395	!-- half of the grid spacing (2.0*ddzw(k)) any more, since the
[1320]	396	!-- comparison with other (LES) models showed that the values of
[1015]	397	!-- the momentum flux becomes too large in this case.
	398	!-- The term containing w(k-1,..) (see above equation) is removed here
	399	!-- because the vertical velocity is assumed to be zero at the surface.
	400	IF ( use_surface_fluxes ) THEN
	401
[1128]	402	DO i = i_left, i_right
	403	DO j = j_south, j_north
[1015]	404
	405	k = nzb_v_inner(j,i)+1
	406	!
	407	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	408	kmzp = 0.25_wp * &
[1015]	409	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
	410
[1320]	411	tend(k,j,i) = tend(k,j,i) &
	412	& + ( kmzp * ( w(k,j,i) - w(k,j-1,i) ) * ddy &
	413	& ) * ddzw(k) &
	414	& + ( kmzp * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
	415	& + vsws(j,i) &
[1015]	416	& ) * ddzw(k)
	417	ENDDO
	418	ENDDO
	419
	420	ENDIF
	421
	422	!
	423	!-- Vertical diffusion at the first gridpoint below the top boundary,
	424	!-- if the momentum flux at the top is prescribed by the user
	425	IF ( use_top_fluxes .AND. constant_top_momentumflux ) THEN
	426
	427	k = nzt
	428
[1128]	429	DO i = i_left, i_right
	430	DO j = j_south, j_north
[1015]	431
	432	!
	433	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	434	kmzm = 0.25_wp * &
[1015]	435	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
	436
[1320]	437	tend(k,j,i) = tend(k,j,i) &
	438	& - ( kmzm * ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
	439	& ) * ddzw(k) &
	440	& + ( -vswst(j,i) &
	441	& - kmzm * ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
[1015]	442	& ) * ddzw(k)
	443	ENDDO
	444	ENDDO
	445
	446	ENDIF
	447	!$acc end kernels
	448
	449	END SUBROUTINE diffusion_v_acc
	450
	451
	452	!------------------------------------------------------------------------------!
[1682]	453	! Description:
	454	! ------------
	455	!> Call for grid point i,j
[1]	456	!------------------------------------------------------------------------------!
[1001]	457	SUBROUTINE diffusion_v_ij( i, j )
[1]	458
[1320]	459	USE arrays_3d, &
	460	ONLY: ddzu, ddzw, km, tend, u, v, vsws, vswst, w
	461
	462	USE control_parameters, &
	463	ONLY: constant_top_momentumflux, use_surface_fluxes, use_top_fluxes
	464
	465	USE grid_variables, &
	466	ONLY: ddx, ddy, ddy2, fxm, fxp, wall_v
	467
	468	USE indices, &
	469	ONLY: nzb, nzb_diff_v, nzb_v_inner, nzb_v_outer, nzt, nzt_diff
	470
	471	USE kinds
[1]	472
	473	IMPLICIT NONE
	474
[1682]	475	INTEGER(iwp) :: i !<
	476	INTEGER(iwp) :: j !<
	477	INTEGER(iwp) :: k !<
	478	REAL(wp) :: kmxm !<
	479	REAL(wp) :: kmxp !<
	480	REAL(wp) :: kmzm !<
	481	REAL(wp) :: kmzp !<
[1]	482
[1682]	483	REAL(wp), DIMENSION(nzb:nzt+1) :: vsus !<
[1001]	484
[1]	485	!
	486	!-- Compute horizontal diffusion
	487	DO k = nzb_v_outer(j,i)+1, nzt
	488	!
	489	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	490	kmxp = 0.25_wp * ( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
	491	kmxm = 0.25_wp * ( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
[1]	492
[1320]	493	tend(k,j,i) = tend(k,j,i) &
	494	& + ( kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
	495	& + kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
	496	& - kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
	497	& - kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
	498	& ) * ddx &
[1340]	499	& + 2.0_wp * ( &
[1320]	500	& km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
	501	& - km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
[1340]	502	& ) * ddy2
[1]	503	ENDDO
	504
	505	!
	506	!-- Wall functions at the left and right walls, respectively
[1340]	507	IF ( wall_v(j,i) /= 0.0_wp ) THEN
[51]	508
	509	!
	510	!-- Calculate the horizontal momentum flux v'u'
[1320]	511	CALL wall_fluxes( i, j, nzb_v_inner(j,i)+1, nzb_v_outer(j,i), &
	512	vsus, 0.0_wp, 1.0_wp, 0.0_wp, 0.0_wp )
[51]	513
[1]	514	DO k = nzb_v_inner(j,i)+1, nzb_v_outer(j,i)
[1340]	515	kmxp = 0.25_wp * &
[978]	516	( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
[1340]	517	kmxm = 0.25_wp * &
[978]	518	( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
[1]	519
	520	tend(k,j,i) = tend(k,j,i) &
[1340]	521	+ 2.0_wp * ( &
[1]	522	km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
	523	- km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
[1340]	524	) * ddy2 &
[1]	525	+ ( fxp(j,i) * ( &
[978]	526	kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
	527	+ kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
[1]	528	) &
	529	- fxm(j,i) * ( &
[978]	530	kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
	531	+ kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
[1]	532	) &
[51]	533	+ wall_v(j,i) * vsus(k) &
[1]	534	) * ddx
	535	ENDDO
	536	ENDIF
	537
	538	!
	539	!-- Compute vertical diffusion. In case of simulating a Prandtl layer,
	540	!-- index k starts at nzb_v_inner+2.
[102]	541	DO k = nzb_diff_v(j,i), nzt_diff
[1]	542	!
	543	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	544	kmzp = 0.25_wp * ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
	545	kmzm = 0.25_wp * ( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
[1]	546
[1320]	547	tend(k,j,i) = tend(k,j,i) &
	548	& + ( kmzp * ( ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
	549	& + ( w(k,j,i) - w(k,j-1,i) ) * ddy &
	550	& ) &
	551	& - kmzm * ( ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
	552	& + ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
	553	& ) &
[1]	554	& ) * ddzw(k)
	555	ENDDO
	556
	557	!
	558	!-- Vertical diffusion at the first grid point above the surface, if the
	559	!-- momentum flux at the bottom is given by the Prandtl law or if it is
	560	!-- prescribed by the user.
	561	!-- Difference quotient of the momentum flux is not formed over half of
	562	!-- the grid spacing (2.0*ddzw(k)) any more, since the comparison with
[1320]	563	!-- other (LES) models showed that the values of the momentum flux becomes
[1]	564	!-- too large in this case.
	565	!-- The term containing w(k-1,..) (see above equation) is removed here
	566	!-- because the vertical velocity is assumed to be zero at the surface.
	567	IF ( use_surface_fluxes ) THEN
	568	k = nzb_v_inner(j,i)+1
	569	!
	570	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	571	kmzp = 0.25_wp * ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
[1]	572
[1320]	573	tend(k,j,i) = tend(k,j,i) &
	574	& + ( kmzp * ( w(k,j,i) - w(k,j-1,i) ) * ddy &
	575	& ) * ddzw(k) &
	576	& + ( kmzp * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
	577	& + vsws(j,i) &
[1]	578	& ) * ddzw(k)
	579	ENDIF
	580
[102]	581	!
	582	!-- Vertical diffusion at the first gridpoint below the top boundary,
	583	!-- if the momentum flux at the top is prescribed by the user
[103]	584	IF ( use_top_fluxes .AND. constant_top_momentumflux ) THEN
[102]	585	k = nzt
	586	!
	587	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	588	kmzm = 0.25_wp * &
[102]	589	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
	590
[1320]	591	tend(k,j,i) = tend(k,j,i) &
	592	& - ( kmzm * ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
	593	& ) * ddzw(k) &
	594	& + ( -vswst(j,i) &
	595	& - kmzm * ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
[102]	596	& ) * ddzw(k)
	597	ENDIF
	598
[1]	599	END SUBROUTINE diffusion_v_ij
	600
[1321]	601	END MODULE diffusion_v_mod

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

Context Navigation

source: palm/trunk/SOURCE/diffusion_v_mod.f90 @ 1867

Download in other formats: