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.
fractional_step_settings_periodic.h
Go to the documentation of this file.
1 #ifndef KRATOS_FRACTIONAL_STEP_SETTINGS_PERIODIC_H
2 #define KRATOS_FRACTIONAL_STEP_SETTINGS_PERIODIC_H
3 
4 // System includes
5 #include <string>
6 #include <iostream>
7 
8 // External includes
9 
10 // Project includes
11 #include "includes/define.h"
20 #include "processes/process.h"
21 
22 // Application includes
24 #include "custom_utilities/solver_settings.h"
26 
27 
28 namespace Kratos
29 {
32 
35 
39 
43 
47 
51 
53 template< class TSparseSpace,
54  class TDenseSpace,
55  class TLinearSolver
56  >
57 class FractionalStepSettingsPeriodic: public SolverSettings<TSparseSpace,TDenseSpace,TLinearSolver>
58 {
59 public:
62 
65 
70 
73 
77 
80  const unsigned int ThisDomainSize,
81  const unsigned int ThisTimeOrder,
82  const bool UseSlip,
83  const bool MoveMeshFlag,
84  const bool ReformDofSet,
85  const Kratos::Variable<int>& rPeriodicVar):
86  BaseType(rModelPart,ThisDomainSize,ThisTimeOrder,UseSlip,MoveMeshFlag,ReformDofSet),
87  mrPeriodicVar(rPeriodicVar)
88  {}
89 
92 
96 
97 
101 
105 
106  void SetStrategy(StrategyLabel const& rStrategyLabel,
107  typename TLinearSolver::Pointer pLinearSolver,
108  const double Tolerance,
109  const unsigned int MaxIter) override
110  {
111  KRATOS_TRY;
112 
113  // pointer types for solution strategy construcion
114  typedef typename Scheme< TSparseSpace, TDenseSpace >::Pointer SchemePointerType;
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  bool UseSlip = BaseType::UseSlipConditions();
123  // Modification of the DofSet is managed by the fractional step strategy, not the auxiliary velocity and pressure strategies.
124  bool ReformDofSet = false; //BaseType::GetReformDofSet();
125  unsigned int EchoLevel = BaseType::GetEchoLevel();
126  unsigned int StrategyEchoLevel = (EchoLevel > 0) ? (EchoLevel-1) : 0;
127 
128  if ( rStrategyLabel == BaseType::Velocity )
129  {
130  // Velocity Builder and Solver
131  BuilderSolverTypePointer pBuildAndSolver = BuilderSolverTypePointer(new ResidualBasedBlockBuilderAndSolverPeriodic<TSparseSpace, TDenseSpace, TLinearSolver> (pLinearSolver, mrPeriodicVar));
132 
133  SchemePointerType pScheme;
134  //initializing fractional velocity solution step
135  if (UseSlip)
136  {
137  unsigned int DomainSize = BaseType::GetDomainSize();
138  SchemePointerType Temp = SchemePointerType(new ResidualBasedIncrementalUpdateStaticSchemeSlip< TSparseSpace, TDenseSpace > (DomainSize,DomainSize));
139  pScheme.swap(Temp);
140  }
141  else
142  {
143  SchemePointerType Temp = SchemePointerType(new ResidualBasedIncrementalUpdateStaticScheme< TSparseSpace, TDenseSpace > ());
144  pScheme.swap(Temp);
145  }
146 
147  // Strategy
149  rModelPart,
150  pScheme,
151  pBuildAndSolver,
152  CalculateReactions,
153  ReformDofSet,
154  CalculateNormDxFlag));
155 
156  }
157  else if ( rStrategyLabel == BaseType::Pressure )
158  {
159  // Pressure Builder and Solver
160  BuilderSolverTypePointer pBuildAndSolver = BuilderSolverTypePointer(new ResidualBasedBlockBuilderAndSolverPeriodic<TSparseSpace, TDenseSpace, TLinearSolver> (pLinearSolver, mrPeriodicVar));
161 
162  SchemePointerType pScheme = SchemePointerType(new ResidualBasedIncrementalUpdateStaticScheme< TSparseSpace, TDenseSpace > ());
163 
164  // Strategy
166  rModelPart,
167  pScheme,
168  pBuildAndSolver,
169  CalculateReactions,
170  ReformDofSet,
171  CalculateNormDxFlag));
172  }
173  else
174  {
175  KRATOS_THROW_ERROR(std::runtime_error,"Error in FractionalStepSettingsPeriodic: Unknown strategy label.","");
176  }
177 
178  BaseType::mTolerances[rStrategyLabel] = Tolerance;
179 
180  BaseType::mMaxIter[rStrategyLabel] = MaxIter;
181 
182  BaseType::mStrategies[rStrategyLabel]->SetEchoLevel(StrategyEchoLevel);
183 
184  KRATOS_CATCH("");
185  }
186 
187  void SetTurbulenceModel(TurbulenceModelLabel const& rTurbulenceModel,
188  typename TLinearSolver::Pointer pLinearSolver,
189  const double Tolerance,
190  const unsigned int MaxIter) override
191  {
192  KRATOS_TRY;
193 
195 
196  ModelPart& rModelPart = BaseType::GetModelPart();
197  unsigned int DomainSize = BaseType::GetDomainSize();
198  unsigned int TimeOrder = BaseType::GetTimeOrder();
199  bool ReformDofSet = BaseType::GetReformDofSet();
200 
201  if (rTurbulenceModel == BaseType::SpalartAllmaras)
202  {
204  (rModelPart,pLinearSolver,DomainSize,Tolerance,MaxIter,ReformDofSet,TimeOrder));
205  }
206  else
207  {
208  KRATOS_THROW_ERROR(std::runtime_error,"Error in FractionalStepSettingsPeriodic: Unknown turbulence model label.","");
209  }
210 
211  KRATOS_CATCH("");
212  }
213 
214  void SetTurbulenceModel(ProcessPointerType pTurbulenceModel) override
215  {
216  BaseType::SetTurbulenceModel(pTurbulenceModel);
217  }
218 
222 
226 
228  std::string Info() const override
229  {
230  std::stringstream buffer;
231  buffer << "FractionalStepSettingsPeriodic" ;
232  return buffer.str();
233  }
234 
236  void PrintInfo(std::ostream& rOStream) const override {rOStream << "FractionalStepSettingsPeriodic";}
237 
239  void PrintData(std::ostream& rOStream) const override {}
240 
241 
245 
246 
248 
249 protected:
252 
253 
257 
258 
262 
263 
267 
268 
272 
273 
277 
278 
282 
283 
285 
286 private:
289 
290 
294  const Kratos::Variable<int>& mrPeriodicVar;
295 
299 
300 
304 
305 
309 
310 
314 
315 
319 
322 
324  FractionalStepSettingsPeriodic& operator=(FractionalStepSettingsPeriodic const& rOther){}
325 
327  FractionalStepSettingsPeriodic(FractionalStepSettingsPeriodic const& rOther){}
328 
329 
331 
332 }; // Class FractionalStepSettingsPeriodic
333 
335 
338 
339 
343 
344 
346 template< class TDenseSpace, class TSparseSpace, class TLinearSolver >
347 inline std::istream& operator >> (std::istream& rIStream,
349 {
350  return rIStream;
351 }
352 
354 template< class TDenseSpace, class TSparseSpace, class TLinearSolver >
355 inline std::ostream& operator << (std::ostream& rOStream,
357 {
358  rThis.PrintInfo(rOStream);
359  rOStream << std::endl;
360  rThis.PrintData(rOStream);
361 
362  return rOStream;
363 }
365 
367 
368 } // namespace Kratos.
369 
370 #endif // KRATOS_FRACTIONAL_STEP_SETTINGS_PERIODIC_H
Current class provides an implementation for the base builder and solving operations.
Definition: builder_and_solver.h:64
Helper class to define solution strategies for FS_Strategy.
Definition: fractional_step_settings_periodic.h:58
BaseType::StrategyPointerType StrategyPointerType
Definition: fractional_step_settings_periodic.h:68
std::string Info() const override
Turn back information as a string.
Definition: fractional_step_settings_periodic.h:228
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: fractional_step_settings_periodic.h:239
FractionalStepSettingsPeriodic(ModelPart &rModelPart, const unsigned int ThisDomainSize, const unsigned int ThisTimeOrder, const bool UseSlip, const bool MoveMeshFlag, const bool ReformDofSet, const Kratos::Variable< int > &rPeriodicVar)
Constructor.
Definition: fractional_step_settings_periodic.h:79
BaseType::StrategyLabel StrategyLabel
Definition: fractional_step_settings_periodic.h:71
SolverSettings< TSparseSpace, TDenseSpace, TLinearSolver > BaseType
Definition: fractional_step_settings_periodic.h:66
~FractionalStepSettingsPeriodic() override
Destructor.
Definition: fractional_step_settings_periodic.h:91
BaseType::TurbulenceModelLabel TurbulenceModelLabel
Definition: fractional_step_settings_periodic.h:72
BaseType::StrategyType StrategyType
Definition: fractional_step_settings_periodic.h:67
void SetStrategy(StrategyLabel const &rStrategyLabel, typename TLinearSolver::Pointer pLinearSolver, const double Tolerance, const unsigned int MaxIter) override
Definition: fractional_step_settings_periodic.h:106
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: fractional_step_settings_periodic.h:236
BaseType::ProcessPointerType ProcessPointerType
Definition: fractional_step_settings_periodic.h:69
void SetTurbulenceModel(ProcessPointerType pTurbulenceModel) override
Definition: fractional_step_settings_periodic.h:214
void SetTurbulenceModel(TurbulenceModelLabel const &rTurbulenceModel, typename TLinearSolver::Pointer pLinearSolver, const double Tolerance, const unsigned int MaxIter) override
Definition: fractional_step_settings_periodic.h:187
KRATOS_CLASS_POINTER_DEFINITION(FractionalStepSettingsPeriodic)
Pointer definition of FractionalStepSettingsPeriodic.
Implicit solving strategy base class This is the base class from which we will derive all the implici...
Definition: implicit_solving_strategy.h:61
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
Definition: residualbased_block_builder_and_solver_periodic.h:69
This class provides the implementation of a static scheme.
Definition: residualbased_incrementalupdate_static_scheme.h:57
Scheme for the solution of problems involving a slip condition.
Definition: residualbased_incrementalupdate_static_scheme_slip.h:70
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
Helper class to define solution strategies for FS_Strategy.
Definition: solver_settings.h:53
StrategyLabel
Definition: solver_settings.h:66
@ Velocity
Definition: solver_settings.h:66
@ Pressure
Definition: solver_settings.h:66
virtual unsigned int GetDomainSize() const
Definition: solver_settings.h:146
ModelPart & GetModelPart()
Definition: solver_settings.h:285
virtual void SetTurbulenceModel(TurbulenceModelLabel const &rTurbulenceModel, typename TLinearSolver::Pointer pLinearSolver, const double Tolerance, const unsigned int MaxIter)
Definition: solver_settings.h:124
bool GetReformDofSet()
Definition: solver_settings.h:226
StrategyType::Pointer StrategyPointerType
Definition: solver_settings.h:62
virtual bool UseSlipConditions() const
Definition: solver_settings.h:156
std::map< StrategyLabel, double > mTolerances
Definition: solver_settings.h:292
std::map< StrategyLabel, StrategyPointerType > mStrategies
Definition: solver_settings.h:290
std::map< StrategyLabel, unsigned int > mMaxIter
Definition: solver_settings.h:294
Process::Pointer ProcessPointerType
Definition: solver_settings.h:63
virtual unsigned int GetEchoLevel()
Definition: solver_settings.h:221
TurbulenceModelLabel
Definition: solver_settings.h:68
@ SpalartAllmaras
Definition: solver_settings.h:68
bool mHaveTurbulenceModel
Definition: solver_settings.h:298
ProcessPointerType mpTurbulenceModel
Definition: solver_settings.h:296
virtual unsigned int GetTimeOrder() const
Definition: solver_settings.h:151
An impelementation of the Spalart-Allmaras turbulence model for incompressible flows.
Definition: spalart_allmaras_turbulence_model.h:77
#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