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: Andrea Favali, Alessandro Tasora, Radu Serban 13 // ============================================================================= 14 15 #ifndef CHELEMENTGENERIC_H 16 #define CHELEMENTGENERIC_H 17 18 #include "chrono/solver/ChKblockGeneric.h" 19 #include "chrono/solver/ChVariablesNode.h" 20 #include "chrono/fea/ChElementBase.h" 21 22 namespace chrono { 23 namespace fea { 24 25 /// @addtogroup fea_elements 26 /// @{ 27 28 /// Class for all elements whose stiffness matrix can be seen 29 /// as a NxN block-matrix to be splitted between N nodes. 30 /// Helps reducing the complexity of inherited FEA elements because 31 /// it implements some bookkeeping for the interface with solver. 32 /// This means that most FEA elements inherited from ChElementGeneric 33 /// need to implement at most the following two fundamental methods: 34 /// ComputeKRMmatricesGlobal(), ComputeInternalForces(), 35 /// and optionally ComputeGravityForces(). 36 class ChApi ChElementGeneric : public ChElementBase { 37 protected: 38 ChKblockGeneric Kmatr; 39 40 public: ChElementGeneric()41 ChElementGeneric() {} ~ChElementGeneric()42 virtual ~ChElementGeneric() {} 43 44 /// Access the proxy to stiffness, for sparse solver Kstiffness()45 ChKblockGeneric& Kstiffness() { return Kmatr; } 46 47 // 48 // Functions for interfacing to the state bookkeeping 49 // 50 51 52 53 /// (This is a default (a bit unoptimal) book keeping so that in children classes you can avoid 54 /// implementing this EleIntLoadResidual_F function, unless you need faster code) 55 virtual void EleIntLoadResidual_F(ChVectorDynamic<>& R, const double c) override; 56 57 /// (This is a default (VERY UNOPTIMAL) book keeping so that in children classes you can avoid 58 /// implementing this EleIntLoadResidual_Mv function, unless you need faster code.) 59 virtual void EleIntLoadResidual_Mv(ChVectorDynamic<>& R, const ChVectorDynamic<>& w, const double c) override; 60 61 /// (This is a default (VERY UNOPTIMAL) book keeping so that in children classes you can avoid 62 /// implementing this EleIntLoadResidual_F_gravity function, unless you need faster code. 63 /// This fallback implementation uses a temp ChLoaderGravity that applies the load to elements 64 /// only if they are inherited by ChLoadableUVW so it can use GetDensity() and Gauss quadrature. 65 virtual void EleIntLoadResidual_F_gravity(ChVectorDynamic<>& R, const ChVector<>& G_acc, const double c) override; 66 67 68 // 69 // FEM functions 70 // 71 72 /// (This is the default implementation (POTENTIALLY INEFFICIENT) so that in children classes you can avoid 73 /// implementing this ComputeGravityForces() function, unless you need faster code.) 74 /// This fallback implementation uses a temp ChLoaderGravity that applies the load to elements 75 /// only if they are inherited by ChLoadableUVW so it can use GetDensity() and Gauss quadrature. 76 virtual void ComputeGravityForces(ChVectorDynamic<>& Fg, const ChVector<>& G_acc) override; 77 78 /// Returns the global mass matrix. 79 /// This is the default implementation, POTENTIALLY VERY INEFFICIENT. 80 /// Children classes may need to override this with a more efficient version. ComputeMmatrixGlobal(ChMatrixRef M)81 virtual void ComputeMmatrixGlobal(ChMatrixRef M) override { ComputeKRMmatricesGlobal(M, 0, 0, 1.0); } 82 83 // 84 // Functions for interfacing to the solver 85 // 86 87 /// Tell to a system descriptor that there are item(s) of type 88 /// ChKblock in this object (for further passing it to a solver) InjectKRMmatrices(ChSystemDescriptor & mdescriptor)89 virtual void InjectKRMmatrices(ChSystemDescriptor& mdescriptor) override { mdescriptor.InsertKblock(&Kmatr); } 90 91 /// Adds the current stiffness K and damping R and mass M matrices in encapsulated 92 /// ChKblock item(s), if any. The K, R, M matrices are load with scaling 93 /// values Kfactor, Rfactor, Mfactor. KRMmatricesLoad(double Kfactor,double Rfactor,double Mfactor)94 virtual void KRMmatricesLoad(double Kfactor, double Rfactor, double Mfactor) override { 95 this->ComputeKRMmatricesGlobal(this->Kmatr.Get_K(), Kfactor, Rfactor, Mfactor); 96 } 97 98 /// Adds the internal forces, expressed as nodal forces, into the 99 /// encapsulated ChVariables, in the 'fb' part: qf+=forces*factor 100 /// (This is a default (a bit unoptimal) book keeping so that in children classes you can avoid 101 /// implementing this VariablesFbLoadInternalForces function, unless you need faster code) 102 virtual void VariablesFbLoadInternalForces(double factor = 1.) override; 103 104 /// Adds M*q (internal masses multiplied current 'qb') to Fb, ex. if qb is initialized 105 /// with v_old using VariablesQbLoadSpeed, this method can be used in 106 /// timestepping schemes that do: M*v_new = M*v_old + forces*dt 107 /// (This is a default (VERY UNOPTIMAL) book keeping so that in children classes you can avoid 108 /// implementing this VariablesFbIncrementMq function, unless you need faster code.) 109 virtual void VariablesFbIncrementMq() override; 110 }; 111 112 /// @} fea_elements 113 114 } // end namespace fea 115 } // end namespace chrono 116 117 #endif 118