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_smooth_joint_CL.h
Go to the documentation of this file.
1 // Main author: Chengshun Shang (CIMNE)
3 // Email: cshang@cimne.upc.edu, chengshun.shang1996@gmail.com
4 // Date: June 2023
6 
7 #if !defined(DEM_SMOOTH_JOINT_CL_H_INCLUDE)
8 #define DEM_SMOOTH_JOINT_CL_H_INCLUDE
9 
10 // Project includes
13 
14 namespace Kratos{
15 
16  class KRATOS_API(DEM_APPLICATION) DEM_smooth_joint : public DEMContinuumConstitutiveLaw {
17 
19 
20  public:
21 
23 
25 
26  void Initialize(SphericContinuumParticle* element1, SphericContinuumParticle* element2, Properties::Pointer pProps) override;
27  void TransferParametersToProperties(const Parameters& parameters, Properties::Pointer pProp) override;
28  void Check(Properties::Pointer pProp) const override;
29 
31 
32  DEMContinuumConstitutiveLaw::Pointer Clone() const override;
33 
34  std::string GetTypeOfLaw() override;
35 
36  virtual void CalculateContactArea(double radius, double other_radius, double& calculation_area) override;
37  virtual double CalculateContactArea(double radius, double other_radius, Vector& v) override;
38  void GetContactArea(const double radius, const double other_radius, const Vector& vector_of_initial_areas, const int neighbour_position, double& calculation_area) override;
39  void CalculateElasticConstants(double& kn_el, double& kt_el, double initial_dist, double equiv_young,
40  double equiv_poisson, double calculation_area, SphericContinuumParticle* element1, SphericContinuumParticle* element2, double indentation) override;
41  //virtual void InitializeContact(SphericParticle* const element1, SphericParticle* const element2, const double indentation);
42 
43  double LocalMaxSearchDistance(const int i,
44  SphericContinuumParticle* element1,
45  SphericContinuumParticle* element2) override;
46 
47  //virtual double GetYoungModulusForComputingRotationalMoments(const double& equiv_young);
48 
49  virtual void CheckFailure(const int i_neighbour_count,
50  SphericContinuumParticle* element1,
51  SphericContinuumParticle* element2,
52  double& contact_sigma,
53  double& contact_tau,
54  double LocalElasticContactForce[3],
55  double ViscoDampingLocalContactForce[3],
56  double ElasticLocalRotationalMoment[3],
57  double ViscoLocalRotationalMoment[3]) override;
58 
59  void CalculateForces(const ProcessInfo& r_process_info,
60  double OldLocalElasticContactForce[3],
61  double LocalElasticContactForce[3],
62  double LocalElasticExtraContactForce[3],
63  double LocalCoordSystem[3][3],
64  double LocalDeltDisp[3],
65  const double kn_el,
66  const double kt_el,
67  double& contact_sigma,
68  double& contact_tau,
69  double& failure_criterion_state,
70  double equiv_young,
71  double equiv_shear,
72  double indentation,
73  double calculation_area,
74  double& acumulated_damage,
75  SphericContinuumParticle* element1,
76  SphericContinuumParticle* element2,
77  int i_neighbour_count,
78  int time_steps,
79  bool& sliding,
80  double &equiv_visco_damp_coeff_normal,
81  double &equiv_visco_damp_coeff_tangential,
82  double LocalRelVel[3],
83  double ViscoDampingLocalContactForce[3]) override;
84 
85  void CalculateNormalForces(double LocalElasticContactForce[3],
86  const double kn_el,
87  double equiv_young,
88  double indentation,
89  double calculation_area,
90  double& acumulated_damage,
91  SphericContinuumParticle* element1,
92  SphericContinuumParticle* element2,
93  int i_neighbour_count,
94  int time_steps,
95  const ProcessInfo& r_process_info,
96  double& contact_sigma);
97 
98  virtual void CalculateTangentialForces(double OldLocalElasticContactForce[3],
99  double LocalElasticContactForce[3],
100  double LocalElasticExtraContactForce[3],
101  double ViscoDampingLocalContactForce[3],
102  double LocalCoordSystem[3][3],
103  double LocalDeltDisp[3],
104  double LocalRelVel[3],
105  const double kt_el,
106  const double equiv_shear,
107  double& contact_tau,
108  double indentation,
109  double calculation_area,
110  double& failure_criterion_state,
111  SphericContinuumParticle* element1,
112  SphericContinuumParticle* element2,
113  int i_neighbour_count,
114  bool& sliding,
115  const ProcessInfo& r_process_info,
116  int time_steps);
117 
118  void CalculateMoments(SphericContinuumParticle* element,
119  SphericContinuumParticle* neighbor,
120  double equiv_young,
121  double distance,
122  double calculation_area,
123  double LocalCoordSystem[3][3],
124  double ElasticLocalRotationalMoment[3],
125  double ViscoLocalRotationalMoment[3],
126  double equiv_poisson,
127  double indentation,
128  double LocalElasticContactForce[3],
129  double normalLocalContactForce,
130  double GlobalElasticContactForces[3],
131  double LocalCoordSystem_2[3],
132  const int i_neighbor_count) override;
133 
134  void AddContributionOfShearStrainParallelToBond(double OldLocalElasticContactForce[3],
135  double LocalElasticExtraContactForce[3],
136  array_1d<double, 3>& OldElasticExtraContactForce,
137  double LocalCoordSystem[3][3],
138  const double kt_el,
139  const double calculation_area,
140  SphericContinuumParticle* element1,
141  SphericContinuumParticle* element2);
142 
143 
144  double mAccumulatedJointTangentialLocalDisplacement[2] = {0.0};
145  double mKn;
146  double mKt;
147  double mLocalJointNormal[3] = {0.0};
148  double mInitialDistanceJoint = 0.0;
149 
150  protected:
151 
152  private:
155 
156  };
157 
158 } // namespace Kratos
159 
160 #endif /*DEM_SMOOTH_JOINT_CL_H_INCLUDE defined*/
Definition: DEM_smooth_joint_CL.h:16
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)
double mKn
Definition: DEM_smooth_joint_CL.h:145
KRATOS_CLASS_POINTER_DEFINITION(DEM_smooth_joint)
~DEM_smooth_joint()
Definition: DEM_smooth_joint_CL.h:30
DEM_smooth_joint()
Definition: DEM_smooth_joint_CL.h:24
double mKt
Definition: DEM_smooth_joint_CL.h:146
Definition: DEM_continuum_constitutive_law.h:23
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
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