Home

Context Navigation

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

Last change on this file since 4330 was 4329, checked in by motisi, 5 years ago
Renamed wall_flags_0 to wall_flags_static_0
Property svn:keywords set to `Id`
File size: 11.4 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |