Home

Context Navigation

coriolis.f90 @ 4869

Last change on this file since 4869 was 4828, checked in by Giersch, 4 years ago
Copyright updated to year 2021, interface pmc_sort removed to accelarate the nesting code
Property svn:keywords set to `Id`
File size: 12.0 KB

Line
1	!> @file coriolis.f90
2	!--------------------------------------------------------------------------------------------------!
3	! This file is part of the PALM model system.
4	!
5	! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
6	! Public License as published by the Free Software Foundation, either version 3 of the License, or
7	! (at your option) any later version.
8	!
9	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
10	! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
11	! Public License for more details.
12	!
13	! You should have received a copy of the GNU General Public License along with PALM. If not, see
14	! <http://www.gnu.org/licenses/>.
15	!
16	! Copyright 1997-2021 Leibniz Universitaet Hannover
17	!--------------------------------------------------------------------------------------------------!
18	!
19	! Current revisions:
20	! -----------------
21	!
22	!
23	! Former revisions:
24	! -----------------
25	! $Id: coriolis.f90 4828 2021-01-05 11:21:41Z schwenkel $
26	! file re-formatted to follow the PALM coding standard
27	!
28	! 4360 2020-01-07 11:25:50Z suehring
29	! Introduction of wall_flags_total_0, which currently sets bits based on static topography
30	! information used in wall_flags_static_0
31	!
32	! 4329 2019-12-10 15:46:36Z motisi
33	! Renamed wall_flags_0 to wall_flags_static_0
34	!
35	! 4196 2019-08-29 11:02:06Z gronemeier
36	! Consider rotation of model domain
37	!
38	! 4182 2019-08-22 15:20:23Z scharf
39	! Corrected "Former revisions" section
40	!
41	! 3655 2019-01-07 16:51:22Z knoop
42	! OpenACC port for SPEC
43	!
44	! Revision 1.1 1997/08/29 08:57:38 raasch
45	! Initial revision
46	!
47	!
48	! Description:
49	! ------------
50	!> Computation of all Coriolis terms in the equations of motion.
51	!>
52	!> @note In this routine the topography is masked, even though this is again done in
53	!> prognostic_equations. However, omitting the masking here lead to slightly different results.
54	!> Reason unknown.
55	!--------------------------------------------------------------------------------------------------!
56	MODULE coriolis_mod
57
58
59	PRIVATE
60	PUBLIC coriolis
61
62	INTERFACE coriolis
63	MODULE PROCEDURE coriolis
64	MODULE PROCEDURE coriolis_ij
65	END INTERFACE coriolis
66
67	CONTAINS
68
69
70	!--------------------------------------------------------------------------------------------------!
71	! Description:
72	! ------------
73	!> Call for all grid points
74	!--------------------------------------------------------------------------------------------------!
75	SUBROUTINE coriolis( component )
76
77	USE arrays_3d, &
78	ONLY: tend, u, ug, v, vg, w
79
80	USE basic_constants_and_equations_mod, &
81	ONLY: pi
82
83	USE control_parameters, &
84	ONLY: f, fs, message_string, rotation_angle
85
86	USE indices, &
87	ONLY: nxl, nxlu, nxr, nyn, nys, nysv, nzb, nzt, wall_flags_total_0
88
89	USE kinds
90
91	IMPLICIT NONE
92
93	INTEGER(iwp) :: component !< component of momentum equation
94	INTEGER(iwp) :: i !< running index x direction
95	INTEGER(iwp) :: j !< running index y direction
96	INTEGER(iwp) :: k !< running index z direction
97
98	REAL(wp) :: cos_rot_angle !< cosine of model rotation angle
99	REAL(wp) :: flag !< flag to mask topography
100	REAL(wp) :: sin_rot_angle !< sine of model rotation angle
101
102	!
103	!-- Precalculate cosine and sine of rotation angle
104	cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
105	sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
106
107	!
108	!-- Compute Coriolis terms for the three velocity components
109	SELECT CASE ( component )
110
111	!
112	!-- u-component
113	CASE ( 1 )
114	!$ACC PARALLEL LOOP COLLAPSE(3) PRIVATE(i, j, k, flag) &
115	!$ACC PRESENT(wall_flags_total_0) &
116	!$ACC PRESENT(v, w, vg) &
117	!$ACC PRESENT(tend)
118	DO i = nxlu, nxr
119	DO j = nys, nyn
120	DO k = nzb+1, nzt
121	!
122	!-- Predetermine flag to mask topography
123	flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 1 ) )
124
125	tend(k,j,i) = tend(k,j,i) + flag * &
126	( f &
127	* ( 0.25_wp * ( v(k,j,i-1) + v(k,j,i) + v(k,j+1,i-1) + v(k,j+1,i) ) &
128	- vg(k) ) &
129	- fs * cos_rot_angle &
130	* 0.25_wp * ( w(k-1,j,i-1) + w(k-1,j,i) + w(k,j,i-1) + w(k,j,i) ) &
131	)
132	ENDDO
133	ENDDO
134	ENDDO
135
136	!
137	!-- v-component
138	CASE ( 2 )
139	!$ACC PARALLEL LOOP COLLAPSE(3) PRIVATE(i, j, k, flag) &
140	!$ACC PRESENT(wall_flags_total_0) &
141	!$ACC PRESENT(u, w, ug) &
142	!$ACC PRESENT(tend)
143	DO i = nxl, nxr
144	DO j = nysv, nyn
145	DO k = nzb+1, nzt
146	!
147	!-- Predetermine flag to mask topography
148	flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2 ) )
149
150	tend(k,j,i) = tend(k,j,i) - flag * &
151	( f &
152	* ( 0.25_wp * ( u(k,j-1,i) + u(k,j,i) + u(k,j-1,i+1) + u(k,j,i+1) ) &
153	- ug(k) ) &
154	+ fs * sin_rot_angle &
155	* 0.25_wp * ( w(k,j,i) + w(k-1,j,i) + w(k,j-1,i) + w(k-1,j-1,i) ) &
156	)
157	ENDDO
158	ENDDO
159	ENDDO
160
161	!
162	!-- w-component
163	CASE ( 3 )
164	!$ACC PARALLEL LOOP COLLAPSE(3) PRIVATE(i, j, k, flag) &
165	!$ACC PRESENT(wall_flags_total_0) &
166	!$ACC PRESENT(u, v) &
167	!$ACC PRESENT(tend)
168	DO i = nxl, nxr
169	DO j = nys, nyn
170	DO k = nzb+1, nzt
171	!
172	!-- Predetermine flag to mask topography
173	flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ) )
174
175	tend(k,j,i) = tend(k,j,i) &
176	+ fs * 0.25_wp * flag &
177	* ( cos_rot_angle &
178	* ( u(k,j,i) + u(k+1,j,i) + u(k,j,i+1) + u(k+1,j,i+1) ) &
179	+ sin_rot_angle &
180	* ( v(k,j,i) + v(k+1,j,i) + v(k,j+1,i) + v(k+1,j+1,i) ) &
181	)
182	ENDDO
183	ENDDO
184	ENDDO
185
186	CASE DEFAULT
187
188	WRITE( message_string, * ) ' wrong component: ', component
189	CALL message( 'coriolis', 'PA0173', 1, 2, 0, 6, 0 )
190
191	END SELECT
192
193	END SUBROUTINE coriolis
194
195
196	!--------------------------------------------------------------------------------------------------!
197	! Description:
198	! ------------
199	!> Call for grid point i,j
200	!--------------------------------------------------------------------------------------------------!
201	SUBROUTINE coriolis_ij( i, j, component )
202
203	USE arrays_3d, &
204	ONLY: tend, u, ug, v, vg, w
205
206	USE basic_constants_and_equations_mod, &
207	ONLY: pi
208
209	USE control_parameters, &
210	ONLY: f, fs, message_string, rotation_angle
211
212	USE indices, &
213	ONLY: nzb, nzt, wall_flags_total_0
214
215	USE kinds
216
217	IMPLICIT NONE
218
219	INTEGER(iwp) :: component !< component of momentum equation
220	INTEGER(iwp) :: i !< running index x direction
221	INTEGER(iwp) :: j !< running index y direction
222	INTEGER(iwp) :: k !< running index z direction
223
224	REAL(wp) :: cos_rot_angle !< cosine of model rotation angle
225	REAL(wp) :: flag !< flag to mask topography
226	REAL(wp) :: sin_rot_angle !< sine of model rotation angle
227
228	!
229	!-- Precalculate cosine and sine of rotation angle
230	cos_rot_angle = COS( rotation_angle * pi / 180.0_wp )
231	sin_rot_angle = SIN( rotation_angle * pi / 180.0_wp )
232
233	!
234	!-- Compute Coriolis terms for the three velocity components
235	SELECT CASE ( component )
236
237	!
238	!-- u-component
239	CASE ( 1 )
240	DO k = nzb+1, nzt
241	!
242	!-- Predetermine flag to mask topography
243	flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 1 ) )
244
245	tend(k,j,i) = tend(k,j,i) + flag * &
246	( f &
247	* ( 0.25_wp * ( v(k,j,i-1) + v(k,j,i) + v(k,j+1,i-1) + v(k,j+1,i) ) &
248	- vg(k) ) &
249	- fs * cos_rot_angle &
250	* 0.25_wp * ( w(k-1,j,i-1) + w(k-1,j,i) + w(k,j,i-1) + w(k,j,i) ) &
251	)
252	ENDDO
253
254	!
255	!-- v-component
256	CASE ( 2 )
257	DO k = nzb+1, nzt
258	!
259	!-- Predetermine flag to mask topography
260	flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2 ) )
261
262	tend(k,j,i) = tend(k,j,i) - flag * &
263	( f &
264	* ( 0.25_wp * ( u(k,j-1,i) + u(k,j,i) + u(k,j-1,i+1) + u(k,j,i+1) )&
265	- ug(k) ) &
266	+ fs * sin_rot_angle &
267	* 0.25_wp * ( w(k,j,i) + w(k-1,j,i) + w(k,j-1,i) + w(k-1,j-1,i) ) &
268	)
269	ENDDO
270
271	!
272	!-- w-component
273	CASE ( 3 )
274	DO k = nzb+1, nzt
275	!
276	!-- Predetermine flag to mask topography
277	flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ) )
278
279	tend(k,j,i) = tend(k,j,i) &
280	+ fs * 0.25_wp * flag &
281	* ( cos_rot_angle &
282	* ( u(k,j,i) + u(k+1,j,i) + u(k,j,i+1) + u(k+1,j,i+1) ) &
283	+ sin_rot_angle &
284	* ( v(k,j,i) + v(k+1,j,i) + v(k,j+1,i) + v(k+1,j+1,i) ) &
285	)
286	ENDDO
287
288	CASE DEFAULT
289
290	WRITE( message_string, * ) ' wrong component: ', component
291	CALL message( 'coriolis', 'PA0173', 1, 2, 0, 6, 0 )
292
293	END SELECT
294
295	END SUBROUTINE coriolis_ij
296
297	END MODULE coriolis_mod

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

Context Navigation

source: palm/trunk/SOURCE/coriolis.f90 @ 4869

Download in other formats: