Home

Context Navigation

source: palm/trunk/SOURCE/diffusion_w.f90 @ 1001

Last change on this file since 1001 was 1001, checked in by raasch, 12 years ago
leapfrog timestep scheme and upstream-spline advection scheme completely removed from the code, reading of dt_fixed from restart file removed
Property svn:keywords set to `Id`
File size: 12.1 KB

Line
1	MODULE diffusion_w_mod
2
3	!------------------------------------------------------------------------------!
4	! Current revisions:
5	! -----------------
6	! arrays comunicated by module instead of parameter list
7	!
8	! Former revisions:
9	! -----------------
10	! $Id: diffusion_w.f90 1001 2012-09-13 14:08:46Z raasch $
11	!
12	! 978 2012-08-09 08:28:32Z fricke
13	! outflow damping layer removed
14	! kmxm_x/_z and kmxp_x/_z change to kmxm and kmxp
15	! kmym_y/_z and kmyp_y/_z change to kmym and kmyp
16	!
17	! 667 2010-12-23 12:06:00Z suehring/gryschka
18	! nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng
19	!
20	! 366 2009-08-25 08:06:27Z raasch
21	! bc_lr/bc_ns replaced by bc_lr_cyc/bc_ns_cyc
22	!
23	! 75 2007-03-22 09:54:05Z raasch
24	! Wall functions now include diabatic conditions, call of routine wall_fluxes,
25	! z0 removed from argument list
26	!
27	! 20 2007-02-26 00:12:32Z raasch
28	! Bugfix: ddzw dimensioned 1:nzt"+1"
29	!
30	! RCS Log replace by Id keyword, revision history cleaned up
31	!
32	! Revision 1.12 2006/02/23 10:38:03 raasch
33	! nzb_2d replaced by nzb_w_outer, wall functions added for all vertical walls,
34	! +z0 in argument list
35	! WARNING: loops containing the MAX function are still not properly vectorized!
36	!
37	! Revision 1.1 1997/09/12 06:24:11 raasch
38	! Initial revision
39	!
40	!
41	! Description:
42	! ------------
43	! Diffusion term of the w-component
44	!------------------------------------------------------------------------------!
45
46	USE wall_fluxes_mod
47
48	PRIVATE
49	PUBLIC diffusion_w
50
51	INTERFACE diffusion_w
52	MODULE PROCEDURE diffusion_w
53	MODULE PROCEDURE diffusion_w_ij
54	END INTERFACE diffusion_w
55
56	CONTAINS
57
58
59	!------------------------------------------------------------------------------!
60	! Call for all grid points
61	!------------------------------------------------------------------------------!
62	SUBROUTINE diffusion_w
63
64	USE arrays_3d
65	USE control_parameters
66	USE grid_variables
67	USE indices
68
69	IMPLICIT NONE
70
71	INTEGER :: i, j, k
72	REAL :: kmxm, kmxp, kmym, kmyp
73
74	REAL, DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: wsus, wsvs
75
76
77	!
78	!-- First calculate horizontal momentum flux w'u' and/or w'v' at vertical
79	!-- walls, if neccessary
80	IF ( topography /= 'flat' ) THEN
81	CALL wall_fluxes( wsus, 0.0, 0.0, 0.0, 1.0, nzb_w_inner, &
82	nzb_w_outer, wall_w_x )
83	CALL wall_fluxes( wsvs, 0.0, 0.0, 1.0, 0.0, nzb_w_inner, &
84	nzb_w_outer, wall_w_y )
85	ENDIF
86
87	DO i = nxl, nxr
88	DO j = nys, nyn
89	DO k = nzb_w_outer(j,i)+1, nzt-1
90	!
91	!-- Interpolate eddy diffusivities on staggered gridpoints
92	kmxp = 0.25 * &
93	( km(k,j,i)+km(k,j,i+1)+km(k+1,j,i)+km(k+1,j,i+1) )
94	kmxm = 0.25 * &
95	( km(k,j,i)+km(k,j,i-1)+km(k+1,j,i)+km(k+1,j,i-1) )
96	kmyp = 0.25 * &
97	( km(k,j,i)+km(k+1,j,i)+km(k,j+1,i)+km(k+1,j+1,i) )
98	kmym = 0.25 * &
99	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
100
101	tend(k,j,i) = tend(k,j,i) &
102	& + ( kmxp * ( w(k,j,i+1) - w(k,j,i) ) * ddx &
103	& + kmxp * ( u(k+1,j,i+1) - u(k,j,i+1) ) * ddzu(k+1) &
104	& - kmxm * ( w(k,j,i) - w(k,j,i-1) ) * ddx &
105	& - kmxm * ( u(k+1,j,i) - u(k,j,i) ) * ddzu(k+1) &
106	& ) * ddx &
107	& + ( kmyp * ( w(k,j+1,i) - w(k,j,i) ) * ddy &
108	& + kmyp * ( v(k+1,j+1,i) - v(k,j+1,i) ) * ddzu(k+1) &
109	& - kmym * ( w(k,j,i) - w(k,j-1,i) ) * ddy &
110	& - kmym * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
111	& ) * ddy &
112	& + 2.0 * ( &
113	& km(k+1,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1) &
114	& - km(k,j,i) * ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) &
115	& ) * ddzu(k+1)
116	ENDDO
117
118	!
119	!-- Wall functions at all vertical walls, where necessary
120	IF ( wall_w_x(j,i) /= 0.0 .OR. wall_w_y(j,i) /= 0.0 ) THEN
121
122	DO k = nzb_w_inner(j,i)+1, nzb_w_outer(j,i)
123	!
124	!-- Interpolate eddy diffusivities on staggered gridpoints
125	kmxp = 0.25 * &
126	( km(k,j,i)+km(k,j,i+1)+km(k+1,j,i)+km(k+1,j,i+1) )
127	kmxm = 0.25 * &
128	( km(k,j,i)+km(k,j,i-1)+km(k+1,j,i)+km(k+1,j,i-1) )
129	kmyp = 0.25 * &
130	( km(k,j,i)+km(k+1,j,i)+km(k,j+1,i)+km(k+1,j+1,i) )
131	kmym = 0.25 * &
132	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
133
134	tend(k,j,i) = tend(k,j,i) &
135	+ ( fwxp(j,i) * ( &
136	kmxp * ( w(k,j,i+1) - w(k,j,i) ) * ddx &
137	+ kmxp * ( u(k+1,j,i+1) - u(k,j,i+1) ) * ddzu(k+1) &
138	) &
139	- fwxm(j,i) * ( &
140	kmxm * ( w(k,j,i) - w(k,j,i-1) ) * ddx &
141	+ kmxm * ( u(k+1,j,i) - u(k,j,i) ) * ddzu(k+1) &
142	) &
143	+ wall_w_x(j,i) * wsus(k,j,i) &
144	) * ddx &
145	+ ( fwyp(j,i) * ( &
146	kmyp * ( w(k,j+1,i) - w(k,j,i) ) * ddy &
147	+ kmyp * ( v(k+1,j+1,i) - v(k,j+1,i) ) * ddzu(k+1) &
148	) &
149	- fwym(j,i) * ( &
150	kmym * ( w(k,j,i) - w(k,j-1,i) ) * ddy &
151	+ kmym * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
152	) &
153	+ wall_w_y(j,i) * wsvs(k,j,i) &
154	) * ddy &
155	+ 2.0 * ( &
156	km(k+1,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1) &
157	- km(k,j,i) * ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) &
158	) * ddzu(k+1)
159	ENDDO
160	ENDIF
161
162	ENDDO
163	ENDDO
164
165	END SUBROUTINE diffusion_w
166
167
168	!------------------------------------------------------------------------------!
169	! Call for grid point i,j
170	!------------------------------------------------------------------------------!
171	SUBROUTINE diffusion_w_ij( i, j )
172
173	USE arrays_3d
174	USE control_parameters
175	USE grid_variables
176	USE indices
177
178	IMPLICIT NONE
179
180	INTEGER :: i, j, k
181	REAL :: kmxm, kmxp, kmym, kmyp
182
183	REAL, DIMENSION(nzb:nzt+1) :: wsus, wsvs
184
185
186	DO k = nzb_w_outer(j,i)+1, nzt-1
187	!
188	!-- Interpolate eddy diffusivities on staggered gridpoints
189	kmxp = 0.25 * ( km(k,j,i)+km(k,j,i+1)+km(k+1,j,i)+km(k+1,j,i+1) )
190	kmxm = 0.25 * ( km(k,j,i)+km(k,j,i-1)+km(k+1,j,i)+km(k+1,j,i-1) )
191	kmyp = 0.25 * ( km(k,j,i)+km(k+1,j,i)+km(k,j+1,i)+km(k+1,j+1,i) )
192	kmym = 0.25 * ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
193
194	tend(k,j,i) = tend(k,j,i) &
195	& + ( kmxp * ( w(k,j,i+1) - w(k,j,i) ) * ddx &
196	& + kmxp * ( u(k+1,j,i+1) - u(k,j,i+1) ) * ddzu(k+1) &
197	& - kmxm * ( w(k,j,i) - w(k,j,i-1) ) * ddx &
198	& - kmxm * ( u(k+1,j,i) - u(k,j,i) ) * ddzu(k+1) &
199	& ) * ddx &
200	& + ( kmyp * ( w(k,j+1,i) - w(k,j,i) ) * ddy &
201	& + kmyp * ( v(k+1,j+1,i) - v(k,j+1,i) ) * ddzu(k+1) &
202	& - kmym * ( w(k,j,i) - w(k,j-1,i) ) * ddy &
203	& - kmym * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
204	& ) * ddy &
205	& + 2.0 * ( &
206	& km(k+1,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1) &
207	& - km(k,j,i) * ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) &
208	& ) * ddzu(k+1)
209	ENDDO
210
211	!
212	!-- Wall functions at all vertical walls, where necessary
213	IF ( wall_w_x(j,i) /= 0.0 .OR. wall_w_y(j,i) /= 0.0 ) THEN
214
215	!
216	!-- Calculate the horizontal momentum fluxes w'u' and/or w'v'
217	IF ( wall_w_x(j,i) /= 0.0 ) THEN
218	CALL wall_fluxes( i, j, nzb_w_inner(j,i)+1, nzb_w_outer(j,i), &
219	wsus, 0.0, 0.0, 0.0, 1.0 )
220	ELSE
221	wsus = 0.0
222	ENDIF
223
224	IF ( wall_w_y(j,i) /= 0.0 ) THEN
225	CALL wall_fluxes( i, j, nzb_w_inner(j,i)+1, nzb_w_outer(j,i), &
226	wsvs, 0.0, 0.0, 1.0, 0.0 )
227	ELSE
228	wsvs = 0.0
229	ENDIF
230
231	DO k = nzb_w_inner(j,i)+1, nzb_w_outer(j,i)
232	!
233	!-- Interpolate eddy diffusivities on staggered gridpoints
234	kmxp = 0.25 * ( km(k,j,i)+km(k,j,i+1)+km(k+1,j,i)+km(k+1,j,i+1) )
235	kmxm = 0.25 * ( km(k,j,i)+km(k,j,i-1)+km(k+1,j,i)+km(k+1,j,i-1) )
236	kmyp = 0.25 * ( km(k,j,i)+km(k+1,j,i)+km(k,j+1,i)+km(k+1,j+1,i) )
237	kmym = 0.25 * ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
238
239	tend(k,j,i) = tend(k,j,i) &
240	+ ( fwxp(j,i) * ( &
241	kmxp * ( w(k,j,i+1) - w(k,j,i) ) * ddx &
242	+ kmxp * ( u(k+1,j,i+1) - u(k,j,i+1) ) * ddzu(k+1) &
243	) &
244	- fwxm(j,i) * ( &
245	kmxm * ( w(k,j,i) - w(k,j,i-1) ) * ddx &
246	+ kmxm * ( u(k+1,j,i) - u(k,j,i) ) * ddzu(k+1) &
247	) &
248	+ wall_w_x(j,i) * wsus(k) &
249	) * ddx &
250	+ ( fwyp(j,i) * ( &
251	kmyp * ( w(k,j+1,i) - w(k,j,i) ) * ddy &
252	+ kmyp * ( v(k+1,j+1,i) - v(k,j+1,i) ) * ddzu(k+1) &
253	) &
254	- fwym(j,i) * ( &
255	kmym * ( w(k,j,i) - w(k,j-1,i) ) * ddy &
256	+ kmym * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
257	) &
258	+ wall_w_y(j,i) * wsvs(k) &
259	) * ddy &
260	+ 2.0 * ( &
261	km(k+1,j,i) * ( w(k+1,j,i) - w(k,j,i) ) * ddzw(k+1) &
262	- km(k,j,i) * ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) &
263	) * ddzu(k+1)
264	ENDDO
265	ENDIF
266
267	END SUBROUTINE diffusion_w_ij
268
269	END MODULE diffusion_w_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |