Home

Context Navigation

source: palm/trunk/SOURCE/diffusion_s.f90 @ 1520

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

Rev	Line
[1]	1	MODULE diffusion_s_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	!
[1310]	17	! Copyright 1997-2014 Leibniz Universitaet Hannover
[1036]	18	!--------------------------------------------------------------------------------!
	19	!
[484]	20	! Current revisions:
[1001]	21	! ------------------
[1341]	22	!
[1375]	23	!
[1321]	24	! Former revisions:
	25	! -----------------
	26	! $Id: diffusion_s.f90 1375 2014-04-25 13:07:08Z hoffmann $
	27	!
[1375]	28	! 1374 2014-04-25 12:55:07Z raasch
	29	! missing variables added to ONLY list
	30	!
[1341]	31	! 1340 2014-03-25 19:45:13Z kanani
	32	! REAL constants defined as wp-kind
	33	!
[1321]	34	! 1320 2014-03-20 08:40:49Z raasch
[1320]	35	! ONLY-attribute added to USE-statements,
	36	! kind-parameters added to all INTEGER and REAL declaration statements,
	37	! kinds are defined in new module kinds,
	38	! revision history before 2012 removed,
	39	! comment fields (!:) to be used for variable explanations added to
	40	! all variable declaration statements
[1321]	41	!
[1258]	42	! 1257 2013-11-08 15:18:40Z raasch
	43	! openacc loop and loop vector clauses removed
	44	!
[1132]	45	! 1128 2013-04-12 06:19:32Z raasch
	46	! loop index bounds in accelerator version replaced by i_left, i_right, j_south,
	47	! j_north
	48	!
[1093]	49	! 1092 2013-02-02 11:24:22Z raasch
	50	! unused variables removed
	51	!
[1037]	52	! 1036 2012-10-22 13:43:42Z raasch
	53	! code put under GPL (PALM 3.9)
	54	!
[1017]	55	! 1015 2012-09-27 09:23:24Z raasch
	56	! accelerator version (*_acc) added
	57	!
[1011]	58	! 1010 2012-09-20 07:59:54Z raasch
	59	! cpp switch __nopointer added for pointer free version
	60	!
[1002]	61	! 1001 2012-09-13 14:08:46Z raasch
	62	! some arrays comunicated by module instead of parameter list
	63	!
[1]	64	! Revision 1.1 2000/04/13 14:54:02 schroeter
	65	! Initial revision
	66	!
	67	!
	68	! Description:
	69	! ------------
	70	! Diffusion term of scalar quantities (temperature and water content)
	71	!------------------------------------------------------------------------------!
	72
	73	PRIVATE
[1015]	74	PUBLIC diffusion_s, diffusion_s_acc
[1]	75
	76	INTERFACE diffusion_s
	77	MODULE PROCEDURE diffusion_s
	78	MODULE PROCEDURE diffusion_s_ij
	79	END INTERFACE diffusion_s
	80
[1015]	81	INTERFACE diffusion_s_acc
	82	MODULE PROCEDURE diffusion_s_acc
	83	END INTERFACE diffusion_s_acc
	84
[1]	85	CONTAINS
	86
	87
	88	!------------------------------------------------------------------------------!
	89	! Call for all grid points
	90	!------------------------------------------------------------------------------!
[1001]	91	SUBROUTINE diffusion_s( s, s_flux_b, s_flux_t, wall_s_flux )
[1]	92
[1320]	93	USE arrays_3d, &
	94	ONLY: ddzu, ddzw, kh, tend
	95
	96	USE control_parameters, &
	97	ONLY: use_surface_fluxes, use_top_fluxes
	98
	99	USE grid_variables, &
	100	ONLY: ddx2, ddy2, fwxm, fwxp, fwym, fwyp, wall_w_x, wall_w_y
	101
	102	USE indices, &
[1374]	103	ONLY: nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, &
[1320]	104	nzb_diff_s_inner, nzb_s_inner, nzb_s_outer, nzt, nzt_diff
	105
	106	USE kinds
[1]	107
	108	IMPLICIT NONE
	109
[1320]	110	INTEGER(iwp) :: i !:
	111	INTEGER(iwp) :: j !:
	112	INTEGER(iwp) :: k !:
	113	REAL(wp) :: wall_s_flux(0:4) !:
	114	REAL(wp), DIMENSION(nysg:nyng,nxlg:nxrg) :: s_flux_b, s_flux_t !:
[1010]	115	#if defined( __nopointer )
[1320]	116	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: s !:
[1010]	117	#else
[1320]	118	REAL(wp), DIMENSION(:,:,:), POINTER :: s !:
[1010]	119	#endif
[1]	120
	121	DO i = nxl, nxr
	122	DO j = nys,nyn
	123	!
	124	!-- Compute horizontal diffusion
[19]	125	DO k = nzb_s_outer(j,i)+1, nzt
[1]	126
[1320]	127	tend(k,j,i) = tend(k,j,i) &
[1340]	128	+ 0.5_wp * ( &
[1320]	129	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
	130	- ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
[1340]	131	) * ddx2 &
	132	+ 0.5_wp * ( &
[1320]	133	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
	134	- ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
[1340]	135	) * ddy2
[1]	136	ENDDO
	137
	138	!
	139	!-- Apply prescribed horizontal wall heatflux where necessary
[1340]	140	IF ( ( wall_w_x(j,i) .NE. 0.0_wp ) .OR. ( wall_w_y(j,i) .NE. 0.0_wp ) ) &
[1]	141	THEN
	142	DO k = nzb_s_inner(j,i)+1, nzb_s_outer(j,i)
	143
[1320]	144	tend(k,j,i) = tend(k,j,i) &
[1340]	145	+ ( fwxp(j,i) * 0.5_wp * &
[1320]	146	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
[1340]	147	+ ( 1.0_wp - fwxp(j,i) ) * wall_s_flux(1) &
	148	-fwxm(j,i) * 0.5_wp * &
[1320]	149	( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
[1340]	150	+ ( 1.0_wp - fwxm(j,i) ) * wall_s_flux(2) &
[1320]	151	) * ddx2 &
[1340]	152	+ ( fwyp(j,i) * 0.5_wp * &
[1320]	153	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
[1340]	154	+ ( 1.0_wp - fwyp(j,i) ) * wall_s_flux(3) &
	155	-fwym(j,i) * 0.5_wp * &
[1320]	156	( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
[1340]	157	+ ( 1.0_wp - fwym(j,i) ) * wall_s_flux(4) &
[1]	158	) * ddy2
	159	ENDDO
	160	ENDIF
	161
	162	!
	163	!-- Compute vertical diffusion. In case that surface fluxes have been
[19]	164	!-- prescribed or computed at bottom and/or top, index k starts/ends at
	165	!-- nzb+2 or nzt-1, respectively.
	166	DO k = nzb_diff_s_inner(j,i), nzt_diff
[1]	167
[1320]	168	tend(k,j,i) = tend(k,j,i) &
[1340]	169	+ 0.5_wp * ( &
[1320]	170	( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) &
	171	- ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) &
[1340]	172	) * ddzw(k)
[1]	173	ENDDO
	174
	175	!
[19]	176	!-- Vertical diffusion at the first computational gridpoint along
[1]	177	!-- z-direction
	178	IF ( use_surface_fluxes ) THEN
	179
	180	k = nzb_s_inner(j,i)+1
	181
[1320]	182	tend(k,j,i) = tend(k,j,i) &
[1340]	183	+ ( 0.5_wp * ( kh(k,j,i)+kh(k+1,j,i) ) &
	184	* ( s(k+1,j,i)-s(k,j,i) ) &
	185	* ddzu(k+1) &
[1320]	186	+ s_flux_b(j,i) &
[1]	187	) * ddzw(k)
	188
	189	ENDIF
	190
[19]	191	!
	192	!-- Vertical diffusion at the last computational gridpoint along
	193	!-- z-direction
	194	IF ( use_top_fluxes ) THEN
	195
	196	k = nzt
	197
[1320]	198	tend(k,j,i) = tend(k,j,i) &
	199	+ ( - s_flux_t(j,i) &
[1340]	200	- 0.5_wp * ( kh(k-1,j,i)+kh(k,j,i) )&
	201	* ( s(k,j,i)-s(k-1,j,i) ) &
	202	* ddzu(k) &
[19]	203	) * ddzw(k)
	204
	205	ENDIF
	206
[1]	207	ENDDO
	208	ENDDO
	209
	210	END SUBROUTINE diffusion_s
	211
	212
	213	!------------------------------------------------------------------------------!
[1015]	214	! Call for all grid points - accelerator version
	215	!------------------------------------------------------------------------------!
	216	SUBROUTINE diffusion_s_acc( s, s_flux_b, s_flux_t, wall_s_flux )
	217
[1320]	218	USE arrays_3d, &
	219	ONLY: ddzu, ddzw, kh, tend
	220
	221	USE control_parameters, &
	222	ONLY: use_surface_fluxes, use_top_fluxes
	223
	224	USE grid_variables, &
	225	ONLY: ddx2, ddy2, fwxm, fwxp, fwym, fwyp, wall_w_x, wall_w_y
	226
	227	USE indices, &
	228	ONLY: i_left, i_right, j_north, j_south, nxlg, nxrg, nyng, nysg, &
[1374]	229	nzb, nzb_diff_s_inner, nzb_s_inner, nzb_s_outer, nzt, nzt_diff
[1320]	230
	231	USE kinds
[1015]	232
	233	IMPLICIT NONE
	234
[1320]	235	INTEGER(iwp) :: i !:
	236	INTEGER(iwp) :: j !:
	237	INTEGER(iwp) :: k !:
	238	REAL(wp) :: wall_s_flux(0:4) !:
	239	REAL(wp), DIMENSION(nysg:nyng,nxlg:nxrg) :: s_flux_b !:
	240	REAL(wp), DIMENSION(nysg:nyng,nxlg:nxrg) :: s_flux_t !:
[1015]	241	#if defined( __nopointer )
[1320]	242	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: s !:
[1015]	243	#else
[1320]	244	REAL(wp), DIMENSION(:,:,:), POINTER :: s !:
[1015]	245	#endif
	246
	247	!$acc kernels present( ddzu, ddzw, fwxm, fwxp, fwym, fwyp, kh ) &
	248	!$acc present( nzb_diff_s_inner, nzb_s_inner, nzb_s_outer, s ) &
	249	!$acc present( s_flux_b, s_flux_t, tend, wall_s_flux ) &
	250	!$acc present( wall_w_x, wall_w_y )
[1128]	251	DO i = i_left, i_right
	252	DO j = j_south, j_north
[1015]	253	!
	254	!-- Compute horizontal diffusion
	255	DO k = 1, nzt
	256	IF ( k > nzb_s_outer(j,i) ) THEN
	257
[1320]	258	tend(k,j,i) = tend(k,j,i) &
[1340]	259	+ 0.5_wp * ( &
[1320]	260	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
	261	- ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
[1340]	262	) * ddx2 &
	263	+ 0.5_wp * ( &
[1320]	264	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
	265	- ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
[1340]	266	) * ddy2
[1015]	267	ENDIF
	268	ENDDO
	269
	270	!
	271	!-- Apply prescribed horizontal wall heatflux where necessary
	272	DO k = 1, nzt
	273	IF ( k > nzb_s_inner(j,i) .AND. k <= nzb_s_outer(j,i) .AND. &
[1340]	274	( wall_w_x(j,i) /= 0.0_wp .OR. wall_w_y(j,i) /= 0.0_wp ) ) &
[1015]	275	THEN
[1320]	276	tend(k,j,i) = tend(k,j,i) &
[1340]	277	+ ( fwxp(j,i) * 0.5_wp * &
[1320]	278	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
[1340]	279	+ ( 1.0_wp - fwxp(j,i) ) * wall_s_flux(1) &
	280	-fwxm(j,i) * 0.5_wp * &
[1320]	281	( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
[1340]	282	+ ( 1.0_wp - fwxm(j,i) ) * wall_s_flux(2) &
[1320]	283	) * ddx2 &
[1340]	284	+ ( fwyp(j,i) * 0.5_wp * &
[1320]	285	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
[1340]	286	+ ( 1.0_wp - fwyp(j,i) ) * wall_s_flux(3) &
	287	-fwym(j,i) * 0.5_wp * &
[1320]	288	( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
[1340]	289	+ ( 1.0_wp - fwym(j,i) ) * wall_s_flux(4) &
[1015]	290	) * ddy2
	291	ENDIF
	292	ENDDO
	293
	294	!
	295	!-- Compute vertical diffusion. In case that surface fluxes have been
	296	!-- prescribed or computed at bottom and/or top, index k starts/ends at
	297	!-- nzb+2 or nzt-1, respectively.
	298	DO k = 1, nzt_diff
	299	IF ( k >= nzb_diff_s_inner(j,i) ) THEN
[1320]	300	tend(k,j,i) = tend(k,j,i) &
[1340]	301	+ 0.5_wp * ( &
[1320]	302	( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) &
	303	- ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) &
[1340]	304	) * ddzw(k)
[1015]	305	ENDIF
	306	ENDDO
	307
	308	!
	309	!-- Vertical diffusion at the first computational gridpoint along
	310	!-- z-direction
	311	DO k = 1, nzt
	312	IF ( use_surface_fluxes .AND. k == nzb_s_inner(j,i)+1 ) THEN
[1320]	313	tend(k,j,i) = tend(k,j,i) &
[1340]	314	+ ( 0.5_wp * ( kh(k,j,i)+kh(k+1,j,i) )&
	315	* ( s(k+1,j,i)-s(k,j,i) ) &
	316	* ddzu(k+1) &
[1320]	317	+ s_flux_b(j,i) &
[1015]	318	) * ddzw(k)
	319	ENDIF
	320
	321	!
	322	!-- Vertical diffusion at the last computational gridpoint along
	323	!-- z-direction
	324	IF ( use_top_fluxes .AND. k == nzt ) THEN
	325	tend(k,j,i) = tend(k,j,i) &
	326	+ ( - s_flux_t(j,i) &
[1340]	327	- 0.5_wp * ( kh(k-1,j,i)+kh(k,j,i) )&
	328	* ( s(k,j,i)-s(k-1,j,i) ) &
	329	* ddzu(k) &
[1015]	330	) * ddzw(k)
	331	ENDIF
	332	ENDDO
	333
	334	ENDDO
	335	ENDDO
	336	!$acc end kernels
	337
	338	END SUBROUTINE diffusion_s_acc
	339
	340
	341	!------------------------------------------------------------------------------!
[1]	342	! Call for grid point i,j
	343	!------------------------------------------------------------------------------!
[1001]	344	SUBROUTINE diffusion_s_ij( i, j, s, s_flux_b, s_flux_t, wall_s_flux )
[1]	345
[1320]	346	USE arrays_3d, &
	347	ONLY: ddzu, ddzw, kh, tend
	348
	349	USE control_parameters, &
	350	ONLY: use_surface_fluxes, use_top_fluxes
	351
	352	USE grid_variables, &
	353	ONLY: ddx2, ddy2, fwxm, fwxp, fwym, fwyp, wall_w_x, wall_w_y
	354
	355	USE indices, &
[1374]	356	ONLY: nxlg, nxrg, nyng, nysg, nzb, nzb_diff_s_inner, nzb_s_inner, &
[1320]	357	nzb_s_outer, nzt, nzt_diff
	358
	359	USE kinds
[1]	360
	361	IMPLICIT NONE
	362
[1320]	363	INTEGER(iwp) :: i !:
	364	INTEGER(iwp) :: j !:
	365	INTEGER(iwp) :: k !:
	366	REAL(wp) :: wall_s_flux(0:4) !:
	367	REAL(wp), DIMENSION(nysg:nyng,nxlg:nxrg) :: s_flux_b !:
	368	REAL(wp), DIMENSION(nysg:nyng,nxlg:nxrg) :: s_flux_t !:
[1010]	369	#if defined( __nopointer )
[1320]	370	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: s !:
[1010]	371	#else
[1320]	372	REAL(wp), DIMENSION(:,:,:), POINTER :: s !:
[1010]	373	#endif
[1]	374
	375	!
	376	!-- Compute horizontal diffusion
[19]	377	DO k = nzb_s_outer(j,i)+1, nzt
[1]	378
[1320]	379	tend(k,j,i) = tend(k,j,i) &
[1340]	380	+ 0.5_wp * ( &
[1320]	381	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
	382	- ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
[1340]	383	) * ddx2 &
	384	+ 0.5_wp * ( &
[1320]	385	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
	386	- ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
[1340]	387	) * ddy2
[1]	388	ENDDO
	389
	390	!
	391	!-- Apply prescribed horizontal wall heatflux where necessary
[1340]	392	IF ( ( wall_w_x(j,i) .NE. 0.0_wp ) .OR. ( wall_w_y(j,i) .NE. 0.0_wp ) ) &
[1]	393	THEN
	394	DO k = nzb_s_inner(j,i)+1, nzb_s_outer(j,i)
	395
[1320]	396	tend(k,j,i) = tend(k,j,i) &
[1340]	397	+ ( fwxp(j,i) * 0.5_wp * &
[1320]	398	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
[1340]	399	+ ( 1.0_wp - fwxp(j,i) ) * wall_s_flux(1) &
	400	-fwxm(j,i) * 0.5_wp * &
[1320]	401	( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
[1340]	402	+ ( 1.0_wp - fwxm(j,i) ) * wall_s_flux(2) &
[1320]	403	) * ddx2 &
[1340]	404	+ ( fwyp(j,i) * 0.5_wp * &
[1320]	405	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
[1340]	406	+ ( 1.0_wp - fwyp(j,i) ) * wall_s_flux(3) &
	407	-fwym(j,i) * 0.5_wp * &
[1320]	408	( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
[1340]	409	+ ( 1.0_wp - fwym(j,i) ) * wall_s_flux(4) &
[1]	410	) * ddy2
	411	ENDDO
	412	ENDIF
	413
	414	!
	415	!-- Compute vertical diffusion. In case that surface fluxes have been
[19]	416	!-- prescribed or computed at bottom and/or top, index k starts/ends at
	417	!-- nzb+2 or nzt-1, respectively.
	418	DO k = nzb_diff_s_inner(j,i), nzt_diff
[1]	419
[1320]	420	tend(k,j,i) = tend(k,j,i) &
[1340]	421	+ 0.5_wp * ( &
[1320]	422	( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) &
	423	- ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) &
[1340]	424	) * ddzw(k)
[1]	425	ENDDO
	426
	427	!
[19]	428	!-- Vertical diffusion at the first computational gridpoint along z-direction
[1]	429	IF ( use_surface_fluxes ) THEN
	430
	431	k = nzb_s_inner(j,i)+1
	432
[1340]	433	tend(k,j,i) = tend(k,j,i) + ( 0.5_wp * ( kh(k,j,i)+kh(k+1,j,i) ) &
	434	* ( s(k+1,j,i)-s(k,j,i) ) &
	435	* ddzu(k+1) &
[1320]	436	+ s_flux_b(j,i) &
[19]	437	) * ddzw(k)
[1]	438
	439	ENDIF
	440
[19]	441	!
	442	!-- Vertical diffusion at the last computational gridpoint along z-direction
	443	IF ( use_top_fluxes ) THEN
	444
	445	k = nzt
	446
[1320]	447	tend(k,j,i) = tend(k,j,i) + ( - s_flux_t(j,i) &
[1340]	448	- 0.5_wp * ( kh(k-1,j,i)+kh(k,j,i) ) &
	449	* ( s(k,j,i)-s(k-1,j,i) ) &
	450	* ddzu(k) &
[19]	451	) * ddzw(k)
	452
	453	ENDIF
	454
[1]	455	END SUBROUTINE diffusion_s_ij
	456
	457	END MODULE diffusion_s_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |