v0.14.0
Loading...
Searching...
No Matches
GenericClimateModel.hpp
Go to the documentation of this file.
1/*
2 * This file is part of MoFEM.
3 * MoFEM is free software: you can redistribute it and/or modify it under
4 * the terms of the GNU Lesser General Public License as published by the
5 * Free Software Foundation, either version 3 of the License, or (at your
6 * option) any later version.
7 *
8 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
9 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
10 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
11 * License for more details.
12 *
13 * You should have received a copy of the GNU Lesser General Public
14 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
15
16#ifndef __GenericClimateModel_HPP
17#define __GenericClimateModel_HPP
18
20
21 double T0; // reference temperature (K)
22 double e0; // reference saturation vapor pressure (es at a certain temp, usually 0 deg C) (Pa)
23 double Rv; // gas constant for water vapor (J*K/Kg)
24 double Lv; // latent heat of vaporization of water (J)
25
26 double u10; //< wind at high 10m (m/s)
27 double CR; //< cloudiness factor (0–1, dimensionless)
28 double Ta; //< air temperature (C)
29 double Td; //< dew point temperature (C)
30 double P; //< pressure
31 double Rs; //< observed solar radiation (W/m2)
32
33 double zenith; //topocentric zenith angle [degrees]
34 double azimuth; //topocentric azimuth angle (eastward from north) [for navigators and solar radiation]
35
36 /** \brief Clausius-Clapeyron equation
37 */
39 return e0*exp((Lv/Rv)*((1./T0)-(1./(T+T0))));
40 }
41
43 const double b = 17.2694;
44 const double T2 = 35.86;
45 return e0*exp(b*T/(T+T0-T2));
46 }
47
49 const double b = 17.2694;
50 const double T2 = 35.86;
51 return e0*exp(b*T/(T+T0-T2))*(b/(T+T0-T2) - b*T/pow(T+T0-T2,2));
52 }
53
54 double calculateVapourPressure(double T) {
55 //This use Tetent's formula by default
57 }
58
59 double calculateVapourPressure_dT(double T) {
60 //This use Tetent's formula by default
62 }
63
64 double calculateMixingRatio(double T,double P) {
65 const double eps = 0.622;
66 double e = calculateVapourPressure(T);
67 return eps*e/(P-e);
68 }
69
70 double calculateMixingRatio_dT(double T,double P) {
71 const double eps = 0.622;
72 double e = calculateVapourPressure(T);
73 double e_dT = calculateVapourPressure_dT(T);
74 return eps*e_dT/(P-e) + eps*e*e_dT/pow(P-e,2);
75 }
76
77 double calculateAbsoluteVirtualTempertaure(double T,double Td,double P) {
78 double r = calculateMixingRatio(Td,P);
79 return T*(0.61+r);
80 }
81
82 double calculateAbsoluteVirtualTempertaure_dT(double T,double Td,double P) {
83 double r = calculateMixingRatio(Td,P);
84 return (0.61+r);
85 }
86
88
89 T0 = 273.15; // reference temperature (K)
90 e0 = 611.3; // reference saturation vapor pressure (es at a certain temp, usually 0 deg C) (Pa)
91 Rv = 461.5; // gas constant for water vapor (J*K/Kg)
92 Lv = 2.5e6; // latent heat of vaporization of water (J)
93
94 u10 = 2.7; //< wind at high 10m (m/s)
95 CR = 0; //< cloudness factor (0–1, dimensionless)
96 Ta = 10; //< air temperature (C)
97 Td = 5; //< dew point temperature (C)
98 P = 101325; //< pressure
99 Rs = 1361; //< observed solar radiation (W/m2)
100
101 zenith = 0; //topocentric zenith angle [degrees]
102 azimuth = 0; //topocentric azimuth angle (eastward from north) [for navigators and solar radiation]
103
104 }
105
106 virtual MoFEMErrorCode set(double t = 0) = 0;
107};
108
109#endif //__GenericClimateModel_HPP
static const double eps
constexpr double t
plate stiffness
Definition plate.cpp:59
double calculateVapourPressure_dT(double T)
double calculateAbsoluteVirtualTempertaure_dT(double T, double Td, double P)
virtual MoFEMErrorCode set(double t=0)=0
double calculateVapourPressureClausiusClapeyron(double T)
Clausius-Clapeyron equation.
double calculateAbsoluteVirtualTempertaure(double T, double Td, double P)
double calculateVapourPressure(double T)
double calculateMixingRatio(double T, double P)
double calculateVapourPressureTetenFormula_dT(double T)
double calculateVapourPressureTetenFormula(double T)
double calculateMixingRatio_dT(double T, double P)