Home

Context Navigation

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

Last change on this file since 1029 was 1017, checked in by raasch, 12 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 16.9 KB

Line
1	MODULE diffusion_s_mod
2
3	!------------------------------------------------------------------------------!
4	! Current revisions:
5	! ------------------
6	!
7	!
8	! Former revisions:
9	! -----------------
10	! $Id: diffusion_s.f90 1017 2012-09-27 11:28:50Z maronga $
11	!
12	! 1015 2012-09-27 09:23:24Z raasch
13	! accelerator version (*_acc) added
14	!
15	! 1010 2012-09-20 07:59:54Z raasch
16	! cpp switch __nopointer added for pointer free version
17	!
18	! 1001 2012-09-13 14:08:46Z raasch
19	! some arrays comunicated by module instead of parameter list
20	!
21	! 667 2010-12-23 12:06:00Z suehring/gryschka
22	! nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng
23	!
24	! 183 2008-08-04 15:39:12Z letzel
25	! bugfix: calculation of fluxes at vertical surfaces
26	!
27	! 129 2007-10-30 12:12:24Z letzel
28	! replace wall_heatflux by wall_s_flux that is now included in the parameter
29	! list, bugfix for assignment of fluxes at walls
30	!
31	! 20 2007-02-26 00:12:32Z raasch
32	! Bugfix: ddzw dimensioned 1:nzt"+1"
33	! Calculation extended for gridpoint nzt, fluxes can be given at top,
34	! +s_flux_t in parameter list, s_flux renamed s_flux_b
35	!
36	! RCS Log replace by Id keyword, revision history cleaned up
37	!
38	! Revision 1.8 2006/02/23 10:34:17 raasch
39	! nzb_2d replaced by nzb_s_outer in horizontal diffusion and by nzb_s_inner
40	! or nzb_diff_s_inner, respectively, in vertical diffusion, prescribed surface
41	! fluxes at vertically oriented topography
42	!
43	! Revision 1.1 2000/04/13 14:54:02 schroeter
44	! Initial revision
45	!
46	!
47	! Description:
48	! ------------
49	! Diffusion term of scalar quantities (temperature and water content)
50	!------------------------------------------------------------------------------!
51
52	PRIVATE
53	PUBLIC diffusion_s, diffusion_s_acc
54
55	INTERFACE diffusion_s
56	MODULE PROCEDURE diffusion_s
57	MODULE PROCEDURE diffusion_s_ij
58	END INTERFACE diffusion_s
59
60	INTERFACE diffusion_s_acc
61	MODULE PROCEDURE diffusion_s_acc
62	END INTERFACE diffusion_s_acc
63
64	CONTAINS
65
66
67	!------------------------------------------------------------------------------!
68	! Call for all grid points
69	!------------------------------------------------------------------------------!
70	SUBROUTINE diffusion_s( s, s_flux_b, s_flux_t, wall_s_flux )
71
72	USE arrays_3d
73	USE control_parameters
74	USE grid_variables
75	USE indices
76
77	IMPLICIT NONE
78
79	INTEGER :: i, j, k
80	REAL :: vertical_gridspace
81	REAL :: wall_s_flux(0:4)
82	REAL, DIMENSION(nysg:nyng,nxlg:nxrg) :: s_flux_b, s_flux_t
83	#if defined( __nopointer )
84	REAL, DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: s
85	#else
86	REAL, DIMENSION(:,:,:), POINTER :: s
87	#endif
88
89	DO i = nxl, nxr
90	DO j = nys,nyn
91	!
92	!-- Compute horizontal diffusion
93	DO k = nzb_s_outer(j,i)+1, nzt
94
95	tend(k,j,i) = tend(k,j,i) &
96	+ 0.5 * ( &
97	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
98	- ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
99	) * ddx2 &
100	+ 0.5 * ( &
101	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
102	- ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
103	) * ddy2
104	ENDDO
105
106	!
107	!-- Apply prescribed horizontal wall heatflux where necessary
108	IF ( ( wall_w_x(j,i) .NE. 0.0 ) .OR. ( wall_w_y(j,i) .NE. 0.0 ) ) &
109	THEN
110	DO k = nzb_s_inner(j,i)+1, nzb_s_outer(j,i)
111
112	tend(k,j,i) = tend(k,j,i) &
113	+ ( fwxp(j,i) * 0.5 * &
114	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
115	+ ( 1.0 - fwxp(j,i) ) * wall_s_flux(1) &
116	-fwxm(j,i) * 0.5 * &
117	( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
118	+ ( 1.0 - fwxm(j,i) ) * wall_s_flux(2) &
119	) * ddx2 &
120	+ ( fwyp(j,i) * 0.5 * &
121	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
122	+ ( 1.0 - fwyp(j,i) ) * wall_s_flux(3) &
123	-fwym(j,i) * 0.5 * &
124	( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
125	+ ( 1.0 - fwym(j,i) ) * wall_s_flux(4) &
126	) * ddy2
127	ENDDO
128	ENDIF
129
130	!
131	!-- Compute vertical diffusion. In case that surface fluxes have been
132	!-- prescribed or computed at bottom and/or top, index k starts/ends at
133	!-- nzb+2 or nzt-1, respectively.
134	DO k = nzb_diff_s_inner(j,i), nzt_diff
135
136	tend(k,j,i) = tend(k,j,i) &
137	+ 0.5 * ( &
138	( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) &
139	- ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) &
140	) * ddzw(k)
141	ENDDO
142
143	!
144	!-- Vertical diffusion at the first computational gridpoint along
145	!-- z-direction
146	IF ( use_surface_fluxes ) THEN
147
148	k = nzb_s_inner(j,i)+1
149
150	tend(k,j,i) = tend(k,j,i) &
151	+ ( 0.5 * ( kh(k,j,i)+kh(k+1,j,i) ) &
152	* ( s(k+1,j,i)-s(k,j,i) ) &
153	* ddzu(k+1) &
154	+ s_flux_b(j,i) &
155	) * ddzw(k)
156
157	ENDIF
158
159	!
160	!-- Vertical diffusion at the last computational gridpoint along
161	!-- z-direction
162	IF ( use_top_fluxes ) THEN
163
164	k = nzt
165
166	tend(k,j,i) = tend(k,j,i) &
167	+ ( - s_flux_t(j,i) &
168	- 0.5 * ( kh(k-1,j,i)+kh(k,j,i) ) &
169	* ( s(k,j,i)-s(k-1,j,i) ) &
170	* ddzu(k) &
171	) * ddzw(k)
172
173	ENDIF
174
175	ENDDO
176	ENDDO
177
178	END SUBROUTINE diffusion_s
179
180
181	!------------------------------------------------------------------------------!
182	! Call for all grid points - accelerator version
183	!------------------------------------------------------------------------------!
184	SUBROUTINE diffusion_s_acc( s, s_flux_b, s_flux_t, wall_s_flux )
185
186	USE arrays_3d
187	USE control_parameters
188	USE grid_variables
189	USE indices
190
191	IMPLICIT NONE
192
193	INTEGER :: i, j, k
194	REAL :: vertical_gridspace
195	REAL :: wall_s_flux(0:4)
196	REAL, DIMENSION(nysg:nyng,nxlg:nxrg) :: s_flux_b, s_flux_t
197	#if defined( __nopointer )
198	REAL, DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: s
199	#else
200	REAL, DIMENSION(:,:,:), POINTER :: s
201	#endif
202
203	!$acc kernels present( ddzu, ddzw, fwxm, fwxp, fwym, fwyp, kh ) &
204	!$acc present( nzb_diff_s_inner, nzb_s_inner, nzb_s_outer, s ) &
205	!$acc present( s_flux_b, s_flux_t, tend, wall_s_flux ) &
206	!$acc present( wall_w_x, wall_w_y )
207	!$acc loop
208	DO i = nxl, nxr
209	DO j = nys,nyn
210	!
211	!-- Compute horizontal diffusion
212	!$acc loop vector( 32 )
213	DO k = 1, nzt
214	IF ( k > nzb_s_outer(j,i) ) THEN
215
216	tend(k,j,i) = tend(k,j,i) &
217	+ 0.5 * ( &
218	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
219	- ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
220	) * ddx2 &
221	+ 0.5 * ( &
222	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
223	- ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
224	) * ddy2
225	ENDIF
226	ENDDO
227
228	!
229	!-- Apply prescribed horizontal wall heatflux where necessary
230	!$acc loop vector(32)
231	DO k = 1, nzt
232	IF ( k > nzb_s_inner(j,i) .AND. k <= nzb_s_outer(j,i) .AND. &
233	( wall_w_x(j,i) /= 0.0 .OR. wall_w_y(j,i) /= 0.0 ) ) &
234	THEN
235	tend(k,j,i) = tend(k,j,i) &
236	+ ( fwxp(j,i) * 0.5 * &
237	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
238	+ ( 1.0 - fwxp(j,i) ) * wall_s_flux(1) &
239	-fwxm(j,i) * 0.5 * &
240	( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
241	+ ( 1.0 - fwxm(j,i) ) * wall_s_flux(2) &
242	) * ddx2 &
243	+ ( fwyp(j,i) * 0.5 * &
244	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
245	+ ( 1.0 - fwyp(j,i) ) * wall_s_flux(3) &
246	-fwym(j,i) * 0.5 * &
247	( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
248	+ ( 1.0 - fwym(j,i) ) * wall_s_flux(4) &
249	) * ddy2
250	ENDIF
251	ENDDO
252
253	!
254	!-- Compute vertical diffusion. In case that surface fluxes have been
255	!-- prescribed or computed at bottom and/or top, index k starts/ends at
256	!-- nzb+2 or nzt-1, respectively.
257	!$acc loop vector( 32 )
258	DO k = 1, nzt_diff
259	IF ( k >= nzb_diff_s_inner(j,i) ) THEN
260	tend(k,j,i) = tend(k,j,i) &
261	+ 0.5 * ( &
262	( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) &
263	- ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) &
264	) * ddzw(k)
265	ENDIF
266	ENDDO
267
268	!
269	!-- Vertical diffusion at the first computational gridpoint along
270	!-- z-direction
271	!$acc loop vector( 32 )
272	DO k = 1, nzt
273	IF ( use_surface_fluxes .AND. k == nzb_s_inner(j,i)+1 ) THEN
274	tend(k,j,i) = tend(k,j,i) &
275	+ ( 0.5 * ( kh(k,j,i)+kh(k+1,j,i) ) &
276	* ( s(k+1,j,i)-s(k,j,i) ) &
277	* ddzu(k+1) &
278	+ s_flux_b(j,i) &
279	) * ddzw(k)
280	ENDIF
281
282	!
283	!-- Vertical diffusion at the last computational gridpoint along
284	!-- z-direction
285	IF ( use_top_fluxes .AND. k == nzt ) THEN
286	tend(k,j,i) = tend(k,j,i) &
287	+ ( - s_flux_t(j,i) &
288	- 0.5 * ( kh(k-1,j,i)+kh(k,j,i) )&
289	* ( s(k,j,i)-s(k-1,j,i) ) &
290	* ddzu(k) &
291	) * ddzw(k)
292	ENDIF
293	ENDDO
294
295	ENDDO
296	ENDDO
297	!$acc end kernels
298
299	END SUBROUTINE diffusion_s_acc
300
301
302	!------------------------------------------------------------------------------!
303	! Call for grid point i,j
304	!------------------------------------------------------------------------------!
305	SUBROUTINE diffusion_s_ij( i, j, s, s_flux_b, s_flux_t, wall_s_flux )
306
307	USE arrays_3d
308	USE control_parameters
309	USE grid_variables
310	USE indices
311
312	IMPLICIT NONE
313
314	INTEGER :: i, j, k
315	REAL :: vertical_gridspace
316	REAL :: wall_s_flux(0:4)
317	REAL, DIMENSION(nysg:nyng,nxlg:nxrg) :: s_flux_b, s_flux_t
318	#if defined( __nopointer )
319	REAL, DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: s
320	#else
321	REAL, DIMENSION(:,:,:), POINTER :: s
322	#endif
323
324	!
325	!-- Compute horizontal diffusion
326	DO k = nzb_s_outer(j,i)+1, nzt
327
328	tend(k,j,i) = tend(k,j,i) &
329	+ 0.5 * ( &
330	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
331	- ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
332	) * ddx2 &
333	+ 0.5 * ( &
334	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
335	- ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
336	) * ddy2
337	ENDDO
338
339	!
340	!-- Apply prescribed horizontal wall heatflux where necessary
341	IF ( ( wall_w_x(j,i) .NE. 0.0 ) .OR. ( wall_w_y(j,i) .NE. 0.0 ) ) &
342	THEN
343	DO k = nzb_s_inner(j,i)+1, nzb_s_outer(j,i)
344
345	tend(k,j,i) = tend(k,j,i) &
346	+ ( fwxp(j,i) * 0.5 * &
347	( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) &
348	+ ( 1.0 - fwxp(j,i) ) * wall_s_flux(1) &
349	-fwxm(j,i) * 0.5 * &
350	( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) &
351	+ ( 1.0 - fwxm(j,i) ) * wall_s_flux(2) &
352	) * ddx2 &
353	+ ( fwyp(j,i) * 0.5 * &
354	( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) &
355	+ ( 1.0 - fwyp(j,i) ) * wall_s_flux(3) &
356	-fwym(j,i) * 0.5 * &
357	( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) &
358	+ ( 1.0 - fwym(j,i) ) * wall_s_flux(4) &
359	) * ddy2
360	ENDDO
361	ENDIF
362
363	!
364	!-- Compute vertical diffusion. In case that surface fluxes have been
365	!-- prescribed or computed at bottom and/or top, index k starts/ends at
366	!-- nzb+2 or nzt-1, respectively.
367	DO k = nzb_diff_s_inner(j,i), nzt_diff
368
369	tend(k,j,i) = tend(k,j,i) &
370	+ 0.5 * ( &
371	( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) &
372	- ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) &
373	) * ddzw(k)
374	ENDDO
375
376	!
377	!-- Vertical diffusion at the first computational gridpoint along z-direction
378	IF ( use_surface_fluxes ) THEN
379
380	k = nzb_s_inner(j,i)+1
381
382	tend(k,j,i) = tend(k,j,i) + ( 0.5 * ( kh(k,j,i)+kh(k+1,j,i) ) &
383	* ( s(k+1,j,i)-s(k,j,i) ) &
384	* ddzu(k+1) &
385	+ s_flux_b(j,i) &
386	) * ddzw(k)
387
388	ENDIF
389
390	!
391	!-- Vertical diffusion at the last computational gridpoint along z-direction
392	IF ( use_top_fluxes ) THEN
393
394	k = nzt
395
396	tend(k,j,i) = tend(k,j,i) + ( - s_flux_t(j,i) &
397	- 0.5 * ( kh(k-1,j,i)+kh(k,j,i) ) &
398	* ( s(k,j,i)-s(k-1,j,i) ) &
399	* ddzu(k) &
400	) * ddzw(k)
401
402	ENDIF
403
404	END SUBROUTINE diffusion_s_ij
405
406	END MODULE diffusion_s_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |