1 // ============================================================================= 2 // PROJECT CHRONO - http://projectchrono.org 3 // 4 // Copyright (c) 2014 projectchrono.org 5 // All rights reserved. 6 // 7 // Use of this source code is governed by a BSD-style license that can be found 8 // in the LICENSE file at the top level of the distribution and at 9 // http://projectchrono.org/license-chrono.txt. 10 // 11 // ============================================================================= 12 // Authors: Alessandro Tasora, Radu Serban 13 // ============================================================================= 14 15 #ifndef CHMESH_H 16 #define CHMESH_H 17 18 #include <cstdlib> 19 #include <cmath> 20 21 #include "chrono/core/ChTimer.h" 22 #include "chrono/physics/ChIndexedNodes.h" 23 #include "chrono/fea/ChContinuumMaterial.h" 24 #include "chrono/fea/ChContactSurface.h" 25 #include "chrono/fea/ChElementBase.h" 26 #include "chrono/fea/ChMeshSurface.h" 27 #include "chrono/fea/ChNodeFEAbase.h" 28 29 namespace chrono { 30 31 class ChAssembly; 32 33 namespace fea { 34 35 /// @addtogroup chrono_fea 36 /// @{ 37 38 /// Class which defines a mesh of finite elements of class ChElementBase, 39 /// between nodes of class ChNodeFEAbase. 40 class ChApi ChMesh : public ChIndexedNodes { 41 42 private: 43 std::vector<std::shared_ptr<ChNodeFEAbase>> vnodes; ///< nodes 44 std::vector<std::shared_ptr<ChElementBase>> velements; ///< elements 45 46 unsigned int n_dofs; ///< total degrees of freedom 47 unsigned int n_dofs_w; ///< total degrees of freedom, derivative (Lie algebra) 48 49 std::vector<std::shared_ptr<ChContactSurface>> vcontactsurfaces; ///< contact surfaces 50 std::vector<std::shared_ptr<ChMeshSurface>> vmeshsurfaces; ///< mesh surfaces, ex.for loads 51 52 bool automatic_gravity_load; 53 int num_points_gravity; 54 55 ChTimer<> timer_internal_forces; 56 ChTimer<> timer_KRMload; 57 int ncalls_internal_forces; 58 int ncalls_KRMload; 59 60 public: ChMesh()61 ChMesh() 62 : n_dofs(0), 63 n_dofs_w(0), 64 automatic_gravity_load(true), 65 num_points_gravity(1), 66 ncalls_internal_forces(0), 67 ncalls_KRMload(0) {} 68 ChMesh(const ChMesh& other); ~ChMesh()69 ~ChMesh() {} 70 71 /// "Virtual" copy constructor (covariant return type). Clone()72 virtual ChMesh* Clone() const override { return new ChMesh(*this); } 73 74 void AddNode(std::shared_ptr<ChNodeFEAbase> m_node); 75 void AddElement(std::shared_ptr<ChElementBase> m_elem); 76 void ClearNodes(); 77 void ClearElements(); 78 79 /// Get the array of nodes of this mesh. GetNodes()80 const std::vector<std::shared_ptr<ChNodeFEAbase>>& GetNodes() const { return vnodes; } 81 82 /// Get the array of elements of this mesh. GetElements()83 const std::vector<std::shared_ptr<ChElementBase>>& GetElements() const { return velements; } 84 85 /// Access the N-th node GetNode(unsigned int n)86 virtual std::shared_ptr<ChNodeBase> GetNode(unsigned int n) override { return vnodes[n]; } 87 88 /// Access the N-th element GetElement(unsigned int n)89 std::shared_ptr<ChElementBase> GetElement(unsigned int n) { return velements[n]; } 90 91 /// Get the number of nodes in the mesh. GetNnodes()92 virtual unsigned int GetNnodes() const override { return (unsigned int)vnodes.size(); } 93 94 /// Get the number of elements in the mesh. GetNelements()95 unsigned int GetNelements() { return (unsigned int)velements.size(); } 96 GetDOF()97 virtual int GetDOF() override { return n_dofs; } GetDOF_w()98 virtual int GetDOF_w() override { return n_dofs_w; } 99 100 /// Override default in ChPhysicsItem. GetCollide()101 virtual bool GetCollide() const override { return true; } 102 103 /// Reset counters for internal force and Jacobian evaluations. ResetCounters()104 void ResetCounters() { 105 ncalls_internal_forces = 0; 106 ncalls_KRMload = 0; 107 } 108 /// Get cumulative number of calls to internal forces evaluation. GetNumCallsInternalForces()109 int GetNumCallsInternalForces() { return ncalls_internal_forces; } 110 /// Get cumulative number of calls to load Jacobian information. GetNumCallsJacobianLoad()111 int GetNumCallsJacobianLoad() { return ncalls_KRMload; } 112 113 /// Reset timers for internal force and Jacobian evaluations. ResetTimers()114 void ResetTimers() { 115 timer_internal_forces.reset(); 116 timer_KRMload.reset(); 117 } 118 /// Get cumulative time for internal force evaluation. GetTimeInternalForces()119 double GetTimeInternalForces() { return timer_internal_forces(); } 120 /// Get cumulative time for Jacobian load calls. GetTimeJacobianLoad()121 double GetTimeJacobianLoad() { return timer_KRMload(); } 122 123 /// Add a contact surface. 124 void AddContactSurface(std::shared_ptr<ChContactSurface> m_surf); 125 126 /// Access the N-th contact surface. GetContactSurface(unsigned int n)127 std::shared_ptr<ChContactSurface> GetContactSurface(unsigned int n) { return vcontactsurfaces[n]; } 128 129 /// Get number of added contact surfaces. GetNcontactSurfaces()130 unsigned int GetNcontactSurfaces() { return (unsigned int)vcontactsurfaces.size(); } 131 132 /// Remove all contact surfaces. 133 void ClearContactSurfaces(); 134 135 /// Add a mesh surface (set of ChLoadableUV items that support area loads such as pressure). 136 void AddMeshSurface(std::shared_ptr<ChMeshSurface> m_surf); 137 138 /// Access the N-th mesh surface. GetMeshSurface(unsigned int n)139 std::shared_ptr<ChMeshSurface> GetMeshSurface(unsigned int n) { return vmeshsurfaces[n]; } 140 141 /// Get number of added mesh surfaces GetNmeshSurfaces()142 unsigned int GetNmeshSurfaces() { return (unsigned int)vmeshsurfaces.size(); } 143 144 /// Remove all mesh surfaces. ClearMeshSurfaces()145 void ClearMeshSurfaces() { vmeshsurfaces.clear(); } 146 147 /// Set reference position of nodes as current position, for all nodes. 148 void Relax(); 149 150 /// Set no speed and no accelerations in nodes (but does not change reference positions). 151 void SetNoSpeedNoAcceleration() override; 152 153 /// This recomputes the number of DOFs, constraints, 154 /// as well as state offsets of contained items. 155 virtual void Setup() override; 156 157 /// Update time dependent data, for all elements. 158 /// Updates all [A] coord.systems for all (corotational) elements. 159 virtual void Update(double m_time, bool update_assets = true) override; 160 161 // Functions to interface this with ChPhysicsItem container 162 virtual void SyncCollisionModels() override; 163 virtual void AddCollisionModelsToSystem() override; 164 virtual void RemoveCollisionModelsFromSystem() override; 165 166 /// If true, as by default, this mesh will add automatically a gravity load 167 /// to all contained elements (that support gravity) using the G value from the ChSystem. 168 /// So this saves you from adding many ChLoad<ChLoaderGravity> to all elements. 169 void SetAutomaticGravity(bool mg, int num_points = 1) { 170 automatic_gravity_load = mg; 171 num_points_gravity = num_points; 172 } 173 /// Tell if this mesh will add automatically a gravity load to all contained elements. GetAutomaticGravity()174 bool GetAutomaticGravity() { return automatic_gravity_load; } 175 176 /// Get ChMesh mass properties 177 void ComputeMassProperties(double& mass, ///< ChMesh object mass 178 ChVector<>& com, ///< ChMesh center of gravity 179 ChMatrix33<>& inertia ///< ChMesh inertia tensor 180 ); 181 182 // 183 // STATE FUNCTIONS 184 // 185 186 // (override/implement interfaces for global state vectors, see ChPhysicsItem for comments.) 187 virtual void IntStateGather(const unsigned int off_x, 188 ChState& x, 189 const unsigned int off_v, 190 ChStateDelta& v, 191 double& T) override; 192 virtual void IntStateScatter(const unsigned int off_x, 193 const ChState& x, 194 const unsigned int off_v, 195 const ChStateDelta& v, 196 const double T, 197 bool full_update) override; 198 virtual void IntStateGatherAcceleration(const unsigned int off_a, ChStateDelta& a) override; 199 virtual void IntStateScatterAcceleration(const unsigned int off_a, const ChStateDelta& a) override; 200 virtual void IntStateIncrement(const unsigned int off_x, 201 ChState& x_new, 202 const ChState& x, 203 const unsigned int off_v, 204 const ChStateDelta& Dv) override; 205 virtual void IntLoadResidual_F(const unsigned int off, ChVectorDynamic<>& R, const double c) override; 206 virtual void IntLoadResidual_Mv(const unsigned int off, 207 ChVectorDynamic<>& R, 208 const ChVectorDynamic<>& w, 209 const double c) override; 210 virtual void IntToDescriptor(const unsigned int off_v, 211 const ChStateDelta& v, 212 const ChVectorDynamic<>& R, 213 const unsigned int off_L, 214 const ChVectorDynamic<>& L, 215 const ChVectorDynamic<>& Qc) override; 216 virtual void IntFromDescriptor(const unsigned int off_v, 217 ChStateDelta& v, 218 const unsigned int off_L, 219 ChVectorDynamic<>& L) override; 220 221 // 222 // SYSTEM FUNCTIONS for interfacing all elements with solver 223 // 224 225 /// Tell to a system descriptor that there are items of type 226 /// ChKblock in this object (for further passing it to a solver) 227 /// Basically does nothing, but maybe that inherited classes may specialize this. 228 virtual void InjectKRMmatrices(ChSystemDescriptor& mdescriptor) override; 229 230 /// Adds the current stiffness K and damping R and mass M matrices in encapsulated 231 /// ChKblock item(s), if any. The K, R, M matrices are added with scaling 232 /// values Kfactor, Rfactor, Mfactor. 233 virtual void KRMmatricesLoad(double Kfactor, double Rfactor, double Mfactor) override; 234 235 /// Sets the 'fb' part (the known term) of the encapsulated ChVariables to zero. 236 virtual void VariablesFbReset() override; 237 238 /// Adds the current forces (applied to item) into the 239 /// encapsulated ChVariables, in the 'fb' part: qf+=forces*factor. 240 virtual void VariablesFbLoadForces(double factor = 1) override; 241 242 /// Initialize the 'qb' part of the ChVariables with the 243 /// current value of speeds. Note: since 'qb' is the unknown, this 244 /// function seems unnecessary, unless used before VariablesFbIncrementMq(). 245 virtual void VariablesQbLoadSpeed() override; 246 247 /// Adds M*q (masses multiplied current 'qb') to Fb, ex. if qb is initialized 248 /// with v_old using VariablesQbLoadSpeed, this method can be used in 249 /// timestepping schemes that do: M*v_new = M*v_old + forces*dt. 250 virtual void VariablesFbIncrementMq() override; 251 252 /// Fetches the item speed (ex. linear and angular vel.in rigid bodies) from the 253 /// 'qb' part of the ChVariables and sets it as the current item speed. 254 /// If 'step' is not 0, also should compute the approximate acceleration of 255 /// the item using backward differences, that is accel=(new_speed-old_speed)/step. 256 /// Mostly used after the solver provided the solution in ChVariables. 257 virtual void VariablesQbSetSpeed(double step = 0) override; 258 259 /// Increment item positions by the 'qb' part of the ChVariables, 260 /// multiplied by a 'step' factor. 261 /// <pre> 262 /// pos += qb * step 263 /// </pre> 264 /// If qb is a speed, this behaves like a single step of 1-st order 265 /// numerical integration (Eulero integration). 266 virtual void VariablesQbIncrementPosition(double step) override; 267 268 /// Tell to a system descriptor that there are variables of type 269 /// ChVariables in this object (for further passing it to a solver) 270 /// Basically does nothing, but maybe that inherited classes may specialize this. 271 virtual void InjectVariables(ChSystemDescriptor& mdescriptor) override; 272 273 private: 274 /// Initial setup (before analysis). 275 /// This function is called from ChSystem::SetupInitial, marking a point where system 276 /// construction is completed. 277 /// <pre> 278 /// - Computes the total number of degrees of freedom 279 /// - Precompute auxiliary data, such as (local) stiffness matrices Kl, if any, for each element. 280 /// </pre> 281 virtual void SetupInitial() override; 282 283 friend class chrono::ChSystem; 284 friend class chrono::ChAssembly; 285 }; 286 287 /// @} chrono_fea 288 289 } // end namespace fea 290 } // end namespace chrono 291 292 #endif 293