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_parallel_bond_CL.h
Go to the documentation of this file.
1 // Main author: Chengshun Shang (CIMNE)
3 // Email: chengshun.shang1996@gmail.com
4 // Date: June 2022
6 
7 #if !defined(DEM_PARALLEL_BOND_CL_H_INCLUDE)
8 #define DEM_PARALLEL_BOND_CL_H_INCLUDE
9 
10 // Project includes
13 
14 namespace Kratos{
15 
16  class KRATOS_API(DEM_APPLICATION) DEM_parallel_bond : public DEMContinuumConstitutiveLaw {
17 
19 
20  public:
21 
22  // TODO: what is the function?
24 
26 
27  void TransferParametersToProperties(const Parameters& parameters, Properties::Pointer pProp) override;
28  std::string GetTypeOfLaw() override;
29  void Check(Properties::Pointer pProp) const override;
30 
32 
33  DEMContinuumConstitutiveLaw::Pointer Clone() const override;
34 
35  virtual void CalculateContactArea(double radius, double other_radius, double& calculation_area) override;
36  virtual double CalculateContactArea(double radius, double other_radius, Vector& v) override;
37  void GetContactArea(const double radius, const double other_radius, const Vector& vector_of_initial_areas, const int neighbour_position, double& calculation_area) override;
38  void CalculateElasticConstants(double& kn_el, double& kt_el, double initial_dist, double equiv_young,
39  double equiv_poisson, double calculation_area, SphericContinuumParticle* element1, SphericContinuumParticle* element2, double indentation) override;
40  virtual void InitializeContact(SphericParticle* const element1, SphericParticle* const element2, const double indentation);
41 
42  virtual void CalculateUnbondedViscoDampingForce(double LocalRelVel[3],
43  double UnbondedViscoDampingLocalContactForce[3],
44  SphericParticle* const element1,
45  SphericParticle* const element2);
46 
47  // TODO: check whether it is necessary
48  double LocalMaxSearchDistance(const int i,
49  SphericContinuumParticle* element1,
50  SphericContinuumParticle* element2) override;
51  //double GetContactSigmaMax();
52 
53  //TODO:CHECK
54  virtual double GetYoungModulusForComputingRotationalMoments(const double& equiv_young);
55 
56  virtual void CheckFailure(const int i_neighbour_count,
57  SphericContinuumParticle* element1,
58  SphericContinuumParticle* element2,
59  double& contact_sigma,
60  double& contact_tau,
61  double LocalElasticContactForce[3],
62  double ViscoDampingLocalContactForce[3],
63  double ElasticLocalRotationalMoment[3],
64  double ViscoLocalRotationalMoment[3]) override;
65 
66  void CalculateForces(const ProcessInfo& r_process_info,
67  double OldLocalElasticContactForce[3],
68  double LocalElasticContactForce[3],
69  double LocalElasticExtraContactForce[3],
70  double LocalCoordSystem[3][3],
71  double LocalDeltDisp[3],
72  const double kn_el,
73  const double kt_el,
74  double& contact_sigma,
75  double& contact_tau,
76  double& failure_criterion_state,
77  double equiv_young,
78  double equiv_shear,
79  double indentation,
80  double calculation_area,
81  double& acumulated_damage,
82  SphericContinuumParticle* element1,
83  SphericContinuumParticle* element2,
84  int i_neighbour_count,
85  int time_steps,
86  bool& sliding,
87  double &equiv_visco_damp_coeff_normal,
88  double &equiv_visco_damp_coeff_tangential,
89  double LocalRelVel[3],
90  double ViscoDampingLocalContactForce[3]) override;
91 
92  virtual double ComputeNormalUnbondedForce(double unbonded_indentation);
93 
94  void CalculateNormalForces(double LocalElasticContactForce[3],
95  const double kn_el,
96  double equiv_young,
97  double indentation,
98  double calculation_area,
99  double& acumulated_damage,
100  SphericContinuumParticle* element1,
101  SphericContinuumParticle* element2,
102  int i_neighbour_count,
103  int time_steps,
104  const ProcessInfo& r_process_info,
105  double& contact_sigma);
106 
107  virtual void CalculateViscoDampingCoeff(double &equiv_visco_damp_coeff_normal,
108  double &equiv_visco_damp_coeff_tangential,
109  SphericContinuumParticle* element1,
110  SphericContinuumParticle* element2,
111  const double kn_el,
112  const double kt_el) override;
113 
114  void CalculateViscoDamping(double LocalRelVel[3],
115  double ViscoDampingLocalContactForce[3],
116  double indentation,
117  double equiv_visco_damp_coeff_normal,
118  double equiv_visco_damp_coeff_tangential,
119  bool& sliding,
120  int failure_id,
121  int i_neighbour_count,
122  SphericContinuumParticle* element1,
123  SphericContinuumParticle* element2);
124 
125  virtual void CalculateTangentialForces(double OldLocalElasticContactForce[3],
126  double LocalElasticContactForce[3],
127  double LocalElasticExtraContactForce[3],
128  double ViscoDampingLocalContactForce[3],
129  double LocalCoordSystem[3][3],
130  double LocalDeltDisp[3],
131  double LocalRelVel[3],
132  const double kt_el,
133  const double equiv_shear,
134  double& contact_tau,
135  double indentation,
136  double calculation_area,
137  double& failure_criterion_state,
138  SphericContinuumParticle* element1,
139  SphericContinuumParticle* element2,
140  int i_neighbour_count,
141  bool& sliding,
142  const ProcessInfo& r_process_info);
143 
144  void CalculateMoments(SphericContinuumParticle* element,
145  SphericContinuumParticle* neighbor,
146  double equiv_young,
147  double distance,
148  double calculation_area,
149  double LocalCoordSystem[3][3],
150  double ElasticLocalRotationalMoment[3],
151  double ViscoLocalRotationalMoment[3],
152  double equiv_poisson,
153  double indentation,
154  double LocalElasticContactForce[3],
155  double normalLocalContactForce,
156  double GlobalElasticContactForces[3],
157  double LocalCoordSystem_2[3],
158  const int i_neighbor_count) override;
159 
160  void ComputeParticleRotationalMoments(SphericContinuumParticle* element,
161  SphericContinuumParticle* neighbor,
162  double equiv_young,
163  double distance,
164  double calculation_area,
165  double LocalCoordSystem[3][3],
166  double ElasticLocalRotationalMoment[3],
167  double ViscoLocalRotationalMoment[3],
168  double equiv_poisson,
169  double indentation,
170  double LocalElasticContactForce[3]) override;
171 
172  void AddContributionOfShearStrainParallelToBond(double OldLocalElasticContactForce[3],
173  double LocalElasticExtraContactForce[3],
174  array_1d<double, 3>& OldElasticExtraContactForce,
175  double LocalCoordSystem[3][3],
176  const double kt_el,
177  const double calculation_area,
178  SphericContinuumParticle* element1,
179  SphericContinuumParticle* element2);
180 
181 
182  double mUnbondedLocalElasticContactForce2 = 0.0;
183  double mUnbondedNormalElasticConstant = 0.0;
184  double mUnbondedTangentialElasticConstant = 0.0;
185  double mUnbondedViscoDampingLocalContactForce[3] = {0.0};
186  double mBondedViscoDampingLocalContactForce[3] = {0.0};
187  double mBondedScalingFactor[3] = {0.0};
188  double mUnbondedEquivViscoDampCoeffTangential = 0.0;
189  double mUnbondedEquivViscoDampCoeffNormal = 0.0;
190  double mInitialIndentationForBondedPart = 0.0;
191  double mAccumulatedBondedTangentialLocalDisplacement[2] = {0.0};
192  double mBondedLocalContactNormalTorque[3] = {0.0};
193  double mBondedLocalContactTangentTorque[3] = {0.0};
194  double mKn;
195  double mKt;
196 
197  protected:
198 
199  private:
203 
204  };
205 
206 } // namespace Kratos
207 
208 #endif /*DEM_PARALLEL_BOND_CL_H_INCLUDE defined*/
Definition: DEM_parallel_bond_CL.h:16
~DEM_parallel_bond()
Definition: DEM_parallel_bond_CL.h:31
double mKt
Definition: DEM_parallel_bond_CL.h:195
void AddContributionOfShearStrainParallelToBond(double OldLocalElasticContactForce[3], double LocalElasticExtraContactForce[3], array_1d< double, 3 > &OldElasticExtraContactForce, double LocalCoordSystem[3][3], const double kt_el, const double calculation_area, SphericContinuumParticle *element1, SphericContinuumParticle *element2)
KRATOS_CLASS_POINTER_DEFINITION(DEM_parallel_bond)
DEM_parallel_bond()
Definition: DEM_parallel_bond_CL.h:25
double mKn
Definition: DEM_parallel_bond_CL.h:194
Definition: DEM_continuum_constitutive_law.h:23
virtual void CalculateViscoDamping(double LocalRelVel[3], double ViscoDampingLocalContactForce[3], double indentation, double equiv_visco_damp_coeff_normal, double equiv_visco_damp_coeff_tangential, bool &sliding, int failure_id)
Definition: DEM_continuum_constitutive_law.cpp:50
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 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
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
Definition: spheric_continuum_particle.h:26
Definition: spheric_particle.h:31
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
integer i
Definition: TensorModule.f:17
Definition: mesh_converter.cpp:33