Home

Context Navigation

source: palm/trunk/SOURCE/diffusion_v.f90 @ 1271

Last change on this file since 1271 was 1258, checked in by raasch, 11 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 24.3 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |