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.
DEM_continuum_constitutive_law.h
Go to the documentation of this file.
1 
2 #if !defined(DEM_CONTINUUM_CONSTITUTIVE_LAW_H_INCLUDED)
3 #define DEM_CONTINUUM_CONSTITUTIVE_LAW_H_INCLUDED
4 
5 /* Project includes */
6 #include "includes/define.h"
7 #include "custom_utilities/AuxiliaryFunctions.h"
8 #include "includes/serializer.h"
9 //#include "includes/properties.h"
10 #include "containers/flags.h"
11 
15 #include "containers/array_1d.h"
16 
17 
18 namespace Kratos {
19 
20  class Properties; //forward declaration
21  class SphericContinuumParticle; // forward declaration of spheric cont particle
22 
23  class KRATOS_API(DEM_APPLICATION) DEMContinuumConstitutiveLaw : public Flags {
24 
25  public:
26 
28 
30 
31  DEMContinuumConstitutiveLaw(const DEMContinuumConstitutiveLaw& rReferenceContinuumConstitutiveLaw);
32 
33  virtual void Initialize(SphericContinuumParticle* element1, SphericContinuumParticle* element2, Properties::Pointer pProps);
34 
35  virtual void SetConstitutiveLawInProperties(Properties::Pointer pProp, bool verbose = true);
36 
37  virtual void SetConstitutiveLawInPropertiesWithParameters(Properties::Pointer pProp, const Parameters& parameters, bool verbose = true);
38 
39  virtual void TransferParametersToProperties(const Parameters& parameters, Properties::Pointer pProp);
40 
41  virtual void Check(Properties::Pointer pProp) const;
42 
43  virtual std::string GetTypeOfLaw();
44 
45  virtual ~DEMContinuumConstitutiveLaw();
46 
47  virtual DEMContinuumConstitutiveLaw::Pointer Clone() const;
48 
49  virtual void CalculateViscoDamping(double LocalRelVel[3],
50  double ViscoDampingLocalContactForce[3],
51  double indentation,
52  double equiv_visco_damp_coeff_normal,
53  double equiv_visco_damp_coeff_tangential,
54  bool& sliding,
55  int failure_id);
56 
57  virtual void CalculateContactArea(double radius,
58  double other_radius,
59  double& calculation_area) {
60  KRATOS_ERROR << "This function (DEMContinuumConstitutiveLaw::CalculateContactArea) shouldn't be accessed, use derived class instead"<<std::endl;
61  };
62 
63  virtual double CalculateContactArea(double radius,
64  double other_radius,
65  Vector& v) {
66  return 0.0;
67  }
68 
69  virtual void GetContactArea(const double radius,
70  const double other_radius,
71  const Vector& vector_of_initial_areas,
72  const int neighbour_position,
73  double& calculation_area){
74 
75  CalculateContactArea(radius, other_radius, calculation_area);
76 
77  }
78 
79  virtual void CalculateElasticConstants(double &kn_el,
80  double &kt_el,
81  double initial_dist,
82  double equiv_young,
83  double equiv_poisson,
84  double calculation_area,
85  SphericContinuumParticle* element1,
86  SphericContinuumParticle* element2, double indentation) {
87  KRATOS_ERROR << "This function (DEMContinuumConstitutiveLaw::CalculateElasticConstants) shouldn't be accessed, use derived class instead"<<std::endl;
88  };
89 
90  virtual void CalculateViscoDampingCoeff(double &equiv_visco_damp_coeff_normal,
91  double &equiv_visco_damp_coeff_tangential,
92  SphericContinuumParticle* element1,
93  SphericContinuumParticle* element2,
94  const double kn_el,
95  const double kt_el) {
96  KRATOS_ERROR << "This function (DEMContinuumConstitutiveLaw::CalculateViscoDampingCoeff) shouldn't be accessed, use derived class instead"<<std::endl;
97  };
98 
99  /*
100  virtual void CheckFailure(const int i_neighbour_count, SphericContinuumParticle* element1, SphericContinuumParticle* element2) {
101 
102  };
103  */
104 
105  virtual void CheckFailure(const int i_neighbour_count,
106  SphericContinuumParticle* element1,
107  SphericContinuumParticle* element2,
108  double& contact_sigma,
109  double& contact_tau,
110  double LocalElasticContactForce[3],
111  double ViscoDampingLocalContactForce[3],
112  double ElasticLocalRotationalMoment[3],
113  double ViscoLocalRotationalMoment[3]) {
114 
115  };
116 
117  virtual void CalculateForces(const ProcessInfo& r_process_info,
118  double OldLocalElasticContactForce[3],
119  double LocalElasticContactForce[3],
120  double LocalElasticExtraContactForce[3],
121  double LocalCoordSystem[3][3],
122  double LocalDeltDisp[3],
123  const double kn_el,
124  const double kt_el,
125  double& contact_sigma,
126  double& contact_tau,
127  double& failure_criterion_state,
128  double equiv_young,
129  double equiv_shear,
130  double indentation,
131  double calculation_area,
132  double& acumulated_damage,
133  SphericContinuumParticle* element1,
134  SphericContinuumParticle* element2,
135  int i_neighbour_count,
136  int time_steps,
137  bool& sliding,
138  double &equiv_visco_damp_coeff_normal,
139  double &equiv_visco_damp_coeff_tangential,
140  double LocalRelVel[3],
141  double ViscoDampingLocalContactForce[3]) {
142  KRATOS_ERROR << "This function (DEMContinuumConstitutiveLaw::CalculateForces) shouldn't be accessed, use derived class instead"<<std::endl;
143  };
144 
145  virtual void CalculateNormalForces(double LocalElasticContactForce[3],
146  const double kn_el,
147  double equiv_young,
148  double indentation,
149  double calculation_area,
150  double& acumulated_damage,
151  SphericContinuumParticle* element1,
152  SphericContinuumParticle* element2,
153  int i_neighbour_count,
154  int time_steps,
155  const ProcessInfo& r_process_info) {
156  KRATOS_ERROR << "This function (DEMContinuumConstitutiveLaw::CalculateNormalForces) shouldn't be accessed, use derived class instead"<<std::endl;
157  }
158 
159  virtual void CalculateTangentialForces(double OldLocalElasticContactForce[3],
160  double LocalElasticContactForce[3],
161  double LocalElasticExtraContactForce[3],
162  double ViscoDampingLocalContactForce[3],
163  double LocalCoordSystem[3][3],
164  double LocalDeltDisp[3],
165  double LocalRelVel[3],
166  const double kt_el,
167  const double equiv_shear,
168  double& contact_sigma,
169  double& contact_tau,
170  double indentation,
171  double calculation_area,
172  double& failure_criterion_state,
173  SphericContinuumParticle* element1,
174  SphericContinuumParticle* element2,
175  int i_neighbour_count,
176  bool& sliding,
177  const ProcessInfo& r_process_info) {
178  KRATOS_ERROR << "This function (DEMContinuumConstitutiveLaw::CalculateTangentialForces) shouldn't be accessed, use derived class instead"<<std::endl;
179  };
180 
182  SphericContinuumParticle* neighbor,
183  double equiv_young,
184  double distance,
185  double calculation_area,
186  double LocalCoordSystem[3][3],
187  double ElasticLocalRotationalMoment[3],
188  double ViscoLocalRotationalMoment[3],
189  double equiv_poisson,
190  double indentation,
191  double LocalElasticContactForce[3],
192  double normalLocalContactForce,
193  double GlobalElasticContactForces[3],
194  double LocalCoordSystem_2[3],
195  const int i_neighbor_count){
196  KRATOS_ERROR << "This function (DEMContinuumConstitutiveLaw::CalculateMoments) shouldn't be accessed, use derived class instead"<<std::endl;
197  };
198 
199  virtual void ComputeParticleRotationalMoments(SphericContinuumParticle* element,
200  SphericContinuumParticle* neighbor,
201  double equiv_young,
202  double distance,
203  double calculation_area,
204  double LocalCoordSystem[3][3],
205  double ElasticLocalRotationalMoment[3],
206  double ViscoLocalRotationalMoment[3],
207  double equiv_poisson,
208  double indentation,
209  double LocalElasticContactForce[3]);
210 
211  virtual void AddPoissonContribution(const double equiv_poisson,
212  double LocalCoordSystem[3][3],
213  double& normal_force,
214  double calculation_area,
215  BoundedMatrix<double, 3, 3>* mSymmStressTensor,
216  SphericContinuumParticle* element1,
217  SphericContinuumParticle* element2,
218  const ProcessInfo& r_process_info,
219  const int i_neighbor_count,
220  const double indentation);
221 
222  virtual double LocalMaxSearchDistance(const int i,
223  SphericContinuumParticle* element1,
224  SphericContinuumParticle* element2);
225 
226  virtual bool CheckRequirementsOfStressTensor();
227 
228  protected:
229  Properties::Pointer mpProperties;
230 
231  private:
232 
233  friend class Serializer;
234 
235  virtual void save(Serializer& rSerializer) const override {
237  //rSerializer.save("MyMemberName",myMember);
238  }
239 
240  virtual void load(Serializer& rSerializer) override {
242  //rSerializer.load("MyMemberName",myMember);
243  }
244  };
245 
246  //This definition is done here to avoid recursive inclusion of header files
247  KRATOS_DEFINE_APPLICATION_VARIABLE(DEM_APPLICATION, DEMContinuumConstitutiveLaw::Pointer, DEM_CONTINUUM_CONSTITUTIVE_LAW_POINTER)
248 
249 } /* namespace Kratos.*/
250 #endif /* DEM_CONSTITUTIVE_LAW_H_INCLUDED defined */
251 
Definition: DEM_continuum_constitutive_law.h:23
virtual void CalculateMoments(SphericContinuumParticle *element, SphericContinuumParticle *neighbor, double equiv_young, double distance, double calculation_area, double LocalCoordSystem[3][3], double ElasticLocalRotationalMoment[3], double ViscoLocalRotationalMoment[3], double equiv_poisson, double indentation, double LocalElasticContactForce[3], double normalLocalContactForce, double GlobalElasticContactForces[3], double LocalCoordSystem_2[3], const int i_neighbor_count)
Definition: DEM_continuum_constitutive_law.h:181
virtual void CalculateForces(const ProcessInfo &r_process_info, double OldLocalElasticContactForce[3], double LocalElasticContactForce[3], double LocalElasticExtraContactForce[3], double LocalCoordSystem[3][3], double LocalDeltDisp[3], const double kn_el, const double kt_el, double &contact_sigma, double &contact_tau, double &failure_criterion_state, double equiv_young, double equiv_shear, double indentation, double calculation_area, double &acumulated_damage, SphericContinuumParticle *element1, SphericContinuumParticle *element2, int i_neighbour_count, int time_steps, bool &sliding, double &equiv_visco_damp_coeff_normal, double &equiv_visco_damp_coeff_tangential, double LocalRelVel[3], double ViscoDampingLocalContactForce[3])
Definition: DEM_continuum_constitutive_law.h:117
virtual void CalculateContactArea(double radius, double other_radius, double &calculation_area)
Definition: DEM_continuum_constitutive_law.h:57
virtual void CheckFailure(const int i_neighbour_count, SphericContinuumParticle *element1, SphericContinuumParticle *element2, double &contact_sigma, double &contact_tau, double LocalElasticContactForce[3], double ViscoDampingLocalContactForce[3], double ElasticLocalRotationalMoment[3], double ViscoLocalRotationalMoment[3])
Definition: DEM_continuum_constitutive_law.h:105
KRATOS_CLASS_POINTER_DEFINITION(DEMContinuumConstitutiveLaw)
virtual void GetContactArea(const double radius, const double other_radius, const Vector &vector_of_initial_areas, const int neighbour_position, double &calculation_area)
Definition: DEM_continuum_constitutive_law.h:69
Properties::Pointer mpProperties
Definition: DEM_continuum_constitutive_law.h:229
virtual double CalculateContactArea(double radius, double other_radius, Vector &v)
Definition: DEM_continuum_constitutive_law.h:63
virtual void CalculateViscoDampingCoeff(double &equiv_visco_damp_coeff_normal, double &equiv_visco_damp_coeff_tangential, SphericContinuumParticle *element1, SphericContinuumParticle *element2, const double kn_el, const double kt_el)
Definition: DEM_continuum_constitutive_law.h:90
virtual void CalculateNormalForces(double LocalElasticContactForce[3], const double kn_el, double equiv_young, double indentation, double calculation_area, double &acumulated_damage, SphericContinuumParticle *element1, SphericContinuumParticle *element2, int i_neighbour_count, int time_steps, const ProcessInfo &r_process_info)
Definition: DEM_continuum_constitutive_law.h:145
virtual void CalculateElasticConstants(double &kn_el, double &kt_el, double initial_dist, double equiv_young, double equiv_poisson, double calculation_area, SphericContinuumParticle *element1, SphericContinuumParticle *element2, double indentation)
Definition: DEM_continuum_constitutive_law.h:79
virtual void CalculateTangentialForces(double OldLocalElasticContactForce[3], double LocalElasticContactForce[3], double LocalElasticExtraContactForce[3], double ViscoDampingLocalContactForce[3], double LocalCoordSystem[3][3], double LocalDeltDisp[3], double LocalRelVel[3], const double kt_el, const double equiv_shear, double &contact_sigma, double &contact_tau, double indentation, double calculation_area, double &failure_criterion_state, SphericContinuumParticle *element1, SphericContinuumParticle *element2, int i_neighbour_count, bool &sliding, const ProcessInfo &r_process_info)
Definition: DEM_continuum_constitutive_law.h:159
Definition: flags.h:58
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
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
Definition: spheric_continuum_particle.h:26
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
#define KRATOS_ERROR
Definition: exception.h:161
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
KRATOS_DEFINE_APPLICATION_VARIABLE(CHIMERA_APPLICATION, double, CHIMERA_DISTANCE)
parameters
Definition: fluid_chimera_analysis.py:35
v
Definition: generate_convection_diffusion_explicit_element.py:114
float radius
Definition: mesh_to_mdpa_converter.py:18
def load(f)
Definition: ode_solve.py:307
integer i
Definition: TensorModule.f:17
Definition: mesh_converter.cpp:33