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_KDEM_CL.h
Go to the documentation of this file.
1 
2 #if !defined(DEM_KDEM_CL_H_INCLUDED)
3 #define DEM_KDEM_CL_H_INCLUDED
4 
5 /* Project includes */
7 
8 
9 namespace Kratos {
10 
11  class KRATOS_API(DEM_APPLICATION) DEM_KDEM : public DEMContinuumConstitutiveLaw {
12 
14 
15  public:
16 
18 
19  DEM_KDEM() {}
20 
21  void TransferParametersToProperties(const Parameters& parameters, Properties::Pointer pProp) override;
22  void Check(Properties::Pointer pProp) const override;
23 
24  ~DEM_KDEM() {}
25 
26  DEMContinuumConstitutiveLaw::Pointer Clone() const override;
27 
28  virtual void CalculateContactArea(double radius, double other_radius, double& calculation_area) override;
29  virtual double CalculateContactArea(double radius, double other_radius, Vector& v) override;
30  void GetContactArea(const double radius, const double other_radius, const Vector& vector_of_initial_areas, const int neighbour_position, double& calculation_area) override;
31  void CalculateElasticConstants(double& kn_el, double& kt_el, double initial_dist, double equiv_young,
32  double equiv_poisson, double calculation_area, SphericContinuumParticle* element1, SphericContinuumParticle* element2, double indentation) override;
33 
34  void CalculateViscoDampingCoeff(double &equiv_visco_damp_coeff_normal,
35  double &equiv_visco_damp_coeff_tangential,
36  SphericContinuumParticle* element1,
37  SphericContinuumParticle* element2,
38  const double kn_el,
39  const double kt_el) override;
40 
41  double LocalMaxSearchDistance(const int i,
42  SphericContinuumParticle* element1,
43  SphericContinuumParticle* element2) override;
44 
45  void CalculateForces(const ProcessInfo& r_process_info,
46  double OldLocalElasticContactForce[3],
47  double LocalElasticContactForce[3],
48  double LocalElasticExtraContactForce[3],
49  double LocalCoordSystem[3][3],
50  double LocalDeltDisp[3],
51  const double kn_el,
52  const double kt_el,
53  double& contact_sigma,
54  double& contact_tau,
55  double& failure_criterion_state,
56  double equiv_young,
57  double equiv_shear,
58  double indentation,
59  double calculation_area,
60  double& acumulated_damage,
61  SphericContinuumParticle* element1,
62  SphericContinuumParticle* element2,
63  int i_neighbour_count,
64  int time_steps,
65  bool& sliding,
66  double &equiv_visco_damp_coeff_normal,
67  double &equiv_visco_damp_coeff_tangential,
68  double LocalRelVel[3],
69  double ViscoDampingLocalContactForce[3]) override;
70 
71  void CalculateNormalForces(double LocalElasticContactForce[3],
72  const double kn_el,
73  double equiv_young,
74  double indentation,
75  double calculation_area,
76  double& acumulated_damage,
77  SphericContinuumParticle* element1,
78  SphericContinuumParticle* element2,
79  int i_neighbour_count,
80  int time_steps,
81  const ProcessInfo& r_process_info) override;
82 
83  double GetContactSigmaMax();
84 
85  virtual double GetYoungModulusForComputingRotationalMoments(const double& equiv_young);
86 
87  void CalculateTangentialForces(double OldLocalElasticContactForce[3],
88  double LocalElasticContactForce[3],
89  double LocalElasticExtraContactForce[3],
90  double ViscoDampingLocalContactForce[3],
91  double LocalCoordSystem[3][3],
92  double LocalDeltDisp[3],
93  double LocalRelVel[3],
94  const double kt_el,
95  const double equiv_shear,
96  double& contact_sigma,
97  double& contact_tau,
98  double indentation,
99  double calculation_area,
100  double& failure_criterion_state,
101  SphericContinuumParticle* element1,
102  SphericContinuumParticle* element2,
103  int i_neighbour_count,
104  bool& sliding,
105  const ProcessInfo& r_process_info) override;
106 
107  void AddContributionOfShearStrainParallelToBond(double OldLocalElasticContactForce[3],
108  double LocalElasticExtraContactForce[3],
109  array_1d<double, 3>& OldElasticExtraContactForce,
110  double LocalCoordSystem[3][3],
111  const double kt_el,
112  const double calculation_area,
113  SphericContinuumParticle* element1,
114  SphericContinuumParticle* element2);
115 
116  void CalculateViscoDamping(double LocalRelVel[3],
117  double ViscoDampingLocalContactForce[3],
118  double indentation,
119  double equiv_visco_damp_coeff_normal,
120  double equiv_visco_damp_coeff_tangential,
121  bool& sliding,
122  int failure_id) override;
123 
124  virtual void CalculateMoments(SphericContinuumParticle* element,
125  SphericContinuumParticle* neighbor,
126  double equiv_young,
127  double distance,
128  double calculation_area,
129  double LocalCoordSystem[3][3],
130  double ElasticLocalRotationalMoment[3],
131  double ViscoLocalRotationalMoment[3],
132  double equiv_poisson,
133  double indentation,
134  double LocalElasticContactForce[3],
135  double normalLocalContactForce,
136  double GlobalElasticContactForces[3],
137  double LocalCoordSystem_2[3],
138  const int i_neighbor_count) override;
139 
140  virtual void ComputeParticleRotationalMoments(SphericContinuumParticle* element,
141  SphericContinuumParticle* neighbor,
142  double equiv_young,
143  double distance,
144  double calculation_area,
145  double LocalCoordSystem[3][3],
146  double ElasticLocalRotationalMoment[3],
147  double ViscoLocalRotationalMoment[3],
148  double equiv_poisson,
149  double indentation,
150  double LocalElasticContactForce[3]) override;
151 
152  void AddPoissonContribution(const double equiv_poisson,
153  double LocalCoordSystem[3][3],
154  double& normal_force,
155  double calculation_area, BoundedMatrix<double, 3, 3>* mSymmStressTensor, SphericContinuumParticle* element1,
156  SphericContinuumParticle* element2, const ProcessInfo& r_process_info, const int i_neighbor_count, const double indentation) override;
157 
158  protected:
159 
160  virtual double GetTauZero(SphericContinuumParticle* element1);
161 
162  virtual double GetInternalFricc(SphericContinuumParticle* element1);
163 
164  private:
165 
166  friend class Serializer;
167 
168  virtual void save(Serializer& rSerializer) const override{
170  //rSerializer.save("MyMemberName",myMember);
171  }
172 
173  virtual void load(Serializer& rSerializer) override{
175  //rSerializer.load("MyMemberName",myMember);
176  }
177  };
178 
179 } /* namespace Kratos.*/
180 #endif /* DEM_KDEM_H_INCLUDED defined */
Definition: DEM_KDEM_CL.h:11
~DEM_KDEM()
Definition: DEM_KDEM_CL.h:24
DEM_KDEM()
Definition: DEM_KDEM_CL.h:19
KRATOS_CLASS_POINTER_DEFINITION(DEM_KDEM)
Definition: DEM_continuum_constitutive_law.h:23
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
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
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