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.
temporal_min_method.h
Go to the documentation of this file.
1 // | / |
2 // ' / __| _` | __| _ \ __|
3 // . \ | ( | | ( |\__ `
4 // _|\_\_| \__,_|\__|\___/ ____/
5 // Multi-Physics
6 //
7 // License: BSD License
8 // Kratos default license: kratos/license.txt
9 //
10 // Main authors: Suneth Warnakulasuriya (https://github.com/sunethwarna)
11 //
12 
13 #if !defined(KRATOS_TEMPORAL_MIN_METHOD_H_INCLUDED)
14 #define KRATOS_TEMPORAL_MIN_METHOD_H_INCLUDED
15 
16 // System includes
17 #include <limits>
18 #include <string>
19 
20 // External includes
21 
22 // Project includes
23 #include "includes/define.h"
24 #include "includes/model_part.h"
25 
26 // Application includes
30 
31 namespace Kratos
32 {
35 
38 
39 namespace TemporalMethods
40 {
41 template <class TContainerType, class TContainerItemType, template <class T> class TDataRetrievalFunctor, template <class T> class TDataStorageFunctor>
43 {
44 public:
45  template <class TDataType>
46  class NormMethod : public TemporalMethod
47  {
48  public:
50 
52  ModelPart& rModelPart,
53  const std::string& rNormType,
54  const Variable<TDataType>& rInputVariable,
55  const int EchoLevel,
56  const Variable<double>& rOutputVariable,
57  const Variable<double>& rMinTimeValueVariable)
58  : TemporalMethod(rModelPart, EchoLevel),
59  mNormType(rNormType),
60  mrInputVariable(rInputVariable),
61  mrOutputVariable(rOutputVariable),
62  mrMinTimeValueVariable(rMinTimeValueVariable)
63  {
65 
66  KRATOS_ERROR_IF(rOutputVariable == rMinTimeValueVariable)
67  << "Same variable is given for minimum and its time value "
68  "variable in norm min method for input variable "
69  << rInputVariable.Name()
70  << ". Please provide two different "
71  "variables. [ variable = "
72  << rOutputVariable.Name() << " ].\n";
73 
74  KRATOS_CATCH("");
75  }
76 
77  void CalculateStatistics() override
78  {
79  TContainerType& r_container =
80  MethodUtilities::GetDataContainer<TContainerType>(this->GetModelPart());
81 
82  const auto& norm_method =
83  MethodUtilities::GetNormMethod(mrInputVariable, mNormType);
84 
85  const double total_time = this->GetTotalTime();
86 
87  const int number_of_items = r_container.size();
88 #pragma omp parallel for
89  for (int i = 0; i < number_of_items; ++i)
90  {
91  TContainerItemType& r_item = *(r_container.begin() + i);
92  const TDataType& r_input_value =
93  TDataRetrievalFunctor<TContainerItemType>()(r_item, mrInputVariable);
94  const double input_norm_value = norm_method(r_input_value);
95 
96  double& r_output_value =
97  TDataStorageFunctor<TContainerItemType>()(r_item, mrOutputVariable);
98  double& r_min_time_value = TDataStorageFunctor<TContainerItemType>()(
99  r_item, mrMinTimeValueVariable);
100 
101  if (input_norm_value < r_output_value)
102  {
103  r_output_value = input_norm_value;
104  r_min_time_value = total_time;
105  }
106  }
107 
108  KRATOS_INFO_IF("TemporalNormMinMethod", this->GetEchoLevel() > 1)
109  << "Calculated temporal norm min for " << mrInputVariable.Name()
110  << " input variable with " << mrOutputVariable.Name()
111  << " min variable and " << mrMinTimeValueVariable.Name()
112  << " time value variable for " << this->GetModelPart().Name() << ".\n";
113  }
114 
116  {
117  TContainerType& r_container =
118  MethodUtilities::GetDataContainer<TContainerType>(this->GetModelPart());
119 
120  auto& initializer_method =
121  TemporalMethodUtilities::InitializeVariables<TContainerType, TContainerItemType, TDataStorageFunctor>;
122  initializer_method(
123  r_container, mrOutputVariable, std::numeric_limits<double>::max());
124  initializer_method(r_container, mrMinTimeValueVariable, 0.0);
125 
126  KRATOS_INFO_IF("TemporalNormMinMethod", this->GetEchoLevel() > 0)
127  << "Initialized temporal norm min method for "
128  << mrInputVariable.Name() << " input variable with "
129  << mrOutputVariable.Name() << " min variable and "
130  << mrMinTimeValueVariable.Name() << " time value variable for "
131  << this->GetModelPart().Name() << ".\n";
132  }
133 
134  private:
135  const std::string mNormType;
136  const Variable<TDataType>& mrInputVariable;
137  const Variable<double>& mrOutputVariable;
138  const Variable<double>& mrMinTimeValueVariable;
139  };
140 
141  std::vector<TemporalMethod::Pointer> static CreateTemporalMethodObject(
142  ModelPart& rModelPart, const std::string& rNormType, const int EchoLevel, Parameters Params)
143  {
144  KRATOS_TRY
145 
146  Parameters default_parameters = Parameters(R"(
147  {
148  "input_variables" : [],
149  "output_variables" : [],
150  "output_time_step_variables" : []
151  })");
152  Params.RecursivelyValidateAndAssignDefaults(default_parameters);
153 
154  const std::vector<std::string>& input_variable_names_list =
155  Params["input_variables"].GetStringArray();
156  const std::vector<std::string>& output_variable_1_names_list =
157  Params["output_variables"].GetStringArray();
158  const std::vector<std::string>& output_variable_2_names_list =
159  Params["output_time_step_variables"].GetStringArray();
160 
161  std::vector<TemporalMethod::Pointer> method_list;
162  if (rNormType == "none") // for non norm types
163  {
164  KRATOS_ERROR << "none norm type is not defined for Min method.\n";
165  }
166  else // for values with norms
167  {
168  MethodUtilities::CheckVariableType<double>(output_variable_1_names_list);
169  MethodUtilities::CheckVariableType<double>(output_variable_2_names_list);
170 
171  const int number_of_variables = input_variable_names_list.size();
172  for (int i = 0; i < number_of_variables; ++i)
173  {
174  const std::string& r_variable_input_name = input_variable_names_list[i];
175  const std::string& r_variable_1_output_name =
176  output_variable_1_names_list[i];
177  const std::string& r_variable_2_output_name =
178  output_variable_2_names_list[i];
180  rModelPart, rNormType, r_variable_input_name, EchoLevel,
181  r_variable_1_output_name, r_variable_2_output_name, method_list, NormMethod)
182  }
183  }
184 
185  return method_list;
186 
187  KRATOS_CATCH("");
188  }
189 };
190 } // namespace TemporalMethods
191 } // namespace Kratos
192 
193 #endif // KRATOS_TEMPORAL_MIN_METHOD_H_INCLUDED
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
std::string & Name()
Definition: model_part.h:1811
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
std::vector< std::string > GetStringArray() const
This method returns the array of strings in the current Parameter.
Definition: kratos_parameters.cpp:693
void RecursivelyValidateAndAssignDefaults(const Parameters &rDefaultParameters)
This function is designed to verify that the parameters under testing match the form prescribed by th...
Definition: kratos_parameters.cpp:1389
Definition: temporal_method.h:40
Definition: temporal_min_method.h:47
Definition: temporal_min_method.h:43
const std::string & Name() const
Definition: variable_data.h:201
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_TRY
Definition: define.h:109
#define KRATOS_ERROR
Definition: exception.h:161
#define KRATOS_ERROR_IF(conditional)
Definition: exception.h:162
#define KRATOS_INFO_IF(label, conditional)
Definition: logger.h:251
int GetEchoLevel() const
Definition: temporal_method.h:112
ModelPart & GetModelPart() const
Definition: temporal_method.h:82
double GetTotalTime() const
Definition: temporal_method.h:87
NormMethod(ModelPart &rModelPart, const std::string &rNormType, const Variable< TDataType > &rInputVariable, const int EchoLevel, const Variable< double > &rOutputVariable, const Variable< double > &rMinTimeValueVariable)
Definition: temporal_min_method.h:51
void InitializeStatisticsVariables() override
Definition: temporal_min_method.h:115
static std::vector< TemporalMethod::Pointer > CreateTemporalMethodObject(ModelPart &rModelPart, const std::string &rNormType, const int EchoLevel, Parameters Params)
Definition: temporal_min_method.h:141
void CalculateStatistics() override
Definition: temporal_min_method.h:77
#define ADD_TEMPORAL_NORM_METHOD_TWO_OUTPUT_VARIABLE_OBJECT( model_part, norm_type, input_variable, echo_level, output_variable_1, output_variable_2, object_list, method)
Definition: method_utilities.h:190
static int EchoLevel
Definition: co_sim_EMPIRE_API.h:42
static double max(double a, double b)
Definition: GeometryFunctions.h:79
const std::function< double(const TDataType &)> GetNormMethod(const Variable< TDataType > &rVariable, const std::string &rNormType)
Definition: method_utilities.cpp:313
template void CheckVariableType< double >(const std::vector< std::string > &rVariableNamesList)
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
integer i
Definition: TensorModule.f:17