Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |