Home

Context Navigation

source: palm/trunk/SOURCE/calc_radiation.f90 @ 1845

Last change on this file since 1845 was 1845, checked in by raasch, 9 years ago
nzb_2d replaced by nzb_..._inner, Kessler precipitation stored on surface grid point
Property svn:keywords set to `Id`
File size: 13.0 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |