Home

Context Navigation

source: palm/trunk/SOURCE/ls_forcing_mod.f90 @ 2234

Last change on this file since 2234 was 2233, checked in by suehring, 8 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 24.9 KB

Rev	Line
[1850]	1	!> @file ls_forcing_mod.f90
[2000]	2	!------------------------------------------------------------------------------!
[1239]	3	! This file is part of PALM.
	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.
[1239]	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	!
[2101]	17	! Copyright 1997-2017 Leibniz Universitaet Hannover
[2000]	18	!------------------------------------------------------------------------------!
[1239]	19	!
	20	! Current revisions:
	21	! ------------------
[2233]	22	!
	23	!
	24	! Former revisions:
	25	! -----------------
	26	! $Id: ls_forcing_mod.f90 2233 2017-05-30 18:08:54Z suehring $
	27	!
	28	! 2232 2017-05-30 17:47:52Z suehring
[2232]	29	! Adopt to new topography structure, even though no well-conceived topography
	30	! concept concerning nudging and large-scale for exist so far.
	31	!
	32	! Also adopt to new surface-structure, i.e. fluxes are obtained from data-types
[1383]	33	!
[2105]	34	! 2104 2017-01-06 16:01:15Z knoop
	35	! Bugfix for approximation related flux input conversion
	36	!
[2072]	37	! 2071 2016-11-17 11:22:14Z maronga
	38	! Small bugfix (Resler)
	39	!
[2038]	40	! 2037 2016-10-26 11:15:40Z knoop
	41	! Anelastic approximation implemented
	42	!
[2001]	43	! 2000 2016-08-20 18:09:15Z knoop
	44	! Forced header and separation lines into 80 columns
	45	!
[1851]	46	! 1850 2016-04-08 13:29:27Z maronga
	47	! Module renamed
	48	!
	49	!
[1683]	50	! 1682 2015-10-07 23:56:08Z knoop
	51	! Code annotations made doxygen readable
	52	!
[1603]	53	! 1602 2015-06-22 07:50:56Z heinze
	54	! PA0370 changed to PA0363
	55	!
[1383]	56	! 1382 2014-04-30 12:15:41Z boeske
[1382]	57	! Renamed variables which store large scale forcing tendencies
	58	! pt_lsa -> td_lsa_lpt, pt_subs -> td_sub_lpt,
	59	! q_lsa -> td_lsa_q, q_subs -> td_sub_q,
	60	! high\|lowpt_lsa -> high\|low_td_lsa_lpt, ...
[1366]	61	!
	62	! 1365 2014-04-22 15:03:56Z boeske
[1365]	63	! Usage of large scale forcing for pt and q enabled:
	64	! Added new subroutine ls_advec for horizontal large scale advection and large
	65	! scale subsidence,
	66	! error message in init_ls_forcing specified,
	67	! variable t renamed nt
[1354]	68	!
	69	! 1353 2014-04-08 15:21:23Z heinze
	70	! REAL constants provided with KIND-attribute
	71	!
[1321]	72	! 1320 2014-03-20 08:40:49Z raasch
[1320]	73	! ONLY-attribute added to USE-statements,
	74	! kind-parameters added to all INTEGER and REAL declaration statements,
	75	! kinds are defined in new module kinds,
	76	! comment fields (!:) to be used for variable explanations added to
	77	! all variable declaration statements
[1239]	78	!
[1319]	79	! 1318 2014-03-17 13:35:16Z raasch
	80	! module interfaces removed
	81	!
[1300]	82	! 1299 2014-03-06 13:15:21Z heinze
	83	! Ensure a zero large scale vertical velocity at the surface
	84	! Bugfix: typo in case of boundary condition in if-clause
	85	!
[1277]	86	! 1276 2014-01-15 13:40:41Z heinze
	87	! Use LSF_DATA also in case of Dirichlet bottom boundary condition for scalars
	88	!
[1250]	89	! 1249 2013-11-06 10:45:47Z heinze
	90	! remove call of user module
	91	! reformatting
	92	!
[1242]	93	! 1241 2013-10-30 11:36:58Z heinze
	94	! Initial revision
[1239]	95	!
	96	! Description:
	97	! ------------
[1682]	98	!> Calculates large scale forcings (geostrophic wind and subsidence velocity) as
	99	!> well as surfaces fluxes dependent on time given in an external file (LSF_DATA).
	100	!> Code is based in parts on DALES and UCLA-LES.
[1239]	101	!--------------------------------------------------------------------------------!
[1682]	102	MODULE ls_forcing_mod
	103
[1239]	104
	105	PRIVATE
[1365]	106	PUBLIC init_ls_forcing, ls_forcing_surf, ls_forcing_vert, ls_advec
[1239]	107	SAVE
	108
[1365]	109	INTERFACE ls_advec
	110	MODULE PROCEDURE ls_advec
	111	MODULE PROCEDURE ls_advec_ij
	112	END INTERFACE ls_advec
[1239]	113
	114	CONTAINS
	115
[1682]	116	!------------------------------------------------------------------------------!
	117	! Description:
	118	! ------------
	119	!> @todo Missing subroutine description.
	120	!------------------------------------------------------------------------------!
[1239]	121	SUBROUTINE init_ls_forcing
	122
[1320]	123	USE arrays_3d, &
[1382]	124	ONLY: p_surf, pt_surf, q_surf, qsws_surf, shf_surf, td_lsa_lpt, &
	125	td_lsa_q, td_sub_lpt, td_sub_q, time_surf, time_vert, &
	126	ug_vert, vg_vert, wsubs_vert, zu
[1239]	127
[1320]	128	USE control_parameters, &
	129	ONLY: end_time, lsf_surf, lsf_vert, message_string, nlsf
	130
	131	USE indices, &
[1365]	132	ONLY: ngp_sums_ls, nzb, nz, nzt
[1320]	133
	134	USE kinds
	135
[1365]	136	USE statistics, &
	137	ONLY: sums_ls_l
[1320]	138
[1365]	139
[1239]	140	IMPLICIT NONE
	141
[1682]	142	CHARACTER(100) :: chmess !<
	143	CHARACTER(1) :: hash !<
[1239]	144
[1682]	145	INTEGER(iwp) :: ierrn !<
	146	INTEGER(iwp) :: finput = 90 !<
	147	INTEGER(iwp) :: k !<
	148	INTEGER(iwp) :: nt !<
[1320]	149
[1682]	150	REAL(wp) :: fac !<
	151	REAL(wp) :: highheight !<
	152	REAL(wp) :: highug_vert !<
	153	REAL(wp) :: highvg_vert !<
	154	REAL(wp) :: highwsubs_vert !<
	155	REAL(wp) :: lowheight !<
	156	REAL(wp) :: lowug_vert !<
	157	REAL(wp) :: lowvg_vert !<
	158	REAL(wp) :: lowwsubs_vert !<
	159	REAL(wp) :: high_td_lsa_lpt !<
	160	REAL(wp) :: low_td_lsa_lpt !<
	161	REAL(wp) :: high_td_lsa_q !<
	162	REAL(wp) :: low_td_lsa_q !<
	163	REAL(wp) :: high_td_sub_lpt !<
	164	REAL(wp) :: low_td_sub_lpt !<
	165	REAL(wp) :: high_td_sub_q !<
	166	REAL(wp) :: low_td_sub_q !<
	167	REAL(wp) :: r_dummy !<
[1320]	168
[1382]	169	ALLOCATE( p_surf(0:nlsf), pt_surf(0:nlsf), q_surf(0:nlsf), &
	170	qsws_surf(0:nlsf), shf_surf(0:nlsf), &
	171	td_lsa_lpt(nzb:nzt+1,0:nlsf), td_lsa_q(nzb:nzt+1,0:nlsf), &
	172	td_sub_lpt(nzb:nzt+1,0:nlsf), td_sub_q(nzb:nzt+1,0:nlsf), &
[1365]	173	time_vert(0:nlsf), time_surf(0:nlsf), &
	174	ug_vert(nzb:nzt+1,0:nlsf), vg_vert(nzb:nzt+1,0:nlsf), &
	175	wsubs_vert(nzb:nzt+1,0:nlsf) )
[1239]	176
[1382]	177	p_surf = 0.0_wp; pt_surf = 0.0_wp; q_surf = 0.0_wp; qsws_surf = 0.0_wp
	178	shf_surf = 0.0_wp; time_vert = 0.0_wp; td_lsa_lpt = 0.0_wp
	179	td_lsa_q = 0.0_wp; td_sub_lpt = 0.0_wp; td_sub_q = 0.0_wp
	180	time_surf = 0.0_wp; ug_vert = 0.0_wp; vg_vert = 0.0_wp
	181	wsubs_vert = 0.0_wp
[1239]	182
[1365]	183	!
	184	!-- Array for storing large scale forcing and nudging tendencies at each
	185	!-- timestep for data output
	186	ALLOCATE( sums_ls_l(nzb:nzt+1,0:7) )
	187	sums_ls_l = 0.0_wp
[1239]	188
[1365]	189	ngp_sums_ls = (nz+2)*6
	190
[1249]	191	OPEN ( finput, FILE='LSF_DATA', STATUS='OLD', &
	192	FORM='FORMATTED', IOSTAT=ierrn )
[1239]	193
[1249]	194	IF ( ierrn /= 0 ) THEN
[1239]	195	message_string = 'file LSF_DATA does not exist'
	196	CALL message( 'ls_forcing', 'PA0368', 1, 2, 0, 6, 0 )
	197	ENDIF
	198
	199	ierrn = 0
	200	!
	201	!-- First three lines of LSF_DATA contain header
[1249]	202	READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
	203	READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
	204	READ ( finput, FMT='(a100)', IOSTAT=ierrn ) chmess
[1239]	205
[1249]	206	IF ( ierrn /= 0 ) THEN
[1239]	207	message_string = 'errors in file LSF_DATA'
	208	CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
	209	ENDIF
	210
	211	!
	212	!-- Surface values are read in
[1365]	213	nt = 0
[1239]	214	ierrn = 0
	215
[1365]	216	DO WHILE ( time_surf(nt) < end_time )
	217	nt = nt + 1
	218	READ ( finput, *, IOSTAT = ierrn ) time_surf(nt), shf_surf(nt), &
	219	qsws_surf(nt), pt_surf(nt), &
	220	q_surf(nt), p_surf(nt)
[1239]	221
[1249]	222	IF ( ierrn < 0 ) THEN
[1299]	223	WRITE ( message_string, * ) 'No time dependent surface variables ',&
[1239]	224	'in&LSF_DATA for end of run found'
	225
[1602]	226	CALL message( 'ls_forcing', 'PA0363', 1, 2, 0, 6, 0 )
[1239]	227	ENDIF
	228	ENDDO
	229
[1249]	230	IF ( time_surf(1) > end_time ) THEN
[1299]	231	WRITE ( message_string, * ) 'No time dependent surface variables in ',&
[1365]	232	'&LSF_DATA for end of run found - ', &
[1239]	233	'lsf_surf is set to FALSE'
	234	CALL message( 'ls_forcing', 'PA0371', 0, 0, 0, 6, 0 )
	235	lsf_surf = .FALSE.
	236	ENDIF
	237
	238	!
	239	!-- Go to the end of the list with surface variables
	240	DO WHILE ( ierrn == 0 )
[1249]	241	READ ( finput, *, IOSTAT = ierrn ) r_dummy
[1239]	242	ENDDO
	243
	244	!
	245	!-- Profiles of ug, vg and w_subs are read in (large scale forcing)
	246
[1365]	247	nt = 0
	248	DO WHILE ( time_vert(nt) < end_time )
	249	nt = nt + 1
[1239]	250	hash = "#"
	251	ierrn = 1 ! not zero
	252	!
	253	!-- Search for the next line consisting of "# time",
	254	!-- from there onwards the profiles will be read
	255	DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) )
[1365]	256	READ ( finput, *, IOSTAT=ierrn ) hash, time_vert(nt)
[1249]	257	IF ( ierrn < 0 ) THEN
[1299]	258	WRITE( message_string, * ) 'No time dependent vertical profiles',&
[1239]	259	' in&LSF_DATA for end of run found'
	260	CALL message( 'ls_forcing', 'PA0372', 1, 2, 0, 6, 0 )
	261	ENDIF
	262	ENDDO
	263
[1365]	264	IF ( nt == 1 .AND. time_vert(nt) > end_time ) EXIT
[1239]	265
[1365]	266	READ ( finput, *, IOSTAT=ierrn ) lowheight, lowug_vert, lowvg_vert, &
[1382]	267	lowwsubs_vert, low_td_lsa_lpt, &
	268	low_td_lsa_q, low_td_sub_lpt, &
	269	low_td_sub_q
[1249]	270	IF ( ierrn /= 0 ) THEN
[1239]	271	message_string = 'errors in file LSF_DATA'
[1299]	272	CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
[1239]	273	ENDIF
	274
[1365]	275	READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert, &
	276	highvg_vert, highwsubs_vert, &
[1382]	277	high_td_lsa_lpt, high_td_lsa_q, &
	278	high_td_sub_lpt, high_td_sub_q
[1239]	279
[1249]	280	IF ( ierrn /= 0 ) THEN
[1239]	281	message_string = 'errors in file LSF_DATA'
[1299]	282	CALL message( 'ls_forcing', 'PA0369', 1, 2, 0, 6, 0 )
[1239]	283	ENDIF
	284
[1241]	285
[1239]	286	DO k = nzb, nzt+1
[1249]	287	IF ( highheight < zu(k) ) THEN
[1382]	288	lowheight = highheight
	289	lowug_vert = highug_vert
	290	lowvg_vert = highvg_vert
	291	lowwsubs_vert = highwsubs_vert
	292	low_td_lsa_lpt = high_td_lsa_lpt
	293	low_td_lsa_q = high_td_lsa_q
	294	low_td_sub_lpt = high_td_sub_lpt
	295	low_td_sub_q = high_td_sub_q
[1239]	296
	297	ierrn = 0
[1365]	298	READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert, &
	299	highvg_vert, highwsubs_vert, &
[1382]	300	high_td_lsa_lpt, &
	301	high_td_lsa_q, &
	302	high_td_sub_lpt, high_td_sub_q
[1239]	303
[1249]	304	IF ( ierrn /= 0 ) THEN
[1365]	305	WRITE( message_string, * ) 'zu(nzt+1) = ', zu(nzt+1), 'm ', &
	306	'is higher than the maximum height in LSF_DATA which ',&
	307	'is ', lowheight, 'm. Interpolation on PALM ', &
	308	'grid is not possible.'
	309	CALL message( 'ls_forcing', 'PA0395', 1, 2, 0, 6, 0 )
[1239]	310	ENDIF
	311
	312	ENDIF
	313
	314	!
	315	!-- Interpolation of prescribed profiles in space
	316	fac = (highheight-zu(k))/(highheight - lowheight)
	317
[1382]	318	ug_vert(k,nt) = fac * lowug_vert &
	319	+ ( 1.0_wp - fac ) * highug_vert
	320	vg_vert(k,nt) = fac * lowvg_vert &
	321	+ ( 1.0_wp - fac ) * highvg_vert
	322	wsubs_vert(k,nt) = fac * lowwsubs_vert &
	323	+ ( 1.0_wp - fac ) * highwsubs_vert
[1239]	324
[1382]	325	td_lsa_lpt(k,nt) = fac * low_td_lsa_lpt &
	326	+ ( 1.0_wp - fac ) * high_td_lsa_lpt
	327	td_lsa_q(k,nt) = fac * low_td_lsa_q &
	328	+ ( 1.0_wp - fac ) * high_td_lsa_q
	329	td_sub_lpt(k,nt) = fac * low_td_sub_lpt &
	330	+ ( 1.0_wp - fac ) * high_td_sub_lpt
	331	td_sub_q(k,nt) = fac * low_td_sub_q &
	332	+ ( 1.0_wp - fac ) * high_td_sub_q
[1365]	333
[1239]	334	ENDDO
	335
	336	ENDDO
	337
[1299]	338	!
	339	!-- Large scale vertical velocity has to be zero at the surface
[1353]	340	wsubs_vert(nzb,:) = 0.0_wp
[1299]	341
[1249]	342	IF ( time_vert(1) > end_time ) THEN
[1365]	343	WRITE ( message_string, * ) 'Time dependent large scale profile ', &
	344	'forcing from&LSF_DATA sets in after end of ' , &
[1249]	345	'simulation - lsf_vert is set to FALSE'
[1239]	346	CALL message( 'ls_forcing', 'PA0373', 0, 0, 0, 6, 0 )
	347	lsf_vert = .FALSE.
	348	ENDIF
	349
	350	CLOSE( finput )
	351
	352
	353	END SUBROUTINE init_ls_forcing
	354
	355
[1682]	356	!------------------------------------------------------------------------------!
	357	! Description:
	358	! ------------
	359	!> @todo Missing subroutine description.
	360	!------------------------------------------------------------------------------!
[1239]	361	SUBROUTINE ls_forcing_surf ( time )
	362
[1320]	363	USE arrays_3d, &
[2232]	364	ONLY: p_surf, pt_surf, q_surf, qsws_surf, shf_surf, &
[2104]	365	heatflux_input_conversion, waterflux_input_conversion, &
[1320]	366	time_surf, time_vert, ug, ug_vert, vg, vg_vert
[1239]	367
[1320]	368	USE control_parameters, &
	369	ONLY: bc_q_b, ibc_pt_b, ibc_q_b, pt_surface, q_surface, &
	370	surface_pressure
	371
[2104]	372	USE indices, &
	373	ONLY: nzb
	374
[1320]	375	USE kinds
	376
[2232]	377	USE surface_mod, &
	378	ONLY: surf_def_h, surf_lsm_h, surf_usm_h
	379
[1239]	380	IMPLICIT NONE
	381
[1682]	382	INTEGER(iwp) :: nt !<
[1239]	383
[2232]	384	REAL(wp) :: dum_surf_flux !<
[1682]	385	REAL(wp) :: fac !<
	386	REAL(wp), INTENT(in) :: time !<
[1320]	387
[1239]	388	!
	389	!-- Interpolation in time of LSF_DATA at the surface
[1365]	390	nt = 1
	391	DO WHILE ( time > time_surf(nt) )
	392	nt = nt + 1
[1239]	393	ENDDO
[1365]	394	IF ( time /= time_surf(nt) ) THEN
	395	nt = nt - 1
[1239]	396	ENDIF
	397
[1365]	398	fac = ( time -time_surf(nt) ) / ( time_surf(nt+1) - time_surf(nt) )
[1239]	399
[1276]	400	IF ( ibc_pt_b == 0 ) THEN
	401	!
	402	!-- In case of Dirichlet boundary condition shf must not
	403	!-- be set - it is calculated via MOST in prandtl_fluxes
[1365]	404	pt_surface = pt_surf(nt) + fac * ( pt_surf(nt+1) - pt_surf(nt) )
[1276]	405
	406	ELSEIF ( ibc_pt_b == 1 ) THEN
	407	!
	408	!-- In case of Neumann boundary condition pt_surface is needed for
	409	!-- calculation of reference density
[2232]	410	dum_surf_flux = ( shf_surf(nt) + fac * &
	411	( shf_surf(nt+1) - shf_surf(nt) ) &
	412	) * heatflux_input_conversion(nzb)
	413	!
	414	!-- Save surface sensible heat flux on default, natural and urban surface
	415	!-- type, if required
	416	IF ( surf_def_h(0)%ns >= 1 ) surf_def_h(0)%shf(:) = dum_surf_flux
	417	IF ( surf_lsm_h%ns >= 1 ) surf_lsm_h%shf(:) = dum_surf_flux
	418	IF ( surf_usm_h%ns >= 1 ) surf_usm_h%shf(:) = dum_surf_flux
[1276]	419
[2232]	420	pt_surface = pt_surf(nt) + fac * ( pt_surf(nt+1) - pt_surf(nt) )
	421
[1276]	422	ENDIF
	423
	424	IF ( ibc_q_b == 0 ) THEN
	425	!
	426	!-- In case of Dirichlet boundary condition qsws must not
	427	!-- be set - it is calculated via MOST in prandtl_fluxes
[1365]	428	q_surface = q_surf(nt) + fac * ( q_surf(nt+1) - q_surf(nt) )
[1276]	429
[1299]	430	ELSEIF ( ibc_q_b == 1 ) THEN
[2232]	431	dum_surf_flux = ( qsws_surf(nt) + fac * &
	432	( qsws_surf(nt+1) - qsws_surf(nt) ) &
	433	) * waterflux_input_conversion(nzb)
	434	!
	435	!-- Save surface sensible heat flux on default, natural and urban surface
	436	!-- type, if required
	437	IF ( surf_def_h(0)%ns >= 1 ) surf_def_h(0)%qsws(:) = dum_surf_flux
	438	IF ( surf_lsm_h%ns >= 1 ) surf_lsm_h%qsws(:) = dum_surf_flux
	439	IF ( surf_usm_h%ns >= 1 ) surf_usm_h%qsws(:) = dum_surf_flux
[1276]	440
	441	ENDIF
	442
[1365]	443	surface_pressure = p_surf(nt) + fac * ( p_surf(nt+1) - p_surf(nt) )
[1239]	444
	445	END SUBROUTINE ls_forcing_surf
	446
	447
[1682]	448	!------------------------------------------------------------------------------!
	449	! Description:
	450	! ------------
	451	!> @todo Missing subroutine description.
	452	!------------------------------------------------------------------------------!
[1239]	453	SUBROUTINE ls_forcing_vert ( time )
	454
[1320]	455	USE arrays_3d, &
	456	ONLY: time_vert, ug, ug_vert, vg, vg_vert, w_subs, wsubs_vert
[1239]	457
[1320]	458	USE control_parameters, &
	459	ONLY: large_scale_subsidence
	460
	461	USE kinds
	462
[1239]	463	IMPLICIT NONE
	464
[1682]	465	INTEGER(iwp) :: nt !<
[1239]	466
[1682]	467	REAL(wp) :: fac !<
	468	REAL(wp), INTENT(in) :: time !<
[1320]	469
[1239]	470	!
	471	!-- Interpolation in time of LSF_DATA for ug, vg and w_subs
[1365]	472	nt = 1
	473	DO WHILE ( time > time_vert(nt) )
	474	nt = nt + 1
[1239]	475	ENDDO
[1365]	476	IF ( time /= time_vert(nt) ) THEN
	477	nt = nt - 1
[1239]	478	ENDIF
	479
[1365]	480	fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
[1239]	481
[1365]	482	ug = ug_vert(:,nt) + fac * ( ug_vert(:,nt+1) - ug_vert(:,nt) )
	483	vg = vg_vert(:,nt) + fac * ( vg_vert(:,nt+1) - vg_vert(:,nt) )
[1239]	484
[1249]	485	IF ( large_scale_subsidence ) THEN
[1365]	486	w_subs = wsubs_vert(:,nt) &
	487	+ fac * ( wsubs_vert(:,nt+1) - wsubs_vert(:,nt) )
[1239]	488	ENDIF
	489
	490	END SUBROUTINE ls_forcing_vert
	491
	492
[1365]	493	!------------------------------------------------------------------------------!
[1682]	494	! Description:
	495	! ------------
	496	!> Call for all grid points
[1365]	497	!------------------------------------------------------------------------------!
	498	SUBROUTINE ls_advec ( time, prog_var )
	499
	500	USE arrays_3d, &
[1382]	501	ONLY: td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, tend, time_vert
[1365]	502
	503	USE control_parameters, &
	504	ONLY: large_scale_subsidence, use_subsidence_tendencies
	505
	506	USE indices
	507
	508	USE kinds
	509
	510	IMPLICIT NONE
	511
[1682]	512	CHARACTER (LEN=*) :: prog_var !<
[1365]	513
[1682]	514	REAL(wp), INTENT(in) :: time !<
	515	REAL(wp) :: fac !<
[1365]	516
[1682]	517	INTEGER(iwp) :: i !<
	518	INTEGER(iwp) :: j !<
	519	INTEGER(iwp) :: k !<
	520	INTEGER(iwp) :: nt !<
[1365]	521
	522	!
	523	!-- Interpolation in time of LSF_DATA
	524	nt = 1
	525	DO WHILE ( time > time_vert(nt) )
	526	nt = nt + 1
	527	ENDDO
	528	IF ( time /= time_vert(nt) ) THEN
	529	nt = nt - 1
	530	ENDIF
	531
	532	fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
	533
	534	!
	535	!-- Add horizontal large scale advection tendencies of pt and q
	536	SELECT CASE ( prog_var )
	537
	538	CASE ( 'pt' )
	539
	540	DO i = nxl, nxr
	541	DO j = nys, nyn
[2232]	542	DO k = nzb+1, nzt
[1382]	543	tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt) + fac * &
[2232]	544	( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) ) *&
	545	MERGE( 1.0_wp, 0.0_wp, &
	546	BTEST( wall_flags_0(k,j,i), 0 ) )
[1365]	547	ENDDO
	548	ENDDO
	549	ENDDO
	550
	551	CASE ( 'q' )
	552
	553	DO i = nxl, nxr
	554	DO j = nys, nyn
[2232]	555	DO k = nzb+1, nzt
[1382]	556	tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt) + fac * &
[2232]	557	( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) ) * &
	558	MERGE( 1.0_wp, 0.0_wp, &
	559	BTEST( wall_flags_0(k,j,i), 0 ) )
[1365]	560	ENDDO
	561	ENDDO
	562	ENDDO
	563
	564	END SELECT
	565
	566	!
	567	!-- Subsidence of pt and q with prescribed subsidence tendencies
	568	IF ( large_scale_subsidence .AND. use_subsidence_tendencies ) THEN
	569
	570	SELECT CASE ( prog_var )
	571
	572	CASE ( 'pt' )
	573
	574	DO i = nxl, nxr
	575	DO j = nys, nyn
[2232]	576	DO k = nzb+1, nzt
[1382]	577	tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac * &
[2232]	578	( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) )*&
	579	MERGE( 1.0_wp, 0.0_wp, &
	580	BTEST( wall_flags_0(k,j,i), 0 ) )
[1365]	581	ENDDO
	582	ENDDO
	583	ENDDO
	584
	585	CASE ( 'q' )
	586
	587	DO i = nxl, nxr
	588	DO j = nys, nyn
[2232]	589	DO k = nzb+1, nzt
[1382]	590	tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac * &
[2232]	591	( td_sub_q(k,nt+1) - td_sub_q(k,nt) ) * &
	592	MERGE( 1.0_wp, 0.0_wp, &
	593	BTEST( wall_flags_0(k,j,i), 0 ) )
[1365]	594	ENDDO
	595	ENDDO
	596	ENDDO
	597
	598	END SELECT
	599
	600	ENDIF
	601
	602	END SUBROUTINE ls_advec
	603
	604
	605	!------------------------------------------------------------------------------!
[1682]	606	! Description:
	607	! ------------
	608	!> Call for grid point i,j
[1365]	609	!------------------------------------------------------------------------------!
	610	SUBROUTINE ls_advec_ij ( i, j, time, prog_var )
	611
	612	USE arrays_3d, &
[1382]	613	ONLY: td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, tend, time_vert
[1365]	614
	615	USE control_parameters, &
	616	ONLY: large_scale_subsidence, use_subsidence_tendencies
	617
	618	USE indices
	619
	620	USE kinds
	621
	622	IMPLICIT NONE
	623
[1682]	624	CHARACTER (LEN=*) :: prog_var !<
[1365]	625
[1682]	626	REAL(wp), INTENT(in) :: time !<
	627	REAL(wp) :: fac !<
[1365]	628
[1682]	629	INTEGER(iwp) :: i !<
	630	INTEGER(iwp) :: j !<
	631	INTEGER(iwp) :: k !<
	632	INTEGER(iwp) :: nt !<
[1365]	633
	634	!
	635	!-- Interpolation in time of LSF_DATA
	636	nt = 1
	637	DO WHILE ( time > time_vert(nt) )
	638	nt = nt + 1
	639	ENDDO
	640	IF ( time /= time_vert(nt) ) THEN
	641	nt = nt - 1
	642	ENDIF
	643
	644	fac = ( time-time_vert(nt) ) / ( time_vert(nt+1)-time_vert(nt) )
	645
	646	!
	647	!-- Add horizontal large scale advection tendencies of pt and q
	648	SELECT CASE ( prog_var )
	649
	650	CASE ( 'pt' )
	651
[2232]	652	DO k = nzb+1, nzt
[1382]	653	tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt) &
[2232]	654	+ fac * ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) )*&
	655	MERGE( 1.0_wp, 0.0_wp, &
	656	BTEST( wall_flags_0(k,j,i), 0 ) )
[1365]	657	ENDDO
	658
	659	CASE ( 'q' )
	660
[2232]	661	DO k = nzb+1, nzt
[1382]	662	tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt) &
[2232]	663	+ fac * ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) ) * &
	664	MERGE( 1.0_wp, 0.0_wp, &
	665	BTEST( wall_flags_0(k,j,i), 0 ) )
[1365]	666	ENDDO
	667
	668	END SELECT
	669
	670	!
	671	!-- Subsidence of pt and q with prescribed profiles
	672	IF ( large_scale_subsidence .AND. use_subsidence_tendencies ) THEN
	673
	674	SELECT CASE ( prog_var )
	675
	676	CASE ( 'pt' )
	677
[2232]	678	DO k = nzb+1, nzt
[1382]	679	tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac * &
[2232]	680	( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) ) * &
	681	MERGE( 1.0_wp, 0.0_wp, &
	682	BTEST( wall_flags_0(k,j,i), 0 ) )
[1365]	683	ENDDO
	684
	685	CASE ( 'q' )
	686
[2232]	687	DO k = nzb+1, nzt
[1382]	688	tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac * &
[2232]	689	( td_sub_q(k,nt+1) - td_sub_q(k,nt) ) * &
	690	MERGE( 1.0_wp, 0.0_wp, &
	691	BTEST( wall_flags_0(k,j,i), 0 ) )
[1365]	692	ENDDO
	693
	694	END SELECT
	695
	696	ENDIF
	697
	698	END SUBROUTINE ls_advec_ij
	699
	700
[1239]	701	END MODULE ls_forcing_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |