SpecialCosseratRodModel.h

Go to the documentation of this file.

 
#pragma once
 
#include <jem/base/Object.h>
 
#include <materials/Material.h>
#include <materials/MaterialFactory.h>
 
#include <jem/base/Array.h>
#include <jem/base/IllegalInputException.h>
#include <jem/base/System.h>
#include <jem/base/Float.h>
#include <jem/numeric/Quaternion.h>
#include <jem/numeric/algebra.h>
#include <jem/numeric/algebra/matmul.h>
#include <jem/util/Event.h>
#include <jem/util/Properties.h>
#include <jem/util/StringUtils.h>
 
#include <jive/Array.h>
#include <jive/algebra/AbstractMatrix.h>
#include <jive/algebra/FlexMatrixBuilder.h>
#include <jive/algebra/MatrixBuilder.h>
#include <jive/app/Names.h>
#include <jive/fem/ElementGroup.h>
#include <jive/fem/ElementSet.h>
#include <jive/fem/NodeGroup.h>
#include <jive/fem/NodeSet.h>
#include <jive/geom/ShapeFactory.h>
#include <jive/geom/StdShape.h>
#include <jive/implict/SolverInfo.h>
#include <jive/model/Actions.h>
#include <jive/model/Model.h>
#include <jive/model/ModelFactory.h>
#include <jive/model/StateVector.h>
#include <jive/util/Assignable.h>
#include <jive/util/Globdat.h>
#include <jive/util/Printer.h>
#include <jive/util/XDofSpace.h>
#include <jive/util/XTable.h>
#include <jive/util/utilities.h>
 
#include <math.h>
 
#include "misc/Line3D.h"
#include "utils/helpers.h"
#include "utils/testing.h"
 
using jem::Slice;
using jem::SliceFrom;
using jem::numeric::matmul;
using jem::numeric::MatmulChain;
using jem::numeric::norm2;
using jem::numeric::Quaternion;
using jem::util::Properties;
 
using jive::BoolVector;
using jive::IdxMatrix;
using jive::algebra::AbstractMatrix;
using jive::algebra::FlexMBuilder;
using jive::algebra::MatrixBuilder;
using jive::fem::ElementGroup;
using jive::fem::ElementSet;
using jive::fem::newNodeGroup;
using jive::fem::NodeGroup;
using jive::fem::NodeSet;
using jive::geom::IShape;
using jive::geom::ShapeFactory;
using jive::model::Model;
using jive::util::Assignable;
using jive::util::DofSpace;
using jive::util::Globdat;
using jive::util::joinNames;
using jive::util::XDofSpace;
using jive::util::XTable;
 
using jive_helpers::e3;
using jive_helpers::expVec;
using jive_helpers::eye;
using jive_helpers::skew;
using jive_helpers::unskew;
using jive_helpers::vec2mat;
 
//-----------------------------------------------------------------------
//   class SpecialCosseratRodModel
//-----------------------------------------------------------------------
 
class SpecialCosseratRodModel : public Model
{
 public:
  static const char *TYPE_NAME;         
  static const char *TRANS_DOF_DEFAULT; 
  static const char *ROT_DOF_DEFAULT;   
  static const char *TRANS_DOF_NAMES;   
  static const char *ROT_DOF_NAMES;     
  static const char *SYMMETRIC_ONLY;    
  static const char *MATERIAL_Y_DIR;    
  static const char *GIVEN_NODES;       
  static const char *GIVEN_DIRS;        
  static const char *LUMPED_MASS;       
  static const char *HINGES;            
 
  static const idx_t TRANS_DOF_COUNT; 
  static const idx_t ROT_DOF_COUNT;   
  static const Slice TRANS_PART;      
  static const Slice ROT_PART;        
 
  JEM_DECLARE_CLASS(SpecialCosseratRodModel, Model);
 
  explicit SpecialCosseratRodModel(const String &name,
                                   const Properties &conf,
                                   const Properties &props,
                                   const Properties &globdat);
 
  virtual bool takeAction(const String &action,
                          const Properties &params,
                          const Properties &globdat) override;
 
  static Ref<Model> makeNew(const String &name,
                            const Properties &conf,
                            const Properties &props,
                            const Properties &globdat);
 
  static void declare();
 
 private:
  void assemble_(MatrixBuilder &mbld,
                 const Vector &fint,
                 const Vector &disp,
                 const String &loadCase = "") const;
 
  void assemble_(const Vector &fint,
                 const Vector &disp,
                 const String &loadCase = "") const;
 
  void assembleGyro_(const Vector &fint,
                     const Vector &velo,
                     const Ref<AbstractMatrix> mass) const;
 
  void assembleM_(MatrixBuilder &mbld,
                  Vector &disp) const;
 
  void getStrainTable_(XTable &strain_table,
                       const Vector &weights,
                       const Vector &disp,
                       const bool mat_vals = false);
 
  void getMaterialTable_(XTable &mat_table,
                         const Vector &weights,
                         const String &name);
 
  void getStressTable_(XTable &stress_table,
                       const Vector &weights,
                       const Vector &disp,
                       const bool mat_vals = false);
 
  void initRotation_();
 
  void initShape_();
 
  void initStrain_();
 
  void getTangents_(const Matrix &phiP,
                    const Matrix &nodePhi_0,
                    const Matrix &nodeU,
                    const idx_t ie) const;
 
  void getGeomtericStiffness_(const Cubix &B,
                              const Matrix &stresses,
                              const Matrix &nodePhi_0,
                              const Matrix &nodeU) const;
 
  void getStrains_(const Matrix &strains,
                   const Vector &w,
                   const Matrix &nodePhi_0,
                   const Matrix &nodeU,
                   const Cubix &nodeLambda,
                   const idx_t ie,
                   const bool spatial = true) const;
 
  void getStresses_(const Matrix &stresses,
                    const Vector &w,
                    const Matrix &nodePhi_0,
                    const Matrix &nodeU,
                    const Cubix &nodeLambda,
                    const idx_t ie,
                    const bool spatial = true,
                    const String &loadCase = "") const;
 
  void updNodeDofs_() const;
 
  void invalidateNodeDofs_();
 
  void getDisplacments_(const Matrix &nodePhi_0,
                        const Matrix &nodeU,
                        const Cubix &nodeLambda,
                        const Vector &disp,
                        const IdxVector &inodes) const;
 
  void getEnergies_(double &E_pot,
                    double &E_diss,
                    const Vector &disp) const;
 
  void getPotentialEnergy_(XTable &energy_table,
                           const Vector &table_weights,
                           const Vector &disp) const;
 
  void getDissipatedEnergy_(XTable &energy_table,
                            const Vector &table_weights,
                            const Vector &disp) const;
 
 private:
  Assignable<ElementGroup> rodElems_; 
  IdxVector rodNodes_;                
  IdxVector rodNodeIdx_;              
  Assignable<ElementSet> allElems_;   
  Assignable<NodeSet> allNodes_;      
 
  Ref<DofSpace> dofs_;     
  Ref<Line3D> shapeK_;     
  Ref<Line3D> shapeM_;     
  Ref<Material> material_; 
  Ref<Model> hinges_;      
 
  IdxVector transTypes_; 
  IdxVector rotTypes_;   
  IdxVector jtypes_;     
 
  mutable IdxMatrix nodeDofs_;   
  mutable bool nodeDofsValid_;   
 
#ifndef NDEBUG
  mutable bool energiesStored_ = false; 
#endif
 
  bool symOnly_;        
  Vector thickFact_;    
  Vector materialYDir_; 
 
  IdxVector givenNodes_; 
  Matrix givenDirs_;     
 
  Cubix LambdaN_;    
  Cubix matStrain0_; 
 
  // The reference geometry does not change during a run (as for matStrain0_)
  Cubix elemGrads_;    
  Matrix elemWeights_; 
 
  // Work arrays of getStrains_, getStresses_ and the assembly: created once in
  // initShape_ and overwritten in every call (see the note in utils/helpers.h).
  // Not thread safe. ipPhiP_ is filled by getTangents_ and consumed at once.
  mutable Matrix curPhi_;     
  mutable Matrix ipPhiP_;     
  mutable Cubix ipLambda_;    
  mutable Cubix ipLambdaP_;   
  mutable Matrix lambdaTP_;   
  mutable Matrix strainWork_; 
  mutable MatmulChain<double, 2> matmul_; 
};