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.
two_step_v_p_settings.h
Go to the documentation of this file.
1 //
2 // Project Name: KratosPFEMFluidDynamicsApplication $
3 // Last modified by: $Author: AFranci $
4 // Date: $Date: January 2016 $
5 // Revision: $Revision: 0.0 $
6 //
7 //
8 
9 #ifndef KRATOS_TWO_STEP_V_P_SETTINGS_H
10 #define KRATOS_TWO_STEP_V_P_SETTINGS_H
11 
12 // System includes
13 #include <string>
14 #include <iostream>
15 
16 // External includes
17 
18 // Project includes
19 #include "includes/define.h"
29 #include "processes/process.h"
30 
31 // Application includes
32 #include "custom_utilities/solver_settings.h"
33 
34 namespace Kratos
35 {
38 
41 
45 
49 
53 
57 
59  template <class TSparseSpace,
60  class TDenseSpace,
61  class TLinearSolver>
62  class TwoStepVPSettings : public TwoStepVPSolverSettings<TSparseSpace, TDenseSpace, TLinearSolver>
63  {
64  public:
67 
70 
75 
77 
81 
84  const unsigned int ThisDomainSize,
85  const unsigned int ThisTimeOrder,
86  const bool ReformDofSet) : BaseType(rModelPart, ThisDomainSize, ThisTimeOrder, ReformDofSet)
87  {
88  }
89 
91  virtual ~TwoStepVPSettings() {}
92 
96 
100 
104 
105  virtual void SetStrategy(StrategyLabel const &rStrategyLabel,
106  typename TLinearSolver::Pointer pLinearSolver,
107  const double Tolerance,
108  const unsigned int MaxIter)
109  {
110  KRATOS_TRY;
111 
112  // pointer types for solution strategy construcion
113  typedef typename Scheme<TSparseSpace, TDenseSpace>::Pointer SchemePointerType;
114  //typedef typename ConvergenceCriteria< TSparseSpace, TDenseSpace >::Pointer ConvergenceCriteriaPointerType;
115  typedef typename BuilderAndSolver<TSparseSpace, TDenseSpace, TLinearSolver>::Pointer BuilderSolverTypePointer;
116 
117  // Default, fixed flags
118  bool CalculateReactions = false;
119  bool CalculateNormDxFlag = true;
120 
121  ModelPart &rModelPart = BaseType::GetModelPart();
122  // Modification of the DofSet is managed by the fractional step strategy, not the auxiliary velocity and pressure strategies.
123  bool ReformDofSet = false; //BaseType::GetReformDofSet();
124  unsigned int EchoLevel = BaseType::GetEchoLevel();
125  unsigned int StrategyEchoLevel = (EchoLevel > 0) ? (EchoLevel - 1) : 0;
126 
127  if (rStrategyLabel == BaseType::Velocity)
128  {
129  // Velocity Builder and Solver
130  BuilderSolverTypePointer pBuildAndSolver = BuilderSolverTypePointer(new ResidualBasedBlockBuilderAndSolver<TSparseSpace, TDenseSpace, TLinearSolver>(pLinearSolver));
131 
132  SchemePointerType pScheme;
133  //initializing fractional velocity solution step
134 
135  SchemePointerType Temp = SchemePointerType(new ResidualBasedIncrementalUpdateStaticScheme<TSparseSpace, TDenseSpace>());
136  pScheme.swap(Temp);
137 
138  // Strategy
139  BaseType::mStrategies[rStrategyLabel] = StrategyPointerType(new ResidualBasedLinearStrategy<TSparseSpace, TDenseSpace, TLinearSolver>(rModelPart, pScheme, pLinearSolver, pBuildAndSolver, CalculateReactions, ReformDofSet, CalculateNormDxFlag));
140  }
141  else if (rStrategyLabel == BaseType::Pressure)
142  {
143  // Pressure Builder and Solver
144  // BuilderSolverTypePointer pBuildAndSolver = BuilderSolverTypePointer(new ResidualBasedEliminationBuilderAndSolverComponentwise<TSparseSpace, TDenseSpace, TLinearSolver, Variable<double> >(pLinearSolver, PRESSURE));
145  BuilderSolverTypePointer pBuildAndSolver = BuilderSolverTypePointer(new ResidualBasedBlockBuilderAndSolver<TSparseSpace, TDenseSpace, TLinearSolver>(pLinearSolver));
146  SchemePointerType pScheme = SchemePointerType(new ResidualBasedIncrementalUpdateStaticScheme<TSparseSpace, TDenseSpace>());
147 
148  // Strategy
149  BaseType::mStrategies[rStrategyLabel] = StrategyPointerType(new ResidualBasedLinearStrategy<TSparseSpace, TDenseSpace, TLinearSolver>(rModelPart, pScheme, pLinearSolver, pBuildAndSolver, CalculateReactions, ReformDofSet, CalculateNormDxFlag));
150  }
151 
152  else
153  {
154  KRATOS_THROW_ERROR(std::runtime_error, "Error in TwoStepVPSettings: Unknown strategy label.", "");
155  }
156 
157  BaseType::mTolerances[rStrategyLabel] = Tolerance;
158 
159  BaseType::mMaxIter[rStrategyLabel] = MaxIter;
160 
161  BaseType::mStrategies[rStrategyLabel]->SetEchoLevel(StrategyEchoLevel);
162 
163  KRATOS_CATCH("");
164  }
165 
169 
173 
175  virtual std::string Info() const
176  {
177  std::stringstream buffer;
178  buffer << "TwoStepVPSettings";
179  return buffer.str();
180  }
181 
183  virtual void PrintInfo(std::ostream &rOStream) const { rOStream << "TwoStepVPSettings"; }
184 
186  virtual void PrintData(std::ostream &rOStream) const {}
187 
191 
193 
194  protected:
197 
201 
205 
209 
213 
217 
221 
223 
224  private:
227 
231 
235 
239 
243 
247 
251 
253  TwoStepVPSettings() {}
254 
256  TwoStepVPSettings &operator=(TwoStepVPSettings const &rOther) {}
257 
259  TwoStepVPSettings(TwoStepVPSettings const &rOther) {}
260 
262 
263  }; // Class TwoStepVPSettings
264 
266 
269 
273 
275  template <class TDenseSpace, class TSparseSpace, class TLinearSolver>
276  inline std::istream &operator>>(std::istream &rIStream,
278  {
279  return rIStream;
280  }
281 
283  template <class TDenseSpace, class TSparseSpace, class TLinearSolver>
284  inline std::ostream &operator<<(std::ostream &rOStream,
286  {
287  rThis.PrintInfo(rOStream);
288  rOStream << std::endl;
289  rThis.PrintData(rOStream);
290 
291  return rOStream;
292  }
294 
296 
297 } // namespace Kratos.
298 
299 #endif // KRATOS_V_P_SETTINGS_H
Current class provides an implementation for the base builder and solving operations.
Definition: builder_and_solver.h:64
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
Current class provides an implementation for standard builder and solving operations.
Definition: residualbased_block_builder_and_solver.h:82
This class provides the implementation of a static scheme.
Definition: residualbased_incrementalupdate_static_scheme.h:57
This is a very simple strategy to solve linearly the problem.
Definition: residualbased_linear_strategy.h:64
This class provides the implementation of the basic tasks that are needed by the solution strategy.
Definition: scheme.h:56
Solving strategy base class This is the base class from which we will derive all the strategies (impl...
Definition: solving_strategy.h:64
Helper class to define solution strategies for FS_Strategy.
Definition: two_step_v_p_settings.h:63
virtual void PrintInfo(std::ostream &rOStream) const
Print information about this object.
Definition: two_step_v_p_settings.h:183
BaseType::StrategyLabel StrategyLabel
Definition: two_step_v_p_settings.h:76
BaseType::StrategyPointerType StrategyPointerType
Definition: two_step_v_p_settings.h:73
BaseType::StrategyType StrategyType
Definition: two_step_v_p_settings.h:72
TwoStepVPSolverSettings< TSparseSpace, TDenseSpace, TLinearSolver > BaseType
Definition: two_step_v_p_settings.h:71
BaseType::ProcessPointerType ProcessPointerType
Definition: two_step_v_p_settings.h:74
KRATOS_CLASS_POINTER_DEFINITION(TwoStepVPSettings)
Pointer definition of TwoStepVPSettings.
virtual std::string Info() const
Turn back information as a string.
Definition: two_step_v_p_settings.h:175
virtual void PrintData(std::ostream &rOStream) const
Print object's data.
Definition: two_step_v_p_settings.h:186
TwoStepVPSettings(ModelPart &rModelPart, const unsigned int ThisDomainSize, const unsigned int ThisTimeOrder, const bool ReformDofSet)
Constructor.
Definition: two_step_v_p_settings.h:83
virtual ~TwoStepVPSettings()
Destructor.
Definition: two_step_v_p_settings.h:91
virtual void SetStrategy(StrategyLabel const &rStrategyLabel, typename TLinearSolver::Pointer pLinearSolver, const double Tolerance, const unsigned int MaxIter)
Definition: two_step_v_p_settings.h:105
Helper class to define solution strategies for TwoStepVPStrategy.
Definition: solver_settings.h:57
ModelPart & GetModelPart()
Definition: solver_settings.h:234
StrategyLabel
Definition: solver_settings.h:70
@ Velocity
Definition: solver_settings.h:71
@ Pressure
Definition: solver_settings.h:72
std::map< StrategyLabel, double > mTolerances
Definition: solver_settings.h:241
virtual unsigned int GetEchoLevel()
Definition: solver_settings.h:176
std::map< StrategyLabel, unsigned int > mMaxIter
Definition: solver_settings.h:243
Process::Pointer ProcessPointerType
Definition: solver_settings.h:67
StrategyType::Pointer StrategyPointerType
Definition: solver_settings.h:66
std::map< StrategyLabel, StrategyPointerType > mStrategies
Definition: solver_settings.h:239
#define KRATOS_THROW_ERROR(ExceptionType, ErrorMessage, MoreInfo)
Definition: define.h:77
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_TRY
Definition: define.h:109
static int EchoLevel
Definition: co_sim_EMPIRE_API.h:42
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
Temp
Definition: PecletTest.py:105