Home

Context Navigation

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

Last change on this file since 4276 was 4196, checked in by gronemeier, 5 years ago

Consider rotation of model domain for claculation of Coriolis force:

check_parameters: Overwrite rotation_angle from namelist by value from static driver;
coriolis: Consider rotation of model domain;
header: Write information about rotation angle;
modules: Added rotation_angle;
netcdf_interface_mod: replaced rotation angle from input-netCDF file by namelist parameter 'rotation_angle';
ocean_mod: Consider rotation of model domain for calculating the Stokes drift;
parin: added rotation_angle to initialization_parameters;

Property svn:keywords set to Id

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |