1 | !> @file basic_constants_and_equations_mod.f90 |
---|
2 | !--------------------------------------------------------------------------------------------------! |
---|
3 | ! This file is part of the PALM model system. |
---|
4 | ! |
---|
5 | ! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General |
---|
6 | ! Public License as published by the Free Software Foundation, either version 3 of the License, or |
---|
7 | ! (at your option) any later version. |
---|
8 | ! |
---|
9 | ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the |
---|
10 | ! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General |
---|
11 | ! Public License for more details. |
---|
12 | ! |
---|
13 | ! You should have received a copy of the GNU General Public License along with PALM. If not, see |
---|
14 | ! <http://www.gnu.org/licenses/>. |
---|
15 | ! |
---|
16 | ! Copyright 1997-2020 Leibniz Universitaet Hannover |
---|
17 | !--------------------------------------------------------------------------------------------------! |
---|
18 | ! |
---|
19 | ! Current revisions: |
---|
20 | ! ----------------- |
---|
21 | ! |
---|
22 | ! |
---|
23 | ! Former revisions: |
---|
24 | ! ----------------- |
---|
25 | ! $Id: basic_constants_and_equations_mod.f90 4742 2020-10-14 15:11:02Z moh.hefny $ |
---|
26 | ! Implement snow and graupel (bulk microphysics) |
---|
27 | ! |
---|
28 | ! 4509 2020-04-26 15:57:55Z raasch |
---|
29 | ! file re-formatted to follow the PALM coding standard |
---|
30 | ! |
---|
31 | ! 4502 2020-04-17 16:14:16Z schwenkel |
---|
32 | ! Implementation of ice microphysics |
---|
33 | ! |
---|
34 | ! 4400 2020-02-10 20:32:41Z suehring |
---|
35 | ! Move routine to transform coordinates from netcdf_interface_mod to |
---|
36 | ! basic_constants_and_equations_mod |
---|
37 | ! |
---|
38 | ! 4360 2020-01-07 11:25:50Z suehring |
---|
39 | ! Corrected "Former revisions" section |
---|
40 | ! |
---|
41 | ! 4088 2019-07-11 13:57:56Z Giersch |
---|
42 | ! Comment of barometric formula improved, function for ideal gas law revised |
---|
43 | ! |
---|
44 | ! 4084 2019-07-10 17:09:11Z knoop |
---|
45 | ! Changed precomputed fractions to be variable based |
---|
46 | ! |
---|
47 | ! 4055 2019-06-27 09:47:29Z suehring |
---|
48 | ! Added rgas_univ (universal gas constant) (E.C. Chan) |
---|
49 | ! |
---|
50 | ! |
---|
51 | ! 3655 2019-01-07 16:51:22Z knoop |
---|
52 | ! OpenACC port for SPEC |
---|
53 | ! 3361 2018-10-16 20:39:37Z knoop |
---|
54 | ! New module (introduced with modularization of bulk cloud physics model) |
---|
55 | ! |
---|
56 | ! |
---|
57 | ! |
---|
58 | ! |
---|
59 | ! Description: |
---|
60 | ! ------------ |
---|
61 | !> This module contains all basic (physical) constants and functions for the calculation of |
---|
62 | !> diagnostic quantities. |
---|
63 | !- -----------------------------------------------------------------------------! |
---|
64 | MODULE basic_constants_and_equations_mod |
---|
65 | |
---|
66 | |
---|
67 | USE kinds |
---|
68 | |
---|
69 | IMPLICIT NONE |
---|
70 | |
---|
71 | |
---|
72 | REAL(wp), PARAMETER :: c_p = 1005.0_wp !< heat capacity of dry air (J kg-1 K-1) |
---|
73 | REAL(wp), PARAMETER :: c_w = 4185.0_wp !< heat capacity of water at 0°C (J kg-1 K-1) |
---|
74 | REAL(wp), PARAMETER :: degc_to_k = 273.15_wp !< temperature (in K) of 0 deg C (K) |
---|
75 | REAL(wp), PARAMETER :: g = 9.81_wp !< gravitational acceleration (m s-2) |
---|
76 | REAL(wp), PARAMETER :: kappa = 0.4_wp !< von Karman constant |
---|
77 | REAL(wp), PARAMETER :: l_m = 0.33E+06_wp !< latent heat of water melting (J kg-1) |
---|
78 | REAL(wp), PARAMETER :: l_v = 2.5E+06_wp !< latent heat of water vaporization (J kg-1) |
---|
79 | REAL(wp), PARAMETER :: l_s = l_m + l_v !< latent heat of water sublimation (J kg-1) |
---|
80 | REAL(wp), PARAMETER :: molecular_weight_of_nacl = 0.05844_wp !< mol. m. NaCl (kg mol-1) |
---|
81 | REAL(wp), PARAMETER :: molecular_weight_of_c3h4o4 = 0.10406_wp !< mol. m. malonic acid (kg mol-1) |
---|
82 | REAL(wp), PARAMETER :: molecular_weight_of_nh4no3 = 0.08004_wp !< mol. m. ammonium sulfate (kg mol-1) |
---|
83 | REAL(wp), PARAMETER :: molecular_weight_of_water = 0.01801528_wp !< mol. m. H2O (kg mol-1) |
---|
84 | REAL(wp), PARAMETER :: pi = 3.141592654_wp !< PI |
---|
85 | !$ACC DECLARE COPYIN(pi) |
---|
86 | REAL(wp), PARAMETER :: rgas_univ = 8.31446261815324_wp !< universal gas constant (J K-1 mol-1) |
---|
87 | REAL(wp), PARAMETER :: rho_i = 916.7_wp !> density of pure ice (kg m-3) |
---|
88 | REAL(wp), PARAMETER :: rho_l = 1.0E3_wp !< density of water (kg m-3) |
---|
89 | REAL(wp), PARAMETER :: rho_nacl = 2165.0_wp !< density of NaCl (kg m-3) |
---|
90 | REAL(wp), PARAMETER :: rho_c3h4o4 = 1600.0_wp !< density of malonic acid (kg m-3) |
---|
91 | REAL(wp), PARAMETER :: rho_nh4no3 = 1720.0_wp !< density of ammonium sulfate (kg m-3) |
---|
92 | REAL(wp), PARAMETER :: r_d = 287.0_wp !< sp. gas const. dry air (J kg-1 K-1) |
---|
93 | REAL(wp), PARAMETER :: r_v = 461.51_wp !< sp. gas const. water vapor (J kg-1 K-1) |
---|
94 | REAL(wp), PARAMETER :: sigma_sb = 5.67037E-08_wp !< Stefan-Boltzmann constant |
---|
95 | REAL(wp), PARAMETER :: solar_constant = 1368.0_wp !< solar constant at top of atmosphere |
---|
96 | REAL(wp), PARAMETER :: vanthoff_nacl = 2.0_wp !< van't Hoff factor for NaCl |
---|
97 | REAL(wp), PARAMETER :: vanthoff_c3h4o4 = 1.37_wp !< van't Hoff factor for malonic acid |
---|
98 | REAL(wp), PARAMETER :: vanthoff_nh4no3 = 2.31_wp !< van't Hoff factor for ammonium sulfate |
---|
99 | |
---|
100 | REAL(wp), PARAMETER :: p_0 = 100000.0_wp !< standard pressure reference state |
---|
101 | |
---|
102 | REAL(wp), PARAMETER :: cp_d_rd = c_p / r_d !< precomputed c_p / r_d |
---|
103 | REAL(wp), PARAMETER :: g_d_cp = g / c_p !< precomputed g / c_p |
---|
104 | REAL(wp), PARAMETER :: lv_d_cp = l_v / c_p !< precomputed l_v / c_p |
---|
105 | REAL(wp), PARAMETER :: ls_d_cp = l_s / c_p !< precomputed l_s / c_p |
---|
106 | REAL(wp), PARAMETER :: lv_d_rd = l_v / r_d !< precomputed l_v / r_d |
---|
107 | REAL(wp), PARAMETER :: rd_d_rv = r_d / r_v !< precomputed r_d / r_v |
---|
108 | REAL(wp), PARAMETER :: rd_d_cp = r_d / c_p !< precomputed r_d / c_p |
---|
109 | |
---|
110 | REAL(wp) :: molecular_weight_of_solute = molecular_weight_of_nacl !< mol. m. NaCl (kg mol-1) |
---|
111 | REAL(wp) :: rho_s = rho_nacl !< density of NaCl (kg m-3) |
---|
112 | REAL(wp) :: vanthoff = vanthoff_nacl !< van't Hoff factor for NaCl |
---|
113 | |
---|
114 | SAVE |
---|
115 | |
---|
116 | PRIVATE magnus_0d, & |
---|
117 | magnus_1d, & |
---|
118 | magnus_tl_0d, & |
---|
119 | magnus_tl_1d, & |
---|
120 | magnus_0d_ice, & |
---|
121 | magnus_1d_ice, & |
---|
122 | ideal_gas_law_rho_0d, & |
---|
123 | ideal_gas_law_rho_1d, & |
---|
124 | ideal_gas_law_rho_pt_0d, & |
---|
125 | ideal_gas_law_rho_pt_1d, & |
---|
126 | exner_function_0d, & |
---|
127 | exner_function_1d, & |
---|
128 | exner_function_invers_0d, & |
---|
129 | exner_function_invers_1d, & |
---|
130 | barometric_formula_0d, & |
---|
131 | barometric_formula_1d |
---|
132 | |
---|
133 | |
---|
134 | INTERFACE convert_utm_to_geographic |
---|
135 | MODULE PROCEDURE convert_utm_to_geographic |
---|
136 | END INTERFACE convert_utm_to_geographic |
---|
137 | |
---|
138 | INTERFACE magnus |
---|
139 | MODULE PROCEDURE magnus_0d |
---|
140 | MODULE PROCEDURE magnus_1d |
---|
141 | END INTERFACE magnus |
---|
142 | |
---|
143 | INTERFACE magnus_tl |
---|
144 | MODULE PROCEDURE magnus_tl_0d |
---|
145 | MODULE PROCEDURE magnus_tl_1d |
---|
146 | END INTERFACE magnus_tl |
---|
147 | |
---|
148 | INTERFACE magnus_ice |
---|
149 | MODULE PROCEDURE magnus_0d_ice |
---|
150 | MODULE PROCEDURE magnus_1d_ice |
---|
151 | END INTERFACE magnus_ice |
---|
152 | |
---|
153 | INTERFACE ideal_gas_law_rho |
---|
154 | MODULE PROCEDURE ideal_gas_law_rho_0d |
---|
155 | MODULE PROCEDURE ideal_gas_law_rho_1d |
---|
156 | END INTERFACE ideal_gas_law_rho |
---|
157 | |
---|
158 | INTERFACE ideal_gas_law_rho_pt |
---|
159 | MODULE PROCEDURE ideal_gas_law_rho_pt_0d |
---|
160 | MODULE PROCEDURE ideal_gas_law_rho_pt_1d |
---|
161 | END INTERFACE ideal_gas_law_rho_pt |
---|
162 | |
---|
163 | INTERFACE exner_function |
---|
164 | MODULE PROCEDURE exner_function_0d |
---|
165 | MODULE PROCEDURE exner_function_1d |
---|
166 | END INTERFACE exner_function |
---|
167 | |
---|
168 | INTERFACE exner_function_invers |
---|
169 | MODULE PROCEDURE exner_function_invers_0d |
---|
170 | MODULE PROCEDURE exner_function_invers_1d |
---|
171 | END INTERFACE exner_function_invers |
---|
172 | |
---|
173 | INTERFACE barometric_formula |
---|
174 | MODULE PROCEDURE barometric_formula_0d |
---|
175 | MODULE PROCEDURE barometric_formula_1d |
---|
176 | END INTERFACE barometric_formula |
---|
177 | ! |
---|
178 | !-- Public routines |
---|
179 | PUBLIC convert_utm_to_geographic |
---|
180 | |
---|
181 | CONTAINS |
---|
182 | |
---|
183 | |
---|
184 | !--------------------------------------------------------------------------------------------------! |
---|
185 | ! Description: |
---|
186 | ! ------------ |
---|
187 | !> Convert UTM coordinates into geographic latitude and longitude. Conversion is based on the work |
---|
188 | !> of KrÃŒger (1912) DOI: 10.2312/GFZ.b103-krueger28 and Karney (2013) DOI: 10.1007/s00190-012-0578-z |
---|
189 | !> Based on a JavaScript of the geodesy function library written by chrisveness |
---|
190 | !> https://github.com/chrisveness/geodesy |
---|
191 | !--------------------------------------------------------------------------------------------------! |
---|
192 | SUBROUTINE convert_utm_to_geographic( crs, eutm, nutm, lon, lat ) |
---|
193 | |
---|
194 | INTEGER(iwp) :: j !< loop index |
---|
195 | |
---|
196 | REAL(wp), INTENT(in) :: eutm !< easting (UTM) |
---|
197 | REAL(wp), INTENT(out) :: lat !< geographic latitude in degree |
---|
198 | REAL(wp), INTENT(out) :: lon !< geographic longitude in degree |
---|
199 | REAL(wp), INTENT(in) :: nutm !< northing (UTM) |
---|
200 | |
---|
201 | REAL(wp) :: a !< 2*pi*a is the circumference of a meridian |
---|
202 | REAL(wp) :: cos_eta_s !< cos(eta_s) |
---|
203 | REAL(wp) :: delta_i !< |
---|
204 | REAL(wp) :: delta_tau_i !< |
---|
205 | REAL(wp) :: e !< eccentricity |
---|
206 | REAL(wp) :: eta !< |
---|
207 | REAL(wp) :: eta_s !< |
---|
208 | REAL(wp) :: n !< 3rd flattening |
---|
209 | REAL(wp) :: n2 !< n^2 |
---|
210 | REAL(wp) :: n3 !< n^3 |
---|
211 | REAL(wp) :: n4 !< n^4 |
---|
212 | REAL(wp) :: n5 !< n^5 |
---|
213 | REAL(wp) :: n6 !< n^6 |
---|
214 | REAL(wp) :: nu !< |
---|
215 | REAL(wp) :: nu_s !< |
---|
216 | REAL(wp) :: sin_eta_s !< sin(eta_s) |
---|
217 | REAL(wp) :: sinh_nu_s !< sinush(nu_s) |
---|
218 | REAL(wp) :: tau_i !< |
---|
219 | REAL(wp) :: tau_i_s !< |
---|
220 | REAL(wp) :: tau_s !< |
---|
221 | REAL(wp) :: x !< adjusted easting |
---|
222 | REAL(wp) :: y !< adjusted northing |
---|
223 | |
---|
224 | REAL(wp), DIMENSION(6) :: beta !< 6th order KrÃŒger expressions |
---|
225 | |
---|
226 | REAL(wp), DIMENSION(8), INTENT(in) :: crs !< coordinate reference system, consists of |
---|
227 | !< (/semi_major_axis, |
---|
228 | !< inverse_flattening, |
---|
229 | !< longitude_of_prime_meridian, |
---|
230 | !< longitude_of_central_meridian, |
---|
231 | !< scale_factor_at_central_meridian, |
---|
232 | !< latitude_of_projection_origin, |
---|
233 | !< false_easting, |
---|
234 | !< false_northing /) |
---|
235 | |
---|
236 | x = eutm - crs(7) ! remove false easting |
---|
237 | y = nutm - crs(8) ! remove false northing |
---|
238 | ! |
---|
239 | !-- From Karney 2011 Eq 15-22, 36: |
---|
240 | e = SQRT( 1.0_wp / crs(2) * ( 2.0_wp - 1.0_wp / crs(2) ) ) |
---|
241 | n = 1.0_wp / crs(2) / ( 2.0_wp - 1.0_wp / crs(2) ) |
---|
242 | n2 = n * n |
---|
243 | n3 = n * n2 |
---|
244 | n4 = n * n3 |
---|
245 | n5 = n * n4 |
---|
246 | n6 = n * n5 |
---|
247 | |
---|
248 | a = crs(1) / ( 1.0_wp + n ) * ( 1.0_wp + 0.25_wp * n2 + 0.015625_wp * n4 + 3.90625E-3_wp * n6 ) |
---|
249 | |
---|
250 | nu = x / ( crs(5) * a ) |
---|
251 | eta = y / ( crs(5) * a ) |
---|
252 | |
---|
253 | !-- According to KrÃŒger (1912), eq. 26* |
---|
254 | beta(1) = 0.5_wp * n & |
---|
255 | - 2.0_wp / 3.0_wp * n2 & |
---|
256 | + 37.0_wp / 96.0_wp * n3 & |
---|
257 | - 1.0_wp / 360.0_wp * n4 & |
---|
258 | - 81.0_wp / 512.0_wp * n5 & |
---|
259 | + 96199.0_wp / 604800.0_wp * n6 |
---|
260 | |
---|
261 | beta(2) = 1.0_wp / 48.0_wp * n2 & |
---|
262 | + 1.0_wp / 15.0_wp * n3 & |
---|
263 | - 437.0_wp / 1440.0_wp * n4 & |
---|
264 | + 46.0_wp / 105.0_wp * n5 & |
---|
265 | - 1118711.0_wp / 3870720.0_wp * n6 |
---|
266 | |
---|
267 | beta(3) = 17.0_wp / 480.0_wp * n3 & |
---|
268 | - 37.0_wp / 840.0_wp * n4 & |
---|
269 | - 209.0_wp / 4480.0_wp * n5 & |
---|
270 | + 5569.0_wp / 90720.0_wp * n6 |
---|
271 | |
---|
272 | beta(4) = 4397.0_wp / 161280.0_wp * n4 & |
---|
273 | - 11.0_wp / 504.0_wp * n5 & |
---|
274 | - 830251.0_wp / 7257600.0_wp * n6 |
---|
275 | |
---|
276 | beta(5) = 4583.0_wp / 161280.0_wp * n5 & |
---|
277 | - 108847.0_wp / 3991680.0_wp * n6 |
---|
278 | |
---|
279 | beta(6) = 20648693.0_wp / 638668800.0_wp * n6 |
---|
280 | |
---|
281 | eta_s = eta |
---|
282 | nu_s = nu |
---|
283 | DO j = 1, 6 |
---|
284 | eta_s = eta_s - beta(j) * SIN(2.0_wp * j * eta) * COSH(2.0_wp * j * nu) |
---|
285 | nu_s = nu_s - beta(j) * COS(2.0_wp * j * eta) * SINH(2.0_wp * j * nu) |
---|
286 | ENDDO |
---|
287 | |
---|
288 | sinh_nu_s = SINH( nu_s ) |
---|
289 | sin_eta_s = SIN( eta_s ) |
---|
290 | cos_eta_s = COS( eta_s ) |
---|
291 | |
---|
292 | tau_s = sin_eta_s / SQRT( sinh_nu_s**2 + cos_eta_s**2 ) |
---|
293 | |
---|
294 | tau_i = tau_s |
---|
295 | delta_tau_i = 1.0_wp |
---|
296 | |
---|
297 | DO WHILE ( ABS( delta_tau_i ) > 1.0E-12_wp ) |
---|
298 | |
---|
299 | delta_i = SINH( e * ATANH( e * tau_i / SQRT( 1.0_wp + tau_i**2 ) ) ) |
---|
300 | |
---|
301 | tau_i_s = tau_i * SQRT( 1.0_wp + delta_i**2 ) - delta_i * SQRT( 1.0_wp + tau_i**2 ) |
---|
302 | |
---|
303 | delta_tau_i = ( tau_s - tau_i_s ) / SQRT( 1.0_wp + tau_i_s**2 ) & |
---|
304 | * ( 1.0_wp + ( 1.0_wp - e**2 ) * tau_i**2 ) & |
---|
305 | / ( ( 1.0_wp - e**2 ) * SQRT( 1.0_wp + tau_i**2 ) ) |
---|
306 | |
---|
307 | tau_i = tau_i + delta_tau_i |
---|
308 | |
---|
309 | ENDDO |
---|
310 | |
---|
311 | lat = ATAN( tau_i ) / pi * 180.0_wp |
---|
312 | lon = ATAN2( sinh_nu_s, cos_eta_s ) / pi * 180.0_wp + crs(4) |
---|
313 | |
---|
314 | END SUBROUTINE convert_utm_to_geographic |
---|
315 | |
---|
316 | !--------------------------------------------------------------------------------------------------! |
---|
317 | ! Description: |
---|
318 | ! ------------ |
---|
319 | !> This function computes the magnus formula (Press et al., 1992). |
---|
320 | !> The magnus formula is needed to calculate the saturation vapor pressure. |
---|
321 | !--------------------------------------------------------------------------------------------------! |
---|
322 | FUNCTION magnus_0d( t ) |
---|
323 | |
---|
324 | IMPLICIT NONE |
---|
325 | |
---|
326 | REAL(wp), INTENT(IN) :: t !< temperature (K) |
---|
327 | |
---|
328 | REAL(wp) :: magnus_0d |
---|
329 | |
---|
330 | ! |
---|
331 | !-- Saturation vapor pressure for a specific temperature: |
---|
332 | magnus_0d = 611.2_wp * EXP( 17.62_wp * ( t - degc_to_k ) / ( t - 29.65_wp ) ) |
---|
333 | |
---|
334 | END FUNCTION magnus_0d |
---|
335 | |
---|
336 | !--------------------------------------------------------------------------------------------------! |
---|
337 | ! Description: |
---|
338 | ! ------------ |
---|
339 | !> This function computes the magnus formula (Press et al., 1992). |
---|
340 | !> The magnus formula is needed to calculate the saturation vapor pressure. |
---|
341 | !--------------------------------------------------------------------------------------------------! |
---|
342 | FUNCTION magnus_1d( t ) |
---|
343 | |
---|
344 | IMPLICIT NONE |
---|
345 | |
---|
346 | REAL(wp), INTENT(IN), DIMENSION(:) :: t !< temperature (K) |
---|
347 | |
---|
348 | REAL(wp), DIMENSION(size(t)) :: magnus_1d |
---|
349 | |
---|
350 | ! |
---|
351 | !-- Saturation vapor pressure for a specific temperature: |
---|
352 | magnus_1d = 611.2_wp * EXP( 17.62_wp * ( t - degc_to_k ) / ( t - 29.65_wp ) ) |
---|
353 | |
---|
354 | END FUNCTION magnus_1d |
---|
355 | |
---|
356 | !--------------------------------------------------------------------------------------------------! |
---|
357 | ! Description: |
---|
358 | ! ------------ |
---|
359 | !> This function computes the magnus formula (Press et al., 1992) using the (ice-) liquid water |
---|
360 | !> potential temperature. |
---|
361 | !> The magnus formula is needed to calculate the saturation vapor pressure over a plane liquid water |
---|
362 | !> surface. |
---|
363 | !--------------------------------------------------------------------------------------------------! |
---|
364 | FUNCTION magnus_tl_0d( t_l ) |
---|
365 | |
---|
366 | IMPLICIT NONE |
---|
367 | |
---|
368 | REAL(wp), INTENT(IN) :: t_l !< liquid water temperature (K) |
---|
369 | |
---|
370 | REAL(wp) :: magnus_tl_0d |
---|
371 | |
---|
372 | ! |
---|
373 | !-- Saturation vapor pressure for a specific temperature: |
---|
374 | magnus_tl_0d = 610.78_wp * EXP( 17.269_wp * ( t_l - 273.16_wp ) / ( t_l - 35.86_wp ) ) |
---|
375 | |
---|
376 | END FUNCTION magnus_tl_0d |
---|
377 | |
---|
378 | !--------------------------------------------------------------------------------------------------! |
---|
379 | ! Description: |
---|
380 | ! ------------ |
---|
381 | !> This function computes the magnus formula (Press et al., 1992) using the (ice-) liquid water |
---|
382 | !> potential temperature. |
---|
383 | !> The magnus formula is needed to calculate the saturation vapor pressure over a plane liquid water |
---|
384 | !> surface. |
---|
385 | !--------------------------------------------------------------------------------------------------! |
---|
386 | FUNCTION magnus_tl_1d( t_l ) |
---|
387 | |
---|
388 | IMPLICIT NONE |
---|
389 | |
---|
390 | REAL(wp), INTENT(IN), DIMENSION(:) :: t_l !< liquid water temperature (K) |
---|
391 | |
---|
392 | REAL(wp), DIMENSION(size(t_l)) :: magnus_tl_1d |
---|
393 | ! |
---|
394 | !-- Saturation vapor pressure for a specific temperature: |
---|
395 | magnus_tl_1d = 610.78_wp * EXP( 17.269_wp * ( t_l - 273.16_wp ) / ( t_l - 35.86_wp ) ) |
---|
396 | |
---|
397 | END FUNCTION magnus_tl_1d |
---|
398 | |
---|
399 | !--------------------------------------------------------------------------------------------------! |
---|
400 | ! Description: |
---|
401 | ! ------------ |
---|
402 | !> This function computes the magnus formula (Press et al., 1992). |
---|
403 | !> The magnus formula is needed to calculate the saturation vapor pressure over a plane ice surface. |
---|
404 | !--------------------------------------------------------------------------------------------------! |
---|
405 | FUNCTION magnus_0d_ice( t ) |
---|
406 | |
---|
407 | IMPLICIT NONE |
---|
408 | |
---|
409 | REAL(wp), INTENT(IN) :: t !< temperature (K) |
---|
410 | |
---|
411 | REAL(wp) :: magnus_0d_ice |
---|
412 | |
---|
413 | ! |
---|
414 | !-- Saturation vapor pressure for a specific temperature: |
---|
415 | !magnus_0d_ice = 611.2_wp * EXP( 22.46_wp * ( t - degc_to_k ) / ( t - 0.53_wp ) ) |
---|
416 | magnus_0d_ice = 610.78_wp * EXP( 21.875_wp * ( t - degc_to_k ) / ( t - 7.66_wp ) ) |
---|
417 | |
---|
418 | |
---|
419 | END FUNCTION magnus_0d_ice |
---|
420 | |
---|
421 | !--------------------------------------------------------------------------------------------------! |
---|
422 | ! Description: |
---|
423 | ! ------------ |
---|
424 | !> This function computes the magnus formula (Press et al., 1992). |
---|
425 | !> The magnus formula is needed to calculate the saturation vapor pressure over a plane ice surface. |
---|
426 | !--------------------------------------------------------------------------------------------------! |
---|
427 | FUNCTION magnus_1d_ice( t ) |
---|
428 | |
---|
429 | IMPLICIT NONE |
---|
430 | |
---|
431 | REAL(wp), INTENT(IN), DIMENSION(:) :: t !< temperature (K) |
---|
432 | |
---|
433 | REAL(wp), DIMENSION(size(t)) :: magnus_1d_ice |
---|
434 | |
---|
435 | ! |
---|
436 | !-- Saturation vapor pressure for a specific temperature: |
---|
437 | !magnus_1d_ice = 611.2_wp * EXP( 22.46_wp * ( t - degc_to_k ) / ( t - 0.53_wp ) ) |
---|
438 | magnus_1d_ice = 610.78_wp * EXP( 21.875_wp * ( t - degc_to_k ) / ( t - 7.66_wp ) ) |
---|
439 | |
---|
440 | |
---|
441 | END FUNCTION magnus_1d_ice |
---|
442 | |
---|
443 | !--------------------------------------------------------------------------------------------------! |
---|
444 | ! Description: |
---|
445 | ! ------------ |
---|
446 | !> Compute the ideal gas law for scalar arguments. |
---|
447 | !--------------------------------------------------------------------------------------------------! |
---|
448 | FUNCTION ideal_gas_law_rho_0d( p, t ) |
---|
449 | |
---|
450 | IMPLICIT NONE |
---|
451 | |
---|
452 | REAL(wp), INTENT(IN) :: p !< pressure (Pa) |
---|
453 | REAL(wp), INTENT(IN) :: t !< temperature (K) |
---|
454 | |
---|
455 | REAL(wp) :: ideal_gas_law_rho_0d |
---|
456 | |
---|
457 | ! |
---|
458 | !-- Compute density according to ideal gas law: |
---|
459 | ideal_gas_law_rho_0d = p / (r_d * t) |
---|
460 | |
---|
461 | END FUNCTION ideal_gas_law_rho_0d |
---|
462 | |
---|
463 | !--------------------------------------------------------------------------------------------------! |
---|
464 | ! Description: |
---|
465 | ! ------------ |
---|
466 | !> Compute the ideal gas law for 1-D array arguments. |
---|
467 | !--------------------------------------------------------------------------------------------------! |
---|
468 | FUNCTION ideal_gas_law_rho_1d( p, t ) |
---|
469 | |
---|
470 | IMPLICIT NONE |
---|
471 | |
---|
472 | REAL(wp), INTENT(IN), DIMENSION(:) :: p !< pressure (Pa) |
---|
473 | REAL(wp), INTENT(IN), DIMENSION(:) :: t !< temperature (K) |
---|
474 | |
---|
475 | REAL(wp), DIMENSION(size(p)) :: ideal_gas_law_rho_1d |
---|
476 | |
---|
477 | ! |
---|
478 | !-- Compute density according to ideal gas law: |
---|
479 | ideal_gas_law_rho_1d = p / (r_d * t) |
---|
480 | |
---|
481 | END FUNCTION ideal_gas_law_rho_1d |
---|
482 | |
---|
483 | !--------------------------------------------------------------------------------------------------! |
---|
484 | ! Description: |
---|
485 | ! ------------ |
---|
486 | !> Compute the ideal gas law for scalar arguments. |
---|
487 | !--------------------------------------------------------------------------------------------------! |
---|
488 | FUNCTION ideal_gas_law_rho_pt_0d( p, t ) |
---|
489 | |
---|
490 | IMPLICIT NONE |
---|
491 | |
---|
492 | REAL(wp), INTENT(IN) :: p !< pressure (Pa) |
---|
493 | REAL(wp), INTENT(IN) :: t !< temperature (K) |
---|
494 | |
---|
495 | REAL(wp) :: ideal_gas_law_rho_pt_0d |
---|
496 | |
---|
497 | ! |
---|
498 | !-- Compute density according to ideal gas law: |
---|
499 | ideal_gas_law_rho_pt_0d = p / (r_d * exner_function(p) * t) |
---|
500 | |
---|
501 | END FUNCTION ideal_gas_law_rho_pt_0d |
---|
502 | |
---|
503 | !--------------------------------------------------------------------------------------------------! |
---|
504 | ! Description: |
---|
505 | ! ------------ |
---|
506 | !> Compute the ideal gas law for 1-D array arguments. |
---|
507 | !--------------------------------------------------------------------------------------------------! |
---|
508 | FUNCTION ideal_gas_law_rho_pt_1d( p, t ) |
---|
509 | |
---|
510 | IMPLICIT NONE |
---|
511 | |
---|
512 | REAL(wp), INTENT(IN), DIMENSION(:) :: p !< pressure (Pa) |
---|
513 | REAL(wp), INTENT(IN), DIMENSION(:) :: t !< temperature (K) |
---|
514 | |
---|
515 | REAL(wp), DIMENSION(size(p)) :: ideal_gas_law_rho_pt_1d |
---|
516 | |
---|
517 | ! |
---|
518 | !-- Compute density according to ideal gas law: |
---|
519 | ideal_gas_law_rho_pt_1d = p / (r_d * exner_function(p) * t) |
---|
520 | |
---|
521 | END FUNCTION ideal_gas_law_rho_pt_1d |
---|
522 | |
---|
523 | !--------------------------------------------------------------------------------------------------! |
---|
524 | ! Description: |
---|
525 | ! ------------ |
---|
526 | !> Compute the exner function for scalar arguments. |
---|
527 | !--------------------------------------------------------------------------------------------------! |
---|
528 | FUNCTION exner_function_0d( p ) |
---|
529 | |
---|
530 | IMPLICIT NONE |
---|
531 | |
---|
532 | REAL(wp), INTENT(IN) :: p !< pressure (Pa) |
---|
533 | |
---|
534 | REAL(wp) :: exner_function_0d |
---|
535 | |
---|
536 | ! |
---|
537 | !-- Compute exner function: |
---|
538 | exner_function_0d = ( p / p_0 )**( rd_d_cp ) |
---|
539 | |
---|
540 | END FUNCTION exner_function_0d |
---|
541 | |
---|
542 | !--------------------------------------------------------------------------------------------------! |
---|
543 | ! Description: |
---|
544 | ! ------------ |
---|
545 | !> Compute the exner function for 1-D array arguments. |
---|
546 | !--------------------------------------------------------------------------------------------------! |
---|
547 | FUNCTION exner_function_1d( p ) |
---|
548 | |
---|
549 | IMPLICIT NONE |
---|
550 | |
---|
551 | REAL(wp), INTENT(IN), DIMENSION(:) :: p !< pressure (Pa) |
---|
552 | |
---|
553 | REAL(wp), DIMENSION(size(p)) :: exner_function_1d |
---|
554 | |
---|
555 | ! |
---|
556 | !-- Compute exner function: |
---|
557 | exner_function_1d = ( p / p_0 )**( rd_d_cp ) |
---|
558 | |
---|
559 | END FUNCTION exner_function_1d |
---|
560 | |
---|
561 | !--------------------------------------------------------------------------------------------------! |
---|
562 | ! Description: |
---|
563 | ! ------------ |
---|
564 | !> Compute the exner function for scalar arguments. |
---|
565 | !--------------------------------------------------------------------------------------------------! |
---|
566 | FUNCTION exner_function_invers_0d( p ) |
---|
567 | |
---|
568 | IMPLICIT NONE |
---|
569 | |
---|
570 | REAL(wp), INTENT(IN) :: p !< pressure (Pa) |
---|
571 | |
---|
572 | REAL(wp) :: exner_function_invers_0d |
---|
573 | |
---|
574 | ! |
---|
575 | !-- Compute exner function: |
---|
576 | exner_function_invers_0d = ( p_0 / p )**( rd_d_cp ) |
---|
577 | |
---|
578 | END FUNCTION exner_function_invers_0d |
---|
579 | |
---|
580 | !--------------------------------------------------------------------------------------------------! |
---|
581 | ! Description: |
---|
582 | ! ------------ |
---|
583 | !> Compute the exner function for 1-D array arguments. |
---|
584 | !--------------------------------------------------------------------------------------------------! |
---|
585 | FUNCTION exner_function_invers_1d( p ) |
---|
586 | |
---|
587 | IMPLICIT NONE |
---|
588 | |
---|
589 | REAL(wp), INTENT(IN), DIMENSION(:) :: p !< pressure (Pa) |
---|
590 | |
---|
591 | REAL(wp), DIMENSION(size(p)) :: exner_function_invers_1d |
---|
592 | |
---|
593 | ! |
---|
594 | !-- Compute exner function: |
---|
595 | exner_function_invers_1d = ( p_0 / p )**( rd_d_cp ) |
---|
596 | |
---|
597 | END FUNCTION exner_function_invers_1d |
---|
598 | |
---|
599 | !--------------------------------------------------------------------------------------------------! |
---|
600 | ! Description: |
---|
601 | ! ------------ |
---|
602 | !> Compute the barometric formula for scalar arguments. The calculation is based on the assumption |
---|
603 | !> of a polytropic atmosphere and neutral stratification, where the temperature lapse rate is g/cp. |
---|
604 | !--------------------------------------------------------------------------------------------------! |
---|
605 | FUNCTION barometric_formula_0d( z, t_0, p_0) |
---|
606 | |
---|
607 | IMPLICIT NONE |
---|
608 | |
---|
609 | REAL(wp), INTENT(IN) :: z !< height (m) |
---|
610 | REAL(wp), INTENT(IN) :: t_0 !< temperature reference state (K) |
---|
611 | REAL(wp), INTENT(IN) :: p_0 !< surface pressure (Pa) |
---|
612 | |
---|
613 | REAL(wp) :: barometric_formula_0d |
---|
614 | |
---|
615 | ! |
---|
616 | !-- Compute barometric formula: |
---|
617 | barometric_formula_0d = p_0 * ( (t_0 - g_d_cp * z) / t_0 )**( cp_d_rd ) |
---|
618 | |
---|
619 | END FUNCTION barometric_formula_0d |
---|
620 | |
---|
621 | !--------------------------------------------------------------------------------------------------! |
---|
622 | ! Description: |
---|
623 | ! ------------ |
---|
624 | !> Compute the barometric formula for 1-D array arguments. The calculation is based on the |
---|
625 | !> assumption of a polytropic atmosphere and neutral stratification, where the temperature lapse |
---|
626 | !> rate is g/cp. |
---|
627 | !--------------------------------------------------------------------------------------------------! |
---|
628 | FUNCTION barometric_formula_1d( z, t_0, p_0) |
---|
629 | |
---|
630 | IMPLICIT NONE |
---|
631 | |
---|
632 | REAL(wp), INTENT(IN), DIMENSION(:) :: z !< height (m) |
---|
633 | REAL(wp), INTENT(IN) :: t_0 !< temperature reference state (K) |
---|
634 | REAL(wp), INTENT(IN) :: p_0 !< surface pressure (Pa) |
---|
635 | |
---|
636 | REAL(wp), DIMENSION(size(z)) :: barometric_formula_1d |
---|
637 | |
---|
638 | ! |
---|
639 | !-- Compute barometric formula: |
---|
640 | barometric_formula_1d = p_0 * ( (t_0 - g_d_cp * z) / t_0 )**( cp_d_rd ) |
---|
641 | |
---|
642 | END FUNCTION barometric_formula_1d |
---|
643 | |
---|
644 | END MODULE basic_constants_and_equations_mod |
---|