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.
ulf_time_step_dec_process.h
Go to the documentation of this file.
1 /*
2 ==============================================================================
3 KratosULFApplication
4 A library based on:
5 Kratos
6 A General Purpose Software for Multi-Physics Finite Element Analysis
7 Version 1.0 (Released on march 05, 2007).
8 
9 Copyright 2007
10 Pooyan Dadvand, Riccardo Rossi, Pawel Ryzhakov
11 pooyan@cimne.upc.edu
12 rrossi@cimne.upc.edu
13 - CIMNE (International Center for Numerical Methods in Engineering),
14 Gran Capita' s/n, 08034 Barcelona, Spain
15 
16 
17 Permission is hereby granted, free of charge, to any person obtaining
18 a copy of this software and associated documentation files (the
19 "Software"), to deal in the Software without restriction, including
20 without limitation the rights to use, copy, modify, merge, publish,
21 distribute, sublicense and/or sell copies of the Software, and to
22 permit persons to whom the Software is furnished to do so, subject to
23 the following condition:
24 
25 Distribution of this code for any commercial purpose is permissible
26 ONLY BY DIRECT ARRANGEMENT WITH THE COPYRIGHT OWNERS.
27 
28 The above copyright notice and this permission notice shall be
29 included in all copies or substantial portions of the Software.
30 
31 THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
32 EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
33 MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
34 IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
35 CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
36 TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
37 SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
38 
39 ==============================================================================
40 */
41 
42 
43 //
44 // Project Name: Kratos
45 // Last Modified by: $Author: anonymous $
46 // Date: $Date: 2008-10-02 10:47:22 $
47 // Revision: $Revision: 1.7 $
48 //
49 // NOW FOR 2D ONLY!!! WRITE a 3D version of it!!!!
50 
51 #if !defined(KRATOS_ULF_TIME_STEP_DEC_PROCESS_INCLUDED )
52 #define KRATOS_ULF_TIME_STEP_DEC_PROCESS_INCLUDED
53 
54 
55 
56 // System includes
57 #include <string>
58 #include <iostream>
59 #include <algorithm>
60 
61 // External includes
62 
63 // Project includes
64 #include "includes/define.h"
65 #include "processes/process.h"
66 #include "includes/node.h"
67 #include "includes/element.h"
68 #include "includes/model_part.h"
69 #include "utilities/geometry_utilities.h"
70 
71 //to compute the determinant
72 #include "utilities/math_utils.h"
73 
74 namespace Kratos
75 {
76 
79 
83 
84 
88 
92 
96 
98 
105  : public Process
106 {
107 public:
110 
113 
117 
120  : mr_model_part(model_part)
121  {
122  }
123 
126  {
127  }
128 
129 
133 
134  void operator()()
135  {
136 
137  }
138 
139 
143 
144  double EstimateDeltaTime(const double dt_max, const double domain_size)
145  {
146  double estimated_dt = 0.00;
147  if(domain_size == 2)
148  estimated_dt = EstimateDeltaTimeTemplated<2>(dt_max);
149  else
150  estimated_dt = EstimateDeltaTimeTemplated<3>(dt_max);
151  return estimated_dt;
152 
153  }
154 
155  template< int TDim>
156  double EstimateDeltaTimeTemplated(const double dt_max)
157  {
158  KRATOS_TRY
159 
160  double deltatime_new;
161  double deltatime = dt_max;
164  //this is for 2D only
172  double Area;
173  double detJ;
174 
175 
176  for(ModelPart::ElementsContainerType::iterator i = mr_model_part.ElementsBegin();
177  i!=mr_model_part.ElementsEnd(); i++)
178  {
179  //calculating the actual Jacobian
180  for(unsigned int iii = 0; iii < i->GetGeometry().size(); iii++)
181  {
182  const array_1d<double,3>& v = i->GetGeometry()[iii].FastGetSolutionStepValue(VELOCITY);
183  for(unsigned int j=0; j <TDim; j++)
184  aux(j,iii) = v[j];
185  }
186 
187  for(unsigned int iii = 0; iii < i->GetGeometry().size(); iii++)
188  {
189  const array_1d<double,3>& a = i->GetGeometry()[iii].FastGetSolutionStepValue(ACCELERATION);
190  for(unsigned int j=0; j <TDim; j++)
191  aux_ac(j,iii) = a[j];
192  }
193 
194  //std::cout<<"AUX from PK2 stress FCT"<<aux<<std::endl;
195  GeometryUtils::CalculateGeometryData(i->GetGeometry(),DN_DX, N, Area);
196 
197  if (Area<=0)
198  KRATOS_THROW_ERROR(std::logic_error,"negative area at the moment of estimating the time step","");
199 
200  noalias(Dv_dx) = prod(aux,DN_DX);
201  noalias(Da_dx) = prod(aux_ac,DN_DX);
202 
203  noalias(J) = I + deltatime*Dv_dx +0.5*deltatime*deltatime*Da_dx;
204 
205  detJ = MathUtils<double>::Det(J);
206 // std::cout<<"DETJ is "<<detJ<<std::endl;
207 
208  if (detJ<=0)
209  {
210 
211  deltatime_new=deltatime/(1-detJ); //x=-b/k
212  deltatime=deltatime_new;
213  noalias(J) = I + deltatime*Dv_dx+0.5*deltatime*deltatime*Da_dx;
214  detJ = MathUtils<double>::Det(J);
215  if (detJ<=0)
216  {
217  deltatime_new=deltatime/(1-detJ); //x=-b/k
218  deltatime=deltatime_new;
219  noalias(J) = I + deltatime*Dv_dx+0.5*deltatime*deltatime*Da_dx;
220  detJ = MathUtils<double>::Det(J);
221 
222 
223 
224  if (detJ<=0)
225  {
226  for(unsigned int iii = 0; iii < i->GetGeometry().size(); iii++)
227 
228  deltatime_new=deltatime/(1-detJ); //x=-b/k
229  deltatime=deltatime_new;
230 
231  }
232  }
233 
234  }
235  }
236  KRATOS_WATCH("Estimated delta time is")
237  KRATOS_WATCH(deltatime)
238 
239  return deltatime;
240  KRATOS_CATCH("")
241 
242 
243  }
247 
248 
252 
253 
257 
259  std::string Info() const override
260  {
261  return "UlfTimeStepDecProcess";
262  }
263 
265  void PrintInfo(std::ostream& rOStream) const override
266  {
267  rOStream << "UlfTimeStepDecProcess";
268  }
269 
271  void PrintData(std::ostream& rOStream) const override
272  {
273  }
274 
275 
279 
280 
282 
283 protected:
286 
287 
291 
292 
296 
297 
301 
302 
306 
307 
311 
312 
316 
317 
319 
320 private:
323 
324 
328  ModelPart& mr_model_part;
329 
333 
334 
338 
339 
343 
344 
348 
349 
353 
355 // UlfTimeStepDecProcess& operator=(UlfTimeStepDecProcess const& rOther);
356 
358 // UlfTimeStepDecProcess(UlfTimeStepDecProcess const& rOther);
359 
360 
362 
363 }; // Class UlfTimeStepDecProcess
364 
366 
369 
370 
374 
375 
377 inline std::istream& operator >> (std::istream& rIStream,
378  UlfTimeStepDecProcess& rThis);
379 
381 inline std::ostream& operator << (std::ostream& rOStream,
382  const UlfTimeStepDecProcess& rThis)
383 {
384  rThis.PrintInfo(rOStream);
385  rOStream << std::endl;
386  rThis.PrintData(rOStream);
387 
388  return rOStream;
389 }
391 
392 
393 } // namespace Kratos.
394 
395 #endif // KRATOS_ULF_TIME_STEP_DEC_PROCESS_INCLUDED defined
396 
397 
static void CalculateGeometryData(const GeometryType &rGeometry, BoundedMatrix< double, 4, 3 > &rDN_DX, array_1d< double, 4 > &rN, double &rVolume)
This function is designed to compute the shape function derivatives, shape functions and volume in 3D...
Definition: geometry_utilities.h:176
Definition: amatrix_interface.h:41
static double Det(const TMatrixType &rA)
Calculates the determinant of a matrix of a square matrix of any size (no check performed on release ...
Definition: math_utils.h:597
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
ElementIterator ElementsBegin(IndexType ThisIndex=0)
Definition: model_part.h:1169
ElementIterator ElementsEnd(IndexType ThisIndex=0)
Definition: model_part.h:1179
The base class for all processes in Kratos.
Definition: process.h:49
Short class definition.
Definition: ulf_time_step_dec_process.h:106
KRATOS_CLASS_POINTER_DEFINITION(UlfTimeStepDecProcess)
Pointer definition of PushStructureProcess.
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: ulf_time_step_dec_process.h:271
~UlfTimeStepDecProcess() override
Destructor.
Definition: ulf_time_step_dec_process.h:125
double EstimateDeltaTime(const double dt_max, const double domain_size)
Definition: ulf_time_step_dec_process.h:144
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: ulf_time_step_dec_process.h:265
double EstimateDeltaTimeTemplated(const double dt_max)
Definition: ulf_time_step_dec_process.h:156
UlfTimeStepDecProcess(ModelPart &model_part)
Default constructor.
Definition: ulf_time_step_dec_process.h:119
void operator()()
Definition: ulf_time_step_dec_process.h:134
std::string Info() const override
Turn back information as a string.
Definition: ulf_time_step_dec_process.h:259
#define KRATOS_THROW_ERROR(ExceptionType, ErrorMessage, MoreInfo)
Definition: define.h:77
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_WATCH(variable)
Definition: define.h:806
#define KRATOS_TRY
Definition: define.h:109
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
AMatrix::IdentityMatrix< double > IdentityMatrix
Definition: amatrix_interface.h:564
KratosZeroMatrix< double > ZeroMatrix
Definition: amatrix_interface.h:559
AMatrix::MatrixProductExpression< TExpression1Type, TExpression2Type > prod(AMatrix::MatrixExpression< TExpression1Type, TCategory1 > const &First, AMatrix::MatrixExpression< TExpression2Type, TCategory2 > const &Second)
Definition: amatrix_interface.h:568
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
T & noalias(T &TheMatrix)
Definition: amatrix_interface.h:484
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
int domain_size
Definition: face_heat.py:4
model_part
Definition: face_heat.py:14
v
Definition: generate_convection_diffusion_explicit_element.py:114
a
Definition: generate_stokes_twofluid_element.py:77
int j
Definition: quadrature.py:648
J
Definition: sensitivityMatrix.py:58
N
Definition: sensitivityMatrix.py:29
integer i
Definition: TensorModule.f:17