Home

Context Navigation

source: palm/trunk/SOURCE/calc_radiation.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: 13.9 KB

Rev	Line
[1873]	1	!> @file calc_radiation.f90
[2000]	2	!------------------------------------------------------------------------------!
[1036]	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.
[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	!
[2101]	17	! Copyright 1997-2017 Leibniz Universitaet Hannover
[2000]	18	!------------------------------------------------------------------------------!
[1036]	19	!
[484]	20	! Current revisions:
[1]	21	! -----------------
[2232]	22	! Adjustments to new topography concept
[1354]	23	!
[1321]	24	! Former revisions:
	25	! -----------------
	26	! $Id: calc_radiation.f90 2232 2017-05-30 17:47:52Z suehring $
	27	!
[2001]	28	! 2000 2016-08-20 18:09:15Z knoop
	29	! Forced header and separation lines into 80 columns
	30	!
[1874]	31	! 1873 2016-04-18 14:50:06Z maronga
	32	! Module renamed (removed _mod)
	33	!
	34	!
[1851]	35	! 1850 2016-04-08 13:29:27Z maronga
	36	! Module renamed
	37	!
	38	!
[1683]	39	! 1682 2015-10-07 23:56:08Z knoop
	40	! Code annotations made doxygen readable
	41	!
[1354]	42	! 1353 2014-04-08 15:21:23Z heinze
	43	! REAL constants provided with KIND-attribute
	44	!
[1323]	45	! 1322 2014-03-20 16:38:49Z raasch
	46	! exponent 4.0 changed to integer
	47	!
[1321]	48	! 1320 2014-03-20 08:40:49Z raasch
[1320]	49	! ONLY-attribute added to USE-statements,
	50	! kind-parameters added to all INTEGER and REAL declaration statements,
	51	! kinds are defined in new module kinds,
	52	! revision history before 2012 removed,
	53	! comment fields (!:) to be used for variable explanations added to
	54	! all variable declaration statements
[1]	55	!
[1037]	56	! 1036 2012-10-22 13:43:42Z raasch
	57	! code put under GPL (PALM 3.9)
	58	!
[1]	59	! Revision 1.1 2000/04/13 14:42:45 schroeter
	60	! Initial revision
	61	!
	62	!
	63	! Description:
	64	! -------------
[1682]	65	!> Calculation of the vertical divergences of the long-wave radiation-fluxes
	66	!> based on the parameterization of the cloud effective emissivity
[1]	67	!------------------------------------------------------------------------------!
[1682]	68	MODULE calc_radiation_mod
	69
[1320]	70	USE kinds
	71
[1]	72	PRIVATE
	73	PUBLIC calc_radiation
	74
[1682]	75	LOGICAL, SAVE :: first_call = .TRUE. !<
	76	REAL(wp), SAVE :: sigma = 5.67E-08_wp !<
[1]	77
[1682]	78	REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE :: lwp_ground !<
	79	REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE :: lwp_top !<
	80	REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE :: blackbody_emission !<
[1]	81
	82	INTERFACE calc_radiation
	83	MODULE PROCEDURE calc_radiation
	84	MODULE PROCEDURE calc_radiation_ij
	85	END INTERFACE calc_radiation
	86
	87	CONTAINS
	88
	89
	90	!------------------------------------------------------------------------------!
[1682]	91	! Description:
	92	! ------------
	93	!> Call for all grid points
[1]	94	!------------------------------------------------------------------------------!
	95	SUBROUTINE calc_radiation
	96
[1320]	97	USE arrays_3d, &
	98	ONLY: dzw, pt, ql, tend
	99
	100	USE cloud_parameters, &
	101	ONLY: cp, l_d_cp, pt_d_t, t_d_pt
	102
	103	USE control_parameters, &
	104	ONLY: rho_surface
	105
	106	USE indices, &
[2232]	107	ONLY: nxl, nxr, nyn, nys, nzb, nzt, wall_flags_0
[1320]	108
	109	USE kinds
	110
[1]	111	USE pegrid
	112
[1320]	113
[1]	114	IMPLICIT NONE
	115
[1682]	116	INTEGER(iwp) :: i !<
	117	INTEGER(iwp) :: j !<
	118	INTEGER(iwp) :: k !<
	119	INTEGER(iwp) :: k_help !<
[1]	120
[1682]	121	REAL(wp) :: df_p !<
	122	REAL(wp) :: df_m !<
	123	REAL(wp) :: effective_emission_up_m !<
	124	REAL(wp) :: effective_emission_up_p !<
	125	REAL(wp) :: effective_emission_down_m !<
	126	REAL(wp) :: effective_emission_down_p !<
	127	REAL(wp) :: f_up_m !<
	128	REAL(wp) :: f_up_p !<
	129	REAL(wp) :: f_down_m !<
	130	REAL(wp) :: f_down_p !<
	131	REAL(wp) :: impinging_flux_at_top !<
	132	REAL(wp) :: temperature !<
[1]	133
	134
	135	!
	136	!-- On first call, allocate temporary arrays
	137	IF ( first_call ) THEN
[1320]	138	ALLOCATE( blackbody_emission(nzb:nzt+1), lwp_ground(nzb:nzt+1), &
[1]	139	lwp_top(nzb:nzt+1) )
	140	first_call = .FALSE.
	141	ENDIF
	142
	143
	144	DO i = nxl, nxr
	145	DO j = nys, nyn
	146	!
	147	!-- Compute the liquid water path (LWP) and blackbody_emission
	148	!-- at all vertical levels
[1353]	149	lwp_ground(nzb) = 0.0_wp
[1]	150	lwp_top(nzt+1) = rho_surface * ql(nzt+1,j,i) * dzw(nzt+1)
	151
	152	temperature = pt(nzb,j,i) * t_d_pt(nzb) + l_d_cp * ql(nzb,j,i)
[1322]	153	blackbody_emission(nzb) = sigma * temperature**4
[1]	154
[2232]	155	DO k = nzb+1, nzt
[1]	156
	157	k_help = ( nzt+nzb+1 ) - k
[1320]	158	lwp_ground(k) = lwp_ground(k-1) + rho_surface * ql(k,j,i) * &
[1]	159	dzw(k)
	160
[1320]	161	lwp_top(k_help) = lwp_top(k_help+1) + &
[1]	162	rho_surface * ql(k_help,j,i) * dzw(k_help)
	163
	164	temperature = pt(k,j,i) * t_d_pt(k) + l_d_cp * ql(k,j,i)
[2232]	165	blackbody_emission(k) = sigma * temperature**4 &
	166	* MERGE( 1.0_wp, 0.0_wp, &
	167	BTEST( wall_flags_0(k,j,i), 0 ) )
[1]	168
	169	ENDDO
	170
[1320]	171	lwp_ground(nzt+1) = lwp_ground(nzt) + &
[1]	172	rho_surface * ql(nzt+1,j,i) * dzw(nzt+1)
	173	lwp_top(nzb) = lwp_top(nzb+1)
	174
[1320]	175	temperature = pt(nzt+1,j,i) * t_d_pt(nzt+1) + l_d_cp * &
[1]	176	ql(nzt+1,j,i)
[1322]	177	blackbody_emission(nzt+1) = sigma * temperature**4
[1]	178
	179	!
	180	!-- See Chlond '92, this is just a first guess
[1353]	181	impinging_flux_at_top = blackbody_emission(nzb) - 100.0_wp
[1]	182
[2232]	183	DO k = nzb+1, nzt
[1]	184	!
	185	!-- Save some computational time, but this may cause load
	186	!-- imbalances if ql is not distributed uniformly
[1353]	187	IF ( ql(k,j,i) /= 0.0_wp ) THEN
[1]	188	!
	189	!-- Compute effective emissivities
[1353]	190	effective_emission_up_p = 1.0_wp - &
	191	EXP( -130.0_wp * lwp_ground(k+1) )
	192	effective_emission_up_m = 1.0_wp - &
	193	EXP( -130.0_wp * lwp_ground(k-1) )
	194	effective_emission_down_p = 1.0_wp - &
	195	EXP( -158.0_wp * lwp_top(k+1) )
	196	effective_emission_down_m = 1.0_wp - &
	197	EXP( -158.0_wp * lwp_top(k-1) )
[1]	198
	199	!
	200	!-- Compute vertical long wave radiation fluxes
[1320]	201	f_up_p = blackbody_emission(nzb) + &
	202	effective_emission_up_p * &
[1]	203	( blackbody_emission(k) - blackbody_emission(nzb) )
	204
[1320]	205	f_up_m = blackbody_emission(nzb) + &
	206	effective_emission_up_m * &
[1]	207	( blackbody_emission(k-1) - blackbody_emission(nzb) )
	208
[1320]	209	f_down_p = impinging_flux_at_top + &
	210	effective_emission_down_p * &
[1]	211	( blackbody_emission(k) - impinging_flux_at_top )
	212
[1320]	213	f_down_m = impinging_flux_at_top + &
	214	effective_emission_down_m * &
[1]	215	( blackbody_emission(k-1) - impinging_flux_at_top )
	216
	217	!
	218	!-- Divergence of vertical long wave radiation fluxes
	219	df_p = f_up_p - f_down_p
	220	df_m = f_up_m - f_down_m
	221
	222	!
	223	!-- Compute tendency term
[1320]	224	tend(k,j,i) = tend(k,j,i) - &
	225	( pt_d_t(k) / ( rho_surface * cp ) * &
[2232]	226	( df_p - df_m ) / dzw(k) ) &
	227	* MERGE( 1.0_wp, 0.0_wp, &
	228	BTEST( wall_flags_0(k,j,i), 0 ) )
[1]	229
	230	ENDIF
	231
	232	ENDDO
	233	ENDDO
	234	ENDDO
	235
	236	END SUBROUTINE calc_radiation
	237
	238
	239	!------------------------------------------------------------------------------!
[1682]	240	! Description:
	241	! ------------
	242	!> Call for grid point i,j
[1]	243	!------------------------------------------------------------------------------!
	244	SUBROUTINE calc_radiation_ij( i, j )
	245
[1320]	246	USE arrays_3d, &
	247	ONLY: dzw, pt, ql, tend
	248
	249	USE cloud_parameters, &
	250	ONLY: cp, l_d_cp, pt_d_t, t_d_pt
	251
	252	USE control_parameters, &
	253	ONLY: rho_surface
	254
	255	USE indices, &
[2232]	256	ONLY: nzb, nzt, wall_flags_0
[1320]	257
	258	USE kinds
	259
[1]	260	USE pegrid
[1320]	261
[1]	262
	263	IMPLICIT NONE
	264
[1682]	265	INTEGER(iwp) :: i !<
	266	INTEGER(iwp) :: j !<
	267	INTEGER(iwp) :: k !<
	268	INTEGER(iwp) :: k_help !<
[1]	269
[1682]	270	REAL(wp) :: df_p !<
	271	REAL(wp) :: df_m !<
	272	REAL(wp) :: effective_emission_up_m !<
	273	REAL(wp) :: effective_emission_up_p !<
	274	REAL(wp) :: effective_emission_down_m !<
	275	REAL(wp) :: effective_emission_down_p !<
	276	REAL(wp) :: f_up_m !<
	277	REAL(wp) :: f_up_p !<
	278	REAL(wp) :: f_down_m !<
	279	REAL(wp) :: f_down_p !<
	280	REAL(wp) :: impinging_flux_at_top !<
	281	REAL(wp) :: temperature !<
[1]	282
[1320]	283
[1]	284	!
	285	!-- On first call, allocate temporary arrays
	286	IF ( first_call ) THEN
[1320]	287	ALLOCATE( blackbody_emission(nzb:nzt+1), lwp_ground(nzb:nzt+1), &
[1]	288	lwp_top(nzb:nzt+1) )
	289	first_call = .FALSE.
	290	ENDIF
	291
	292	!
	293	!-- Compute the liquid water path (LWP) and blackbody_emission
	294	!-- at all vertical levels
[1353]	295	lwp_ground(nzb) = 0.0_wp
[1]	296	lwp_top(nzt+1) = rho_surface * ql(nzt+1,j,i) * dzw(nzt+1)
	297
	298	temperature = pt(nzb,j,i) * t_d_pt(nzb) + l_d_cp * ql(nzb,j,i)
[1322]	299	blackbody_emission(nzb) = sigma * temperature**4
[1]	300
[2232]	301	DO k = nzb+1, nzt
[1]	302	k_help = ( nzt+nzb+1 ) - k
	303	lwp_ground(k) = lwp_ground(k-1) + rho_surface * ql(k,j,i) * dzw(k)
	304
[1320]	305	lwp_top(k_help) = lwp_top(k_help+1) + &
[1]	306	rho_surface * ql(k_help,j,i) * dzw(k_help)
	307
	308	temperature = pt(k,j,i) * t_d_pt(k) + l_d_cp * ql(k,j,i)
[2232]	309	blackbody_emission(k) = sigma * temperature**4 &
	310	* MERGE( 1.0_wp, 0.0_wp, &
	311	BTEST( wall_flags_0(k,j,i), 0 ) )
[1]	312
	313	ENDDO
[1320]	314	lwp_ground(nzt+1) = lwp_ground(nzt) + &
[1]	315	rho_surface * ql(nzt+1,j,i) * dzw(nzt+1)
	316	lwp_top(nzb) = lwp_top(nzb+1)
	317
[1320]	318	temperature = pt(nzt+1,j,i) * t_d_pt(nzt+1) + l_d_cp * &
[1]	319	ql(nzt+1,j,i)
[1322]	320	blackbody_emission(nzt+1) = sigma * temperature**4
[1]	321
	322	!
	323	!-- See Chlond '92, this is just a first guess
[1353]	324	impinging_flux_at_top = blackbody_emission(nzb) - 100.0_wp
[1]	325
[2232]	326	DO k = nzb+1, nzt
[1]	327	!
	328	!-- Store some computational time,
	329	!-- this may cause load imbalances if ql is not distributed uniformly
[1353]	330	IF ( ql(k,j,i) /= 0.0_wp ) THEN
[1]	331	!
	332	!-- Compute effective emissivities
[1353]	333	effective_emission_up_p = 1.0_wp - &
	334	EXP( -130.0_wp * lwp_ground(k+1) )
	335	effective_emission_up_m = 1.0_wp - &
	336	EXP( -130.0_wp * lwp_ground(k-1) )
	337	effective_emission_down_p = 1.0_wp - &
	338	EXP( -158.0_wp * lwp_top(k+1) )
	339	effective_emission_down_m = 1.0_wp - &
	340	EXP( -158.0_wp * lwp_top(k-1) )
[1]	341
	342	!
	343	!-- Compute vertical long wave radiation fluxes
[1320]	344	f_up_p = blackbody_emission(nzb) + effective_emission_up_p * &
[1]	345	( blackbody_emission(k) - blackbody_emission(nzb) )
	346
[1320]	347	f_up_m = blackbody_emission(nzb) + effective_emission_up_m * &
[1]	348	( blackbody_emission(k-1) - blackbody_emission(nzb) )
	349
[1320]	350	f_down_p = impinging_flux_at_top + effective_emission_down_p * &
[1]	351	( blackbody_emission(k) - impinging_flux_at_top )
	352
[1320]	353	f_down_m = impinging_flux_at_top + effective_emission_down_m * &
[1]	354	( blackbody_emission(k-1) - impinging_flux_at_top )
	355
	356	!
	357	!- Divergence of vertical long wave radiation fluxes
	358	df_p = f_up_p - f_down_p
	359	df_m = f_up_m - f_down_m
	360
	361	!
	362	!-- Compute tendency term
[1320]	363	tend(k,j,i) = tend(k,j,i) - ( pt_d_t(k) / ( rho_surface * cp ) * &
[2232]	364	( df_p - df_m ) / dzw(k) ) &
	365	* MERGE( 1.0_wp, 0.0_wp, &
	366	BTEST( wall_flags_0(k,j,i), 0 ) )
[1]	367
	368	ENDIF
	369
	370	ENDDO
	371
	372	END SUBROUTINE calc_radiation_ij
	373
	374	END MODULE calc_radiation_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |