Method getRigidMotionMatrix

tools/Assembler/ref/Assembler.cpp:169–210

}

ZSparseMatrix Assembler::getRigidMotionMatrix() {
    typedef Eigen::Triplet<double> T;
    std::vector<T> triplets;
    size_t dim = m_mesh->getDim();

    // Translation part.
    for (size_t i=0; i<m_mesh->getNbrNodes(); i++) {
        for (size_t j=0; j<dim; j++) {
            triplets.push_back(T(j,i*dim+j, 1.0));
        }
    }

    // Compute centroid
    VectorF c(dim);
    c *= 0;

    // Rotation part.
    size_t rot_degrees = 0;
    if (dim == 2) {
        rot_degrees = 1;
        for (size_t i=0; i<m_mesh->getNbrNodes(); i++) {
            VectorF x = m_mesh->getNode(i);
            triplets.push_back(T(dim  , i*dim  , -x[1]+c[1]));
            triplets.push_back(T(dim  , i*dim+1,  x[0]-c[0]));
        }
    } else {
        rot_degrees = 3;
        for (size_t i=0; i<m_mesh->getNbrNodes(); i++) {
            VectorF x = m_mesh->getNode(i);
            triplets.push_back(T(dim  , i*dim+1, -x[2]+c[2]));
            triplets.push_back(T(dim  , i*dim+2,  x[1]-c[1]));
            triplets.push_back(T(dim+1, i*dim  ,  x[2]-c[2]));
            triplets.push_back(T(dim+1, i*dim+2, -x[0]+c[0]));
            triplets.push_back(T(dim+2, i*dim  , -x[1]+c[1]));
            triplets.push_back(T(dim+2, i*dim+1,  x[0]-c[0]));
        }
    }

    Eigen::SparseMatrix<double> Ru = Eigen::SparseMatrix<double>(dim + rot_degrees, dim * m_mesh->getNbrNodes());
    Ru.setFromTriplets(triplets.begin(), triplets.end());
    return ZSparseMatrix(Ru);
}