KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
T_microclimate_flux_condition.h
Go to the documentation of this file.
1 // KRATOS___
2 // // ) )
3 // // ___ ___
4 // // ____ //___) ) // ) )
5 // // / / // // / /
6 // ((____/ / ((____ ((___/ / MECHANICS
7 //
8 // License: geo_mechanics_application/license.txt
9 //
10 //
11 // Main authors: Mohamed Nabi
12 //
13 
14 #pragma once
15 
16 // Project includes
17 #include "includes/serializer.h"
18 
19 // Application includes
22 #include "custom_utilities/element_utilities.hpp"
24 
25 namespace Kratos
26 {
27 
29 {
30  double precipitation = 0.0;
31  double evaporation = 0.0;
32 };
33 
34 template <unsigned int TDim, unsigned int TNumNodes>
35 class KRATOS_API(GEO_MECHANICS_APPLICATION) GeoTMicroClimateFluxCondition
36  : public GeoTCondition<TDim, TNumNodes>
37 {
38 public:
40 
41  using IndexType = std::size_t;
44 
45  GeoTMicroClimateFluxCondition() : GeoTCondition<TDim, TNumNodes>() {}
46 
47  GeoTMicroClimateFluxCondition(IndexType NewId, GeometryType::Pointer pGeometry)
48  : GeoTCondition<TDim, TNumNodes>(NewId, pGeometry)
49  {
50  }
51 
53  GeometryType::Pointer pGeometry,
54  Properties::Pointer pProperties)
55  : GeoTCondition<TDim, TNumNodes>(NewId, pGeometry, pProperties)
56  {
57  }
58 
59  Condition::Pointer Create(IndexType NewId,
60  const NodesArrayType& rNodes,
61  Properties::Pointer pProperties) const override;
62 
63  void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
64 
65  void InitializeSolutionStep(const ProcessInfo& rCurrentProcessInfo) override;
66 
67  void CalculateLocalSystem(Matrix& rLeftHandSideMatrix,
68  Vector& rRightHandSideVector,
69  const ProcessInfo& rCurrentProcessInfo) override;
70 
71 private:
72  void InitializeProperties();
73 
74  void CalculateAndAddLHS(Matrix& rLeftHandSideMatrix,
76  double IntegrationCoefficient,
77  const array_1d<double, TNumNodes>& rLeftHandSideFluxes);
78 
79  void CalculateAndAddRHS(Vector& rRightHandSideVector,
81  double IntegrationCoefficient,
82  const Vector& rNodalTemperatures,
83  const array_1d<double, TNumNodes>& rLeftHandSideFluxes,
84  const array_1d<double, TNumNodes>& rRightHandSideFluxes);
85 
86  array_1d<double, TNumNodes> CalculateLeftHandSideFluxes() const;
87  array_1d<double, TNumNodes> CalculateRightHandSideFluxes(double TimeStepSize,
88  double PreviousStorage,
89  double PreviousRadiation) const;
90  double CalculateRightHandSideFlux(double NetRadiation,
91  double SurfaceHeatStorage,
92  double ActualEvaporation) const;
93 
94  double CalculateCurrentWaterStorage(double TimeStepSize,
95  double PreviousStorage,
96  double PreviousRadiation) const;
97 
98  double CalculateCurrentNetRadiation() const;
99  double CalculateNetRadiation(unsigned int NodeIndex) const;
100 
101  double CalculateSurfaceHeatStorage(double TimeStepSize,
102  double PreviousRadiation,
103  double NetRadiation) const;
104 
105  WaterFluxes CalculateWaterFluxes(unsigned int NodeIndex,
106  double TimeStepSize,
107  double PreviousStorage,
108  double NetRadiation,
109  double SurfaceHeatStorage) const;
110 
111  double CalculatePotentialEvaporation(unsigned int NodeIndex,
112  double NetRadiation,
113  double SurfaceHeatStorage) const;
114 
115  void CalculateRoughness(const ProcessInfo& rCurrentProcessInfo);
116 
117  double CalculateSurfaceRoughnessFactor(double CurrentAirTemperature,
118  double PreviousRoughnessTemperature,
119  double RichardsonBulkModulus,
120  double FrictionDragCoefficient) const;
121 
122  // Serialization
123  friend class Serializer;
124 
125  void save(Serializer& rSerializer) const override
126  {
128  rSerializer.save("mIsInitialized", mIsInitialized);
129  rSerializer.save("mAlbedoCoefficient", mAlbedoCoefficient);
130  rSerializer.save("mFirstCoverStorageCoefficient", mFirstCoverStorageCoefficient);
131  rSerializer.save("mSecondCoverStorageCoefficient", mSecondCoverStorageCoefficient);
132  rSerializer.save("mThirdCoverStorageCoefficient", mThirdCoverStorageCoefficient);
133  rSerializer.save("mBuildEnvironmentRadiation", mBuildEnvironmentRadiation);
134  rSerializer.save("mMinimalStorage", mMinimalStorage);
135  rSerializer.save("mMaximalStorage", mMaximalStorage);
136  rSerializer.save("mRoughnessTemperature", mRoughnessTemperature);
137  rSerializer.save("mNetRadiation", mNetRadiation);
138  rSerializer.save("mWaterStorage", mWaterStorage);
139  rSerializer.save("mWaterDensity", mWaterDensity);
140  }
141 
142  void load(Serializer& rSerializer) override
143  {
145  rSerializer.load("mIsInitialized", mIsInitialized);
146  rSerializer.load("mAlbedoCoefficient", mAlbedoCoefficient);
147  rSerializer.load("mFirstCoverStorageCoefficient", mFirstCoverStorageCoefficient);
148  rSerializer.load("mSecondCoverStorageCoefficient", mSecondCoverStorageCoefficient);
149  rSerializer.load("mThirdCoverStorageCoefficient", mThirdCoverStorageCoefficient);
150  rSerializer.load("mBuildEnvironmentRadiation", mBuildEnvironmentRadiation);
151  rSerializer.load("mMinimalStorage", mMinimalStorage);
152  rSerializer.load("mMaximalStorage", mMaximalStorage);
153  rSerializer.load("mRoughnessTemperature", mRoughnessTemperature);
154  rSerializer.load("mNetRadiation", mNetRadiation);
155  rSerializer.load("mWaterStorage", mWaterStorage);
156  rSerializer.load("mWaterDensity", mWaterDensity);
157  }
158 
159  bool mIsInitialized = false;
160  double mAlbedoCoefficient = 0.0;
161  double mFirstCoverStorageCoefficient = 0.0;
162  double mSecondCoverStorageCoefficient = 0.0;
163  double mThirdCoverStorageCoefficient = 0.0;
164  double mBuildEnvironmentRadiation = 0.0;
165  double mMinimalStorage = 0.0;
166  double mMaximalStorage = 0.0;
167  double mRoughnessTemperature = 0.0;
168  double mNetRadiation = 0.0;
169  double mWaterStorage = 0.0;
170  double mWaterDensity = 0.0;
171 };
172 
173 // class TMicroClimateFluxCondition.
174 
175 } // namespace Kratos.
Base class for all Conditions.
Definition: condition.h:59
std::size_t IndexType
Definition: flags.h:74
Definition: T_condition.h:27
Definition: T_microclimate_flux_condition.h:37
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(GeoTMicroClimateFluxCondition)
GeoTMicroClimateFluxCondition(IndexType NewId, GeometryType::Pointer pGeometry)
Definition: T_microclimate_flux_condition.h:47
GeoTMicroClimateFluxCondition()
Definition: T_microclimate_flux_condition.h:45
GeoTMicroClimateFluxCondition(IndexType NewId, GeometryType::Pointer pGeometry, Properties::Pointer pProperties)
Definition: T_microclimate_flux_condition.h:52
Geometry base class.
Definition: geometry.h:71
PointerVector< TPointType > PointsArrayType
Definition: geometry.h:118
PointerVector is a container like stl vector but using a vector to store pointers to its data.
Definition: pointer_vector.h:72
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
void load(std::string const &rTag, TDataType &rObject)
Definition: serializer.h:207
void save(std::string const &rTag, std::array< TDataType, TDataSize > const &rObject)
Definition: serializer.h:545
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
std::size_t IndexType
The definition of the index type.
Definition: key_hash.h:35
Modeler::Pointer Create(const std::string &ModelerName, Model &rModel, const Parameters ModelParameters)
Checks if the modeler is registered.
Definition: modeler_factory.cpp:30
void InitializeSolutionStep(ConstructionUtility &rThisUtil, std::string ThermalSubModelPartName, std::string MechanicalSubModelPartName, std::string HeatFluxSubModelPartName, std::string HydraulicPressureSubModelPartName, bool thermal_conditions, bool mechanical_conditions, int phase)
Definition: add_custom_utilities_to_python.cpp:45
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
ModelPart::NodesContainerType NodesArrayType
Definition: gid_gauss_point_container.h:42
def load(f)
Definition: ode_solve.py:307
Definition: T_microclimate_flux_condition.h:29
double precipitation
Definition: T_microclimate_flux_condition.h:30
double evaporation
Definition: T_microclimate_flux_condition.h:31