Home

Context Navigation

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

Last change on this file since 4672 was 4671, checked in by pavelkrc, 4 years ago
Radiative transfer model RTM version 4.1
Property svn:keywords set to `Id`
File size: 25.2 KB

Rev	Line
[1873]	1	!> @file diffusion_v.f90
[4583]	2	!--------------------------------------------------------------------------------------------------!
[2696]	3	! This file is part of the PALM model system.
[1036]	4	!
[4583]	5	! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
	6	! Public License as published by the Free Software Foundation, either version 3 of the License, or
	7	! (at your option) any later version.
[1036]	8	!
[4583]	9	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
	10	! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
	11	! Public License for more details.
[1036]	12	!
[4583]	13	! You should have received a copy of the GNU General Public License along with PALM. If not, see
	14	! <http://www.gnu.org/licenses/>.
[1036]	15	!
[4360]	16	! Copyright 1997-2020 Leibniz Universitaet Hannover
[4583]	17	!--------------------------------------------------------------------------------------------------!
[1036]	18	!
[484]	19	! Current revisions:
[1]	20	! -----------------
[4671]	21	!
	22	!
[1321]	23	! Former revisions:
	24	! -----------------
	25	! $Id: diffusion_v.f90 4671 2020-09-09 20:27:58Z pavelkrc $
[4671]	26	! Implementation of downward facing USM and LSM surfaces
	27	!
	28	! 4583 2020-06-29 12:36:47Z raasch
[4583]	29	! file re-formatted to follow the PALM coding standard
	30	!
	31	! 4360 2020-01-07 11:25:50Z suehring
	32	! Introduction of wall_flags_total_0, which currently sets bits based on static topography
	33	! information used in wall_flags_static_0
	34	!
[4346]	35	! 4329 2019-12-10 15:46:36Z motisi
[4329]	36	! Renamed wall_flags_0 to wall_flags_static_0
[4583]	37	!
[4329]	38	! 4182 2019-08-22 15:20:23Z scharf
[4182]	39	! Corrected "Former revisions" section
[4583]	40	!
[4182]	41	! 3655 2019-01-07 16:51:22Z knoop
[3634]	42	! OpenACC port for SPEC
[2716]	43	!
[4182]	44	! Revision 1.1 1997/09/12 06:24:01 raasch
	45	! Initial revision
	46	!
	47	!
[1]	48	! Description:
	49	! ------------
[1682]	50	!> Diffusion term of the v-component
[4583]	51	!--------------------------------------------------------------------------------------------------!
[1682]	52	MODULE diffusion_v_mod
[1]	53
[4583]	54
[1]	55	PRIVATE
[2118]	56	PUBLIC diffusion_v
[1]	57
	58	INTERFACE diffusion_v
	59	MODULE PROCEDURE diffusion_v
	60	MODULE PROCEDURE diffusion_v_ij
	61	END INTERFACE diffusion_v
	62
	63	CONTAINS
	64
	65
[4583]	66	!--------------------------------------------------------------------------------------------------!
[1682]	67	! Description:
	68	! ------------
	69	!> Call for all grid points
[4583]	70	!--------------------------------------------------------------------------------------------------!
[1001]	71	SUBROUTINE diffusion_v
[1]	72
[4583]	73	USE arrays_3d, &
	74	ONLY: ddzu, ddzw, drho_air, km, rho_air_zw, tend, u, v, w
	75
	76	USE control_parameters, &
	77	ONLY: constant_top_momentumflux, use_surface_fluxes, use_top_fluxes
	78
	79	USE grid_variables, &
[2232]	80	ONLY: ddx, ddy, ddy2
[4583]	81
	82	USE indices, &
[4346]	83	ONLY: nxl, nxr, nyn, nysv, nzb, nzt, wall_flags_total_0
[4583]	84
[1320]	85	USE kinds
[1]	86
[4583]	87	USE surface_mod, &
	88	ONLY : surf_def_h, surf_def_v, surf_lsm_h, surf_lsm_v, surf_usm_h, surf_usm_v
[2232]	89
[1]	90	IMPLICIT NONE
	91
[2232]	92	INTEGER(iwp) :: i !< running index x direction
	93	INTEGER(iwp) :: j !< running index y direction
	94	INTEGER(iwp) :: k !< running index z direction
	95	INTEGER(iwp) :: l !< running index of surface type, south- or north-facing wall
	96	INTEGER(iwp) :: m !< running index surface elements
	97	INTEGER(iwp) :: surf_e !< End index of surface elements at (j,i)-gridpoint
	98	INTEGER(iwp) :: surf_s !< Start index of surface elements at (j,i)-gridpoint
[1001]	99
[2232]	100	REAL(wp) :: flag !< flag to mask topography grid points
[3547]	101	REAL(wp) :: kmxm !< diffusion coefficient on leftward side of the v-gridbox - interpolated onto xu-yv grid
	102	REAL(wp) :: kmxp !< diffusion coefficient on rightward side of the v-gridbox - interpolated onto xu-yv grid
	103	REAL(wp) :: kmzm !< diffusion coefficient on bottom of the gridbox - interpolated onto yv-zw grid
	104	REAL(wp) :: kmzp !< diffusion coefficient on top of the gridbox - interpolated onto yv-zw grid
[4583]	105	REAL(wp) :: mask_bottom !< flag to mask vertical upward-facing surface
	106	REAL(wp) :: mask_east !< flag to mask vertical surface south of the grid point
	107	REAL(wp) :: mask_west !< flag to mask vertical surface north of the grid point
	108	REAL(wp) :: mask_top !< flag to mask vertical downward-facing surface
[1]	109
[3634]	110	!$ACC PARALLEL LOOP COLLAPSE(2) PRIVATE(i, j, k, l, m) &
	111	!$ACC PRIVATE(surf_e, surf_s, flag, kmxm, kmxp, kmzm, kmzp) &
	112	!$ACC PRIVATE(mask_bottom, mask_east, mask_west, mask_top) &
[4346]	113	!$ACC PRESENT(wall_flags_total_0, km) &
[3634]	114	!$ACC PRESENT(u, v, w) &
	115	!$ACC PRESENT(ddzu, ddzw, drho_air, rho_air_zw) &
	116	!$ACC PRESENT(surf_def_h(0:2), surf_def_v(2:3)) &
	117	!$ACC PRESENT(surf_lsm_h, surf_lsm_v(2:3)) &
	118	!$ACC PRESENT(surf_usm_h, surf_usm_v(0:3)) &
	119	!$ACC PRESENT(tend)
[1]	120	DO i = nxl, nxr
[106]	121	DO j = nysv, nyn
[1]	122	!
	123	!-- Compute horizontal diffusion
[2232]	124	DO k = nzb+1, nzt
	125
[1]	126	!
[4583]	127	!-- Predetermine flag to mask topography and wall-bounded grid points.
	128	!-- It is sufficient to masked only east- and west-facing surfaces, which need special
	129	!-- treatment for the v-component.
	130	flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2 ) )
[4346]	131	mask_east = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i+1), 2 ) )
	132	mask_west = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i-1), 2 ) )
[2232]	133	!
[1]	134	!-- Interpolate eddy diffusivities on staggered gridpoints
[2232]	135	kmxp = 0.25_wp * ( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
	136	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]	137
[4583]	138	tend(k,j,i) = tend(k,j,i) + ( &
	139	mask_east * kmxp * ( &
	140	( v(k,j,i+1) - v(k,j,i) ) * ddx &
	141	+ ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
	142	) &
	143	- mask_west * kmxm * ( &
	144	( v(k,j,i) - v(k,j,i-1) ) * ddx &
	145	+ ( u(k,j,i) - u(k,j-1,i) ) * ddy &
	146	) &
	147	) * ddx * flag &
	148	+ 2.0_wp * ( &
	149	km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
	150	- km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
	151	) * ddy2 * flag
[2232]	152
[1]	153	ENDDO
	154
	155	!
[4583]	156	!-- Add horizontal momentum flux v'u' at east- (l=2) and west-facing (l=3) surfaces. Note,
	157	!-- in the the flat case, loops won't be entered as start_index > end_index. Furtermore,
	158	!-- note, no vertical natural surfaces so far.
[2232]	159	!-- Default-type surfaces
	160	DO l = 2, 3
	161	surf_s = surf_def_v(l)%start_index(j,i)
	162	surf_e = surf_def_v(l)%end_index(j,i)
	163	DO m = surf_s, surf_e
	164	k = surf_def_v(l)%k(m)
[4583]	165	tend(k,j,i) = tend(k,j,i) + surf_def_v(l)%mom_flux_uv(m) * ddx
	166	ENDDO
[2232]	167	ENDDO
[1]	168	!
[2232]	169	!-- Natural-type surfaces
	170	DO l = 2, 3
	171	surf_s = surf_lsm_v(l)%start_index(j,i)
	172	surf_e = surf_lsm_v(l)%end_index(j,i)
	173	DO m = surf_s, surf_e
	174	k = surf_lsm_v(l)%k(m)
[4583]	175	tend(k,j,i) = tend(k,j,i) + surf_lsm_v(l)%mom_flux_uv(m) * ddx
	176	ENDDO
[2232]	177	ENDDO
[1]	178	!
[2232]	179	!-- Urban-type surfaces
	180	DO l = 2, 3
	181	surf_s = surf_usm_v(l)%start_index(j,i)
	182	surf_e = surf_usm_v(l)%end_index(j,i)
	183	DO m = surf_s, surf_e
	184	k = surf_usm_v(l)%k(m)
[4583]	185	tend(k,j,i) = tend(k,j,i) + surf_usm_v(l)%mom_flux_uv(m) * ddx
	186	ENDDO
[2232]	187	ENDDO
	188	!
[4583]	189	!-- Compute vertical diffusion. In case of simulating a surface layer, respective grid
	190	!-- diffusive fluxes are masked (flag 10) within this loop, and added further below, else,
	191	!-- simple gradient approach is applied. Model top is also mask if top-momentum flux is
	192	!-- given.
[2232]	193	DO k = nzb+1, nzt
	194	!
[4583]	195	!-- Determine flags to mask topography below and above. Flag 2 is used to mask
	196	!-- topography in general, while flag 8 implies also information about
	197	!-- use_surface_fluxes. Flag 9 is used to control momentum flux at model top.
	198	mask_bottom = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k-1,j,i), 8 ) )
	199	mask_top = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k+1,j,i), 8 ) ) * &
	200	MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k+1,j,i), 9 ) )
	201	flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2 ) )
[2232]	202	!
[1]	203	!-- Interpolate eddy diffusivities on staggered gridpoints
[4583]	204	kmzp = 0.25_wp * ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
	205	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]	206
[4583]	207	tend(k,j,i) = tend(k,j,i) &
	208	& + ( kmzp * ( ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
	209	& + ( w(k,j,i) - w(k,j-1,i) ) * ddy &
	210	& ) * rho_air_zw(k) * mask_top &
	211	& - kmzm * ( ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
	212	& + ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
	213	& ) * rho_air_zw(k-1) * mask_bottom &
	214	& ) * ddzw(k) * drho_air(k) * flag
[1]	215	ENDDO
	216
	217	!
[4583]	218	!-- Vertical diffusion at the first grid point above the surface, if the momentum flux at
	219	!-- the bottom is given by the Prandtl law or if it is prescribed by the user.
	220	!-- Difference quotient of the momentum flux is not formed over half of the grid spacing
	221	!-- (2.0*ddzw(k)) any more, since the comparison with other (LES) models showed that the
	222	!-- values of the momentum flux becomes too large in this case.
[1]	223	IF ( use_surface_fluxes ) THEN
	224	!
[2232]	225	!-- Default-type surfaces, upward-facing
	226	surf_s = surf_def_h(0)%start_index(j,i)
	227	surf_e = surf_def_h(0)%end_index(j,i)
	228	DO m = surf_s, surf_e
	229	k = surf_def_h(0)%k(m)
[1]	230
[4583]	231	tend(k,j,i) = tend(k,j,i) &
	232	+ ( - ( - surf_def_h(0)%vsws(m) ) ) * ddzw(k) * drho_air(k)
[2232]	233	ENDDO
	234	!
	235	!-- Default-type surfaces, dowward-facing
	236	surf_s = surf_def_h(1)%start_index(j,i)
	237	surf_e = surf_def_h(1)%end_index(j,i)
	238	DO m = surf_s, surf_e
	239	k = surf_def_h(1)%k(m)
[1]	240
[4583]	241	tend(k,j,i) = tend(k,j,i) &
	242	+ ( - surf_def_h(1)%vsws(m) ) * ddzw(k) * drho_air(k)
[2232]	243	ENDDO
[102]	244	!
[2232]	245	!-- Natural-type surfaces, upward-facing
[4671]	246	surf_s = surf_lsm_h(0)%start_index(j,i)
	247	surf_e = surf_lsm_h(0)%end_index(j,i)
[2232]	248	DO m = surf_s, surf_e
[4671]	249	k = surf_lsm_h(0)%k(m)
[2232]	250
[4583]	251	tend(k,j,i) = tend(k,j,i) &
[4671]	252	+ ( - ( - surf_lsm_h(0)%vsws(m) ) ) * ddzw(k) * drho_air(k)
[2232]	253
	254	ENDDO
[102]	255	!
[4671]	256	!-- Natural-type surfaces, downward-facing
	257	surf_s = surf_lsm_h(1)%start_index(j,i)
	258	surf_e = surf_lsm_h(1)%end_index(j,i)
	259	DO m = surf_s, surf_e
	260	k = surf_lsm_h(1)%k(m)
	261
	262	tend(k,j,i) = tend(k,j,i) &
	263	+ ( - surf_lsm_h(1)%vsws(m) ) * ddzw(k) * drho_air(k)
	264
	265	ENDDO
	266	!
[2232]	267	!-- Urban-type surfaces, upward-facing
[4671]	268	surf_s = surf_usm_h(0)%start_index(j,i)
	269	surf_e = surf_usm_h(0)%end_index(j,i)
[2232]	270	DO m = surf_s, surf_e
[4671]	271	k = surf_usm_h(0)%k(m)
[102]	272
[4583]	273	tend(k,j,i) = tend(k,j,i) &
[4671]	274	+ ( - ( - surf_usm_h(0)%vsws(m) ) ) * ddzw(k) * drho_air(k)
[2232]	275
	276	ENDDO
[4671]	277	!
	278	!-- Urban-type surfaces, downward-facing
	279	surf_s = surf_usm_h(1)%start_index(j,i)
	280	surf_e = surf_usm_h(1)%end_index(j,i)
	281	DO m = surf_s, surf_e
	282	k = surf_usm_h(1)%k(m)
	283
	284	tend(k,j,i) = tend(k,j,i) &
	285	+ ( - surf_usm_h(1)%vsws(m) ) * ddzw(k) * drho_air(k)
	286
	287	ENDDO
[102]	288	ENDIF
[2232]	289	!
	290	!-- Add momentum flux at model top
[2638]	291	IF ( use_top_fluxes .AND. constant_top_momentumflux ) THEN
[2232]	292	surf_s = surf_def_h(2)%start_index(j,i)
	293	surf_e = surf_def_h(2)%end_index(j,i)
	294	DO m = surf_s, surf_e
[102]	295
[2232]	296	k = surf_def_h(2)%k(m)
	297
[4583]	298	tend(k,j,i) = tend(k,j,i) &
	299	+ ( - surf_def_h(2)%vsws(m) ) * ddzw(k) * drho_air(k)
[2232]	300	ENDDO
	301	ENDIF
	302
[1]	303	ENDDO
	304	ENDDO
	305
	306	END SUBROUTINE diffusion_v
	307
	308
[4583]	309	!--------------------------------------------------------------------------------------------------!
[1682]	310	! Description:
	311	! ------------
	312	!> Call for grid point i,j
[4583]	313	!--------------------------------------------------------------------------------------------------!
[1001]	314	SUBROUTINE diffusion_v_ij( i, j )
[1]	315
[4583]	316	USE arrays_3d, &
	317	ONLY: ddzu, ddzw, drho_air, km, tend, u, v, w, rho_air_zw
	318
	319	USE control_parameters, &
	320	ONLY: constant_top_momentumflux, use_surface_fluxes, use_top_fluxes
	321
	322	USE grid_variables, &
[2232]	323	ONLY: ddx, ddy, ddy2
[4583]	324
	325	USE indices, &
[4346]	326	ONLY: nzb, nzt, wall_flags_total_0
[4583]	327
[1320]	328	USE kinds
[1]	329
[4583]	330	USE surface_mod, &
	331	ONLY : surf_def_h, surf_def_v, surf_lsm_h, surf_lsm_v, surf_usm_h, surf_usm_v
[2232]	332
[1]	333	IMPLICIT NONE
	334
	335
[2232]	336	INTEGER(iwp) :: i !< running index x direction
	337	INTEGER(iwp) :: j !< running index y direction
	338	INTEGER(iwp) :: k !< running index z direction
	339	INTEGER(iwp) :: l !< running index of surface type, south- or north-facing wall
	340	INTEGER(iwp) :: m !< running index surface elements
	341	INTEGER(iwp) :: surf_e !< End index of surface elements at (j,i)-gridpoint
	342	INTEGER(iwp) :: surf_s !< Start index of surface elements at (j,i)-gridpoint
[1001]	343
[2232]	344	REAL(wp) :: flag !< flag to mask topography grid points
[3547]	345	REAL(wp) :: kmxm !< diffusion coefficient on leftward side of the v-gridbox - interpolated onto xu-yv grid
	346	REAL(wp) :: kmxp !< diffusion coefficient on rightward side of the v-gridbox - interpolated onto xu-yv grid
	347	REAL(wp) :: kmzm !< diffusion coefficient on bottom of the gridbox - interpolated onto xu-zw grid
	348	REAL(wp) :: kmzp !< diffusion coefficient on top of the gridbox - interpolated onto xu-zw grid
[4583]	349	REAL(wp) :: mask_bottom !< flag to mask vertical upward-facing surface
	350	REAL(wp) :: mask_east !< flag to mask vertical surface south of the grid point
	351	REAL(wp) :: mask_west !< flag to mask vertical surface north of the grid point
[2232]	352	REAL(wp) :: mask_top !< flag to mask vertical downward-facing surface
	353
[1]	354	!
	355	!-- Compute horizontal diffusion
[2232]	356	DO k = nzb+1, nzt
[1]	357	!
[4583]	358	!-- Predetermine flag to mask topography and wall-bounded grid points.
	359	!-- It is sufficient to masked only east- and west-facing surfaces, which need special
	360	!-- treatment for the v-component.
	361	flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2 ) )
[4346]	362	mask_east = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i+1), 2 ) )
	363	mask_west = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i-1), 2 ) )
[2232]	364	!
[1]	365	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	366	kmxp = 0.25_wp * ( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
	367	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]	368
[4583]	369	tend(k,j,i) = tend(k,j,i) + ( &
	370	mask_east * kmxp * ( &
	371	( v(k,j,i+1) - v(k,j,i) ) * ddx &
	372	+ ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
	373	) &
	374	- mask_west * kmxm * ( &
	375	( v(k,j,i) - v(k,j,i-1) ) * ddx &
	376	+ ( u(k,j,i) - u(k,j-1,i) ) * ddy &
	377	) &
	378	) * ddx * flag &
	379	+ 2.0_wp * ( &
	380	km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
	381	- km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
[2232]	382	) * ddy2 * flag
[1]	383	ENDDO
	384
	385	!
[4583]	386	!-- Add horizontal momentum flux v'u' at east- (l=2) and west-facing (l=3) surfaces. Note, in the
	387	!-- the flat case, loops won't be entered as start_index > end_index. Furtermore, note, no
	388	!-- vertical natural surfaces so far.
[2232]	389	!-- Default-type surfaces
	390	DO l = 2, 3
	391	surf_s = surf_def_v(l)%start_index(j,i)
	392	surf_e = surf_def_v(l)%end_index(j,i)
	393	DO m = surf_s, surf_e
	394	k = surf_def_v(l)%k(m)
	395	tend(k,j,i) = tend(k,j,i) + surf_def_v(l)%mom_flux_uv(m) * ddx
[4583]	396	ENDDO
[2232]	397	ENDDO
[51]	398	!
[2232]	399	!-- Natural-type surfaces
	400	DO l = 2, 3
	401	surf_s = surf_lsm_v(l)%start_index(j,i)
	402	surf_e = surf_lsm_v(l)%end_index(j,i)
	403	DO m = surf_s, surf_e
	404	k = surf_lsm_v(l)%k(m)
	405	tend(k,j,i) = tend(k,j,i) + surf_lsm_v(l)%mom_flux_uv(m) * ddx
[4583]	406	ENDDO
[2232]	407	ENDDO
[1]	408	!
[2232]	409	!-- Urban-type surfaces
	410	DO l = 2, 3
	411	surf_s = surf_usm_v(l)%start_index(j,i)
	412	surf_e = surf_usm_v(l)%end_index(j,i)
	413	DO m = surf_s, surf_e
	414	k = surf_usm_v(l)%k(m)
	415	tend(k,j,i) = tend(k,j,i) + surf_usm_v(l)%mom_flux_uv(m) * ddx
[4583]	416	ENDDO
[2232]	417	ENDDO
[1]	418	!
[4583]	419	!-- Compute vertical diffusion. In case of simulating a surface layer, respective grid diffusive
	420	!-- fluxes are masked (flag 8) within this loop, and added further below, else, simple gradient
	421	!-- approach is applied. Model top is also mask if top-momentum flux is given.
[2232]	422	DO k = nzb+1, nzt
	423	!
[4583]	424	!-- Determine flags to mask topography below and above. Flag 2 is used to mask topography in
	425	!-- general, while flag 10 implies also information about use_surface_fluxes. Flag 9 is used
	426	!-- to control momentum flux at model top.
	427	mask_bottom = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k-1,j,i), 8 ) )
	428	mask_top = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k+1,j,i), 8 ) ) * &
	429	MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k+1,j,i), 9 ) )
	430	flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2 ) )
[2232]	431	!
[1]	432	!-- Interpolate eddy diffusivities on staggered gridpoints
[1340]	433	kmzp = 0.25_wp * ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
	434	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]	435
[4583]	436	tend(k,j,i) = tend(k,j,i) &
	437	& + ( kmzp * ( ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
	438	& + ( w(k,j,i) - w(k,j-1,i) ) * ddy &
	439	& ) * rho_air_zw(k) * mask_top &
	440	& - kmzm * ( ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
	441	& + ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
	442	& ) * rho_air_zw(k-1) * mask_bottom &
	443	& ) * ddzw(k) * drho_air(k) * flag
[1]	444	ENDDO
	445
	446	!
[4583]	447	!-- Vertical diffusion at the first grid point above the surface, if the momentum flux at the
	448	!-- bottom is given by the Prandtl law or if it is prescribed by the user.
	449	!-- Difference quotient of the momentum flux is not formed over half of the grid spacing
	450	!-- (2.0*ddzw(k)) any more, since the comparison with other (LES) models showed that the values
	451	!-- of the momentum flux becomes too large in this case.
[1]	452	IF ( use_surface_fluxes ) THEN
	453	!
[2232]	454	!-- Default-type surfaces, upward-facing
	455	surf_s = surf_def_h(0)%start_index(j,i)
	456	surf_e = surf_def_h(0)%end_index(j,i)
	457	DO m = surf_s, surf_e
	458	k = surf_def_h(0)%k(m)
[1]	459
[4583]	460	tend(k,j,i) = tend(k,j,i) + ( - ( - surf_def_h(0)%vsws(m) ) ) * ddzw(k) * drho_air(k)
[2232]	461	ENDDO
	462	!
	463	!-- Default-type surfaces, dowward-facing
	464	surf_s = surf_def_h(1)%start_index(j,i)
	465	surf_e = surf_def_h(1)%end_index(j,i)
	466	DO m = surf_s, surf_e
	467	k = surf_def_h(1)%k(m)
[1]	468
[4583]	469	tend(k,j,i) = tend(k,j,i) + ( - surf_def_h(1)%vsws(m) ) * ddzw(k) * drho_air(k)
[2232]	470	ENDDO
[102]	471	!
[2232]	472	!-- Natural-type surfaces, upward-facing
[4671]	473	surf_s = surf_lsm_h(0)%start_index(j,i)
	474	surf_e = surf_lsm_h(0)%end_index(j,i)
[2232]	475	DO m = surf_s, surf_e
[4671]	476	k = surf_lsm_h(0)%k(m)
[2232]	477
[4671]	478	tend(k,j,i) = tend(k,j,i) + ( - ( - surf_lsm_h(0)%vsws(m) ) ) * ddzw(k) * drho_air(k)
[2232]	479
	480	ENDDO
[102]	481	!
[4671]	482	!-- Natural-type surfaces, downward-facing
	483	surf_s = surf_lsm_h(1)%start_index(j,i)
	484	surf_e = surf_lsm_h(1)%end_index(j,i)
	485	DO m = surf_s, surf_e
	486	k = surf_lsm_h(1)%k(m)
	487
	488	tend(k,j,i) = tend(k,j,i) + ( - surf_lsm_h(1)%vsws(m) ) * ddzw(k) * drho_air(k)
	489
	490	ENDDO
	491	!
[2232]	492	!-- Urban-type surfaces, upward-facing
[4671]	493	surf_s = surf_usm_h(0)%start_index(j,i)
	494	surf_e = surf_usm_h(0)%end_index(j,i)
[2232]	495	DO m = surf_s, surf_e
[4671]	496	k = surf_usm_h(0)%k(m)
[102]	497
[4671]	498	tend(k,j,i) = tend(k,j,i) + ( - ( - surf_usm_h(0)%vsws(m) ) ) * ddzw(k) * drho_air(k)
[2232]	499
	500	ENDDO
[4671]	501	!
	502	!-- Urban-type surfaces, downward-facing
	503	surf_s = surf_usm_h(1)%start_index(j,i)
	504	surf_e = surf_usm_h(1)%end_index(j,i)
	505	DO m = surf_s, surf_e
	506	k = surf_usm_h(1)%k(m)
	507
	508	tend(k,j,i) = tend(k,j,i) + ( - surf_usm_h(1)%vsws(m) ) * ddzw(k) * drho_air(k)
	509
	510	ENDDO
[102]	511	ENDIF
[2232]	512	!
	513	!-- Add momentum flux at model top
[2638]	514	IF ( use_top_fluxes .AND. constant_top_momentumflux ) THEN
[2232]	515	surf_s = surf_def_h(2)%start_index(j,i)
	516	surf_e = surf_def_h(2)%end_index(j,i)
	517	DO m = surf_s, surf_e
[102]	518
[2232]	519	k = surf_def_h(2)%k(m)
	520
[4583]	521	tend(k,j,i) = tend(k,j,i) + ( - surf_def_h(2)%vsws(m) ) * ddzw(k) * drho_air(k)
[2232]	522	ENDDO
	523	ENDIF
	524
[1]	525	END SUBROUTINE diffusion_v_ij
	526
[1321]	527	END MODULE diffusion_v_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |