Home

Context Navigation

source: palm/trunk/SOURCE/diffusion_u.f90 @ 4185

Last change on this file since 4185 was 4182, checked in by scharf, 5 years ago
corrected "Former revisions" section minor formatting in "Former revisions" section added "Author" section
Property svn:keywords set to `Id`
File size: 24.1 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |