| 117 | } |
| 118 | |
| 119 | ZSparseMatrix Assembler3D::getMassMatrix(bool lumped, int repeats) { |
| 120 | double p = m_density; |
| 121 | |
| 122 | typedef Eigen::Triplet<double> T; |
| 123 | std::vector<T> triplets; |
| 124 | |
| 125 | if (lumped) { |
| 126 | for (size_t i=0; i<m_mesh->getNbrElements(); i++) { |
| 127 | VectorI idx = m_mesh->getElement(i); |
| 128 | assert(idx.size() == 4); |
| 129 | double V = m_mesh->getElementVolume(i); |
| 130 | |
| 131 | for (size_t j=0; j<idx.size(); j++) |
| 132 | for (size_t l=0; l<repeats; l++) |
| 133 | triplets.push_back(T(repeats*idx[j]+l, |
| 134 | repeats*idx[j]+l, p*V/4.0)); |
| 135 | } |
| 136 | } else { |
| 137 | double coeff_jj = 0.1, |
| 138 | coeff_jk = 0.05; |
| 139 | for (size_t i=0; i<m_mesh->getNbrElements(); ++i) { |
| 140 | VectorI idx = m_mesh->getElement(i); |
| 141 | assert(idx.size() == 4); |
| 142 | double V = m_mesh->getElementVolume(i); |
| 143 | |
| 144 | for (size_t j=0; j<4; ++j) { |
| 145 | for (size_t k=0; k<4; ++k) { |
| 146 | if (idx[j] == idx[k]) { |
| 147 | for (size_t l=0; l<repeats; ++l) |
| 148 | triplets.push_back( |
| 149 | T(repeats*idx[j]+l, |
| 150 | repeats*idx[k]+l, p*V*coeff_jj)); |
| 151 | } else { |
| 152 | for (size_t l=0; l<repeats; ++l) |
| 153 | triplets.push_back( |
| 154 | T(repeats*idx[j]+l, |
| 155 | repeats*idx[k]+l, p*V*coeff_jk)); |
| 156 | } |
| 157 | } |
| 158 | } |
| 159 | } |
| 160 | } |
| 161 | |
| 162 | Eigen::SparseMatrix<double> M = Eigen::SparseMatrix<double>( |
| 163 | repeats*m_mesh->getNbrNodes(), repeats*m_mesh->getNbrNodes()); |
| 164 | M.setFromTriplets(triplets.begin(), triplets.end()); |
| 165 | return ZSparseMatrix(M); |
| 166 | } |
| 167 | |
| 168 | ZSparseMatrix Assembler3D::getStiffnessMatrix() { |
| 169 | // Elastic modulii |
nothing calls this directly
no test coverage detected