Class Elasticity_Kugelstadt2021

Inheritance Relationships

Base Type

Class Documentation

class Elasticity_Kugelstadt2021 : public SPH::NonPressureForceBase

This class implements the implicit SPH formulation for incompressible linearly elastic solids introduced by Kugelstadt et al. [KBF+21].

References:

  • [KBF+21] Tassilo Kugelstadt, Jan Bender, José Antonio Fernández-Fernández, Stefan Rhys Jeske, Fabian Löschner, Andreas Longva. Fast Corotated Elastic SPH Solids with Implicit Zero-Energy Mode Control. Proceedings of the ACM on Computer Graphics and Interactive Techniques, 2021. URL: http://dx.doi.org/10.1145/3480142

Public Functions

Elasticity_Kugelstadt2021(FluidModel *model)
virtual ~Elasticity_Kugelstadt2021(void)
inline virtual std::string getMethodName()

returns the name of the method

virtual void step()

Perform a step of the elasticity solver.

virtual void reset()
virtual void performNeighborhoodSearchSort()
virtual void saveState(BinaryFileWriter &binWriter)
virtual void loadState(BinaryFileReader &binReader)
void computeRotations()

Extract rotation matrices from deformation gradients.

Public Static Functions

static inline NonPressureForceBase *creator(FluidModel *model)
static void matrixVecProd(const Real *vec, Real *result, void *userData)

Matrix vector product used by the matrix-free conjugate gradient solver to solve the system in Eq. 30.

Public Static Attributes

static std::string METHOD_NAME = "Kugelstadt et al. 2021"
static int YOUNGS_MODULUS = -1
static int POISSON_RATIO = -1
static int FIXED_BOX_MIN = -1
static int FIXED_BOX_MAX = -1
static int ITERATIONS_V = -1
static int MAX_ITERATIONS_V = -1
static int MAX_ERROR_V = -1
static int ALPHA = -1
static int MAX_NEIGHBORS = -1

Protected Types

typedef Eigen::SimplicialLLT<Eigen::SparseMatrix<double>, Eigen::Lower, Eigen::AMDOrdering<int>> SolverLLT
typedef Eigen::ConjugateGradient<MatrixReplacement, Eigen::Lower | Eigen::Upper, Eigen::IdentityPreconditioner> Solver

Protected Functions

void computeRHS(VectorXr &rhs)

Compute right hand side of the linear system of the volume solver (Eq. 30).

void determineFixedParticles()

Mark all particles in the bounding box as fixed.

std::string computeMD5(const unsigned int objIndex)

Compute an MD4 check sum using the neighborhood structure in order to recognize known particle models (cache).

void initValues()

Initialize the particle neighborhoods in the reference configuration. Fix particles which lie in the user-defined bounding box. Find out if there are multiple separate objects in the phase. Finally, compute kernel gradient correction matrices and factorization.

void initSystem()

Initialize the solver for the linear system by either computing a factorization or loading a factorization from the cache.

void initFactorization(std::shared_ptr<Factorization> factorization, std::vector<unsigned int> &particleIndices, const unsigned int nFixed, const Real dt, const Real mu)

Compute the factorization of the linear system matrix. This is only done once at the beginning of the simulation.

void findObjects()

Find separate objects by object id.

void computeMatrixL()

Compute kernel gradient correction matrices (Eq. 8).

void precomputeValues()

Precompute some values and products to improve the performance of the solvers.

void stepElasticitySolver()

Solve the linear system for the stretching forces including zero energy mode control using the precomputed matrix factorization.

void stepVolumeSolver()

Solver for the volume conservation forces (Eq. 30).

virtual void initParameters()
virtual void deferredInit()

This function is called after the simulation scene is loaded and all parameters are initialized. While reading a scene file several parameters can change. The deferred init function should initialize all values which depend on these parameters.

inline FORCE_INLINE void symMatTimesVec (const Vector6r &M, const Vector3r &v, Vector3r &res)
void rotationMatricesToAVXQuaternions()

Protected Attributes

Real m_youngsModulus
Real m_poissonRatio
Vector3r m_fixedBoxMin
Vector3r m_fixedBoxMax
std::vector<unsigned int> m_current_to_initial_index
std::vector<unsigned int> m_initial_to_current_index
std::vector<std::vector<unsigned int>> m_initialNeighbors
std::vector<Real> m_restVolumes
std::vector<Matrix3r> m_rotations
std::vector<Real> m_stress
std::vector<Matrix3r> m_L
std::vector<Matrix3r> m_F
std::vector<Vector3r> m_vDiff
std::vector<Matrix3r> m_RL
unsigned int m_iterationsV
unsigned int m_maxIterV
Real m_maxErrorV
Real m_alpha
int m_maxNeighbors
unsigned int m_totalNeighbors
std::vector<ElasticObject*> m_objects
Real m_lambda
Real m_mu
std::vector<Vector3r, Eigen::aligned_allocator<Vector3r>> m_precomp_RL_gradW
std::vector<Vector3r, Eigen::aligned_allocator<Vector3r>> m_precomp_L_gradW
std::vector<Vector3r, Eigen::aligned_allocator<Vector3r>> m_precomp_RLj_gradW
std::vector<unsigned int> m_precomputed_indices
Solver m_solver
struct ElasticObject

Public Functions

inline ElasticObject()
inline ~ElasticObject()

Public Members

std::string m_md5
std::vector<unsigned int> m_particleIndices
unsigned int m_nFixed
std::shared_ptr<Factorization> m_factorization
std::vector<Vector3r, Eigen::aligned_allocator<Vector3r>> m_f
std::vector<Vector3r, Eigen::aligned_allocator<Vector3r>> m_sol
std::vector<Vector3r, Eigen::aligned_allocator<Vector3r>> m_RHS
std::vector<Vector3r, Eigen::aligned_allocator<Vector3r>> m_RHS_perm
std::vector<Quaternionr, Eigen::aligned_allocator<Quaternionr>> m_quats
struct Factorization

Public Members

Real m_dt
Real m_mu
Eigen::SparseMatrix<Real, Eigen::RowMajor> m_DT_K
Eigen::SparseMatrix<Real, Eigen::RowMajor> m_D
Eigen::SparseMatrix<Real, Eigen::ColMajor> m_matHTH
Eigen::SparseMatrix<Real, Eigen::ColMajor> m_matL
Eigen::SparseMatrix<Real, Eigen::ColMajor> m_matLT
Eigen::VectorXi m_permInd
Eigen::VectorXi m_permInvInd