Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |