PhysicsBasedAnimationToolkit 0.0.10
Cross-platform C++20 library of algorithms and data structures commonly used in computer graphics research on physically-based simulation.
Loading...
Searching...
No Matches
Mass.h
Go to the documentation of this file.
1
14
15#ifndef PBAT_FEM_MASS_H
16#define PBAT_FEM_MASS_H
17
18#include "Concepts.h"
19#include "ShapeFunctions.h"
20#include "pbat/Aliases.h"
21#include "pbat/common/Eigen.h"
23
24#include <exception>
25#include <fmt/core.h>
26#include <tbb/parallel_for.h>
27
28namespace pbat::fem {
29
41template <typename TElement, typename TN>
42inline auto
43ElementMassMatrix(Eigen::MatrixBase<TN> const& N, typename TN::Scalar w, typename TN::Scalar rho)
44{
45 static_assert(
46 TElement::kNodes == TN::RowsAtCompileTime || TN::RowsAtCompileTime == Eigen::Dynamic,
47 "Shape function vector size must match number of element nodes");
48 return w * rho * N * N.transpose();
49}
50
65template <
66 typename TElement,
67 typename TN,
68 typename TDerivedwg,
69 typename TDerivedrhog,
70 typename TDerivedOut>
72 Eigen::MatrixBase<TN> const& Neg,
73 Eigen::DenseBase<TDerivedwg> const& wg,
74 Eigen::DenseBase<TDerivedrhog> const& rhog,
75 Eigen::PlainObjectBase<TDerivedOut>& Meg)
76{
77 constexpr int kNodes = TElement::kNodes;
78 const auto nQuadPts = Neg.cols();
79 Meg.resize(kNodes, kNodes * nQuadPts);
80 for (Eigen::Index g = 0; g < nQuadPts; ++g)
81 {
82 auto const N = Neg.col(g);
83 auto const w = wg(g);
84 auto const rho = rhog(g);
85 auto Me = w * rho * N * N.transpose();
86 Meg.template block<kNodes, kNodes>(0, g * kNodes) = Me;
87 }
88}
89
102template <typename TElement, typename TN, typename TDerivedwg, typename TDerivedrhog>
104 Eigen::MatrixBase<TN> const& Neg,
105 Eigen::DenseBase<TDerivedwg> const& wg,
106 Eigen::DenseBase<TDerivedrhog> const& rhog)
107{
108 using ScalarType = typename TN::Scalar;
109 constexpr int kNodes = TElement::kNodes;
110 Eigen::Matrix<ScalarType, kNodes, Eigen::Dynamic> Meg;
111 ToElementMassMatrices<TElement>(Neg.derived(), wg.derived(), rhog.derived(), Meg);
112 return Meg;
113}
114
135template <
136 CElement TElement,
137 int Dims,
138 class TDerivedE,
139 class TDerivedeg,
140 class TDerivedMe,
141 class TDerivedIn,
142 class TDerivedOut>
143void GemmMass(
144 Eigen::DenseBase<TDerivedE> const& E,
145 Eigen::Index nNodes,
146 Eigen::DenseBase<TDerivedeg> const& eg,
147 Eigen::MatrixBase<TDerivedMe> const& Me,
148 int dims,
149 Eigen::MatrixBase<TDerivedIn> const& X,
150 Eigen::DenseBase<TDerivedOut>& Y);
151
168template <CMesh TMesh, class TDerivedeg, class TDerivedMe, class TDerivedIn, class TDerivedOut>
169inline void GemmMass(
170 TMesh const& mesh,
171 Eigen::DenseBase<TDerivedeg> const& eg,
172 Eigen::MatrixBase<TDerivedMe> const& Me,
173 int dims,
174 Eigen::MatrixBase<TDerivedIn> const& X,
175 Eigen::DenseBase<TDerivedOut>& Y)
176{
178 mesh.E,
179 static_cast<typename TMesh::IndexType>(mesh.X.cols()),
180 eg.derived(),
181 Me.derived(),
182 dims,
183 X.derived(),
184 Y.derived());
185}
186
206template <
207 CElement TElement,
208 int Dims,
209 Eigen::StorageOptions Options,
210 class TDerivedE,
211 class TDerivedeg,
212 class TDerivedwg,
213 class TDerivedrhog,
214 class TDerivedNeg>
215auto MassMatrix(
216 Eigen::DenseBase<TDerivedE> const& E,
217 Eigen::Index nNodes,
218 Eigen::DenseBase<TDerivedeg> const& eg,
219 Eigen::DenseBase<TDerivedwg> const& wg,
220 Eigen::DenseBase<TDerivedrhog> const& rhog,
221 Eigen::MatrixBase<TDerivedNeg> const& Neg,
222 int dims = 1)
223 -> Eigen::SparseMatrix<typename TDerivedNeg::Scalar, Options, typename TDerivedE::Scalar>;
224
242template <
243 CElement TElement,
244 int Dims,
245 Eigen::StorageOptions Options,
246 class TDerivedE,
247 class TDerivedeg,
248 class TDerivedMe>
249auto MassMatrix(
250 Eigen::DenseBase<TDerivedE> const& E,
251 Eigen::Index nNodes,
252 Eigen::DenseBase<TDerivedeg> const& eg,
253 Eigen::MatrixBase<TDerivedMe> const& Meg,
254 int dims = 1)
255 -> Eigen::SparseMatrix<typename TDerivedMe::Scalar, Options, typename TDerivedE::Scalar>;
256
273template <
274 Eigen::StorageOptions Options,
275 CMesh TMesh,
276 class TDerivedeg,
277 class TDerivedwg,
278 class TDerivedrhog,
279 class TDerivedNeg>
281 TMesh const& mesh,
282 Eigen::DenseBase<TDerivedeg> const& eg,
283 Eigen::DenseBase<TDerivedwg> const& wg,
284 Eigen::DenseBase<TDerivedrhog> const& rhog,
285 Eigen::MatrixBase<TDerivedNeg> const& Neg,
286 int dims = 1)
287 -> Eigen::SparseMatrix<typename TDerivedNeg::Scalar, Options, typename TMesh::IndexType>
288{
290 mesh.E,
291 static_cast<typename TMesh::IndexType>(mesh.X.cols()),
292 eg.derived(),
293 wg.derived(),
294 rhog.derived(),
295 Neg.derived(),
296 dims);
297}
298
313template <Eigen::StorageOptions Options, CMesh TMesh, class TDerivedeg, class TDerivedMe>
315 TMesh const& mesh,
316 Eigen::DenseBase<TDerivedeg> const& eg,
317 Eigen::MatrixBase<TDerivedMe> const& Meg,
318 int dims = 1)
319 -> Eigen::SparseMatrix<typename TDerivedMe::Scalar, Options, typename TMesh::IndexType>
320{
322 mesh.E,
323 static_cast<typename TMesh::IndexType>(mesh.X.cols()),
324 eg.derived(),
325 Meg.derived(),
326 dims);
327}
328
348template <
349 CElement TElement,
350 int Dims,
351 class TDerivedE,
352 class TDerivedeg,
353 class TDerivedwg,
354 class TDerivedrhog,
355 class TDerivedNeg,
356 class TDerivedOut>
358 Eigen::DenseBase<TDerivedE> const& E,
359 Eigen::Index nNodes,
360 Eigen::DenseBase<TDerivedeg> const& eg,
361 Eigen::DenseBase<TDerivedwg> const& wg,
362 Eigen::DenseBase<TDerivedrhog> const& rhog,
363 Eigen::MatrixBase<TDerivedNeg> const& Neg,
364 int dims,
365 Eigen::PlainObjectBase<TDerivedOut>& m);
366
384template <
385 CElement TElement,
386 int Dims,
387 class TDerivedE,
388 class TDerivedeg,
389 class TDerivedMe,
390 class TDerivedOut>
392 Eigen::DenseBase<TDerivedE> const& E,
393 Eigen::Index nNodes,
394 Eigen::DenseBase<TDerivedeg> const& eg,
395 Eigen::MatrixBase<TDerivedMe> const& Meg,
396 int dims,
397 Eigen::PlainObjectBase<TDerivedOut>& m);
398
416template <
417 CMesh TMesh,
418 class TDerivedeg,
419 class TDerivedwg,
420 class TDerivedrhog,
421 class TDerivedNeg,
422 class TDerivedOut>
424 TMesh const& mesh,
425 Eigen::DenseBase<TDerivedeg> const& eg,
426 Eigen::DenseBase<TDerivedwg> const& wg,
427 Eigen::DenseBase<TDerivedrhog> const& rhog,
428 Eigen::MatrixBase<TDerivedNeg> const& Neg,
429 int dims,
430 Eigen::PlainObjectBase<TDerivedOut>& m)
431{
433 mesh.E,
434 static_cast<typename TMesh::IndexType>(mesh.X.cols()),
435 eg.derived(),
436 wg.derived(),
437 rhog.derived(),
438 Neg.derived(),
439 dims,
440 m.derived());
441}
442
457template <CMesh TMesh, class TDerivedeg, class TDerivedMe, class TDerivedOut>
459 TMesh const& mesh,
460 Eigen::DenseBase<TDerivedeg> const& eg,
461 Eigen::MatrixBase<TDerivedMe> const& Meg,
462 int dims,
463 Eigen::PlainObjectBase<TDerivedOut>& m)
464{
466 mesh.E,
467 static_cast<typename TMesh::IndexType>(mesh.X.cols()),
468 eg.derived(),
469 Meg.derived(),
470 dims,
471 m.derived());
472}
473
492template <
493 CElement TElement,
494 int Dims,
495 class TDerivedE,
496 class TDerivedeg,
497 class TDerivedwg,
498 class TDerivedrhog,
499 class TDerivedNeg>
501 Eigen::DenseBase<TDerivedE> const& E,
502 Eigen::Index nNodes,
503 Eigen::DenseBase<TDerivedeg> const& eg,
504 Eigen::DenseBase<TDerivedwg> const& wg,
505 Eigen::DenseBase<TDerivedrhog> const& rhog,
506 Eigen::MatrixBase<TDerivedNeg> const& Neg,
507 int dims = 1) -> Eigen::Vector<typename TDerivedNeg::Scalar, Eigen::Dynamic>
508{
509 using ScalarType = typename TDerivedNeg::Scalar;
510 auto const numberOfDofs = dims * nNodes;
511 Eigen::Vector<ScalarType, Eigen::Dynamic> m(numberOfDofs);
512 ToLumpedMassMatrix<TElement, Dims>(E, nNodes, eg, wg, rhog, Neg, dims, m);
513 return m;
514}
515
532template <CElement TElement, int Dims, class TDerivedE, class TDerivedeg, class TDerivedMe>
534 Eigen::DenseBase<TDerivedE> const& E,
535 Eigen::Index nNodes,
536 Eigen::DenseBase<TDerivedeg> const& eg,
537 Eigen::MatrixBase<TDerivedMe> const& Meg,
538 int dims = 1) -> Eigen::Vector<typename TDerivedMe::Scalar, Eigen::Dynamic>
539{
540 using ScalarType = typename TDerivedMe::Scalar;
541 auto const numberOfDofs = dims * nNodes;
542 Eigen::Vector<ScalarType, Eigen::Dynamic> m(numberOfDofs);
543 ToLumpedMassMatrix<TElement, Dims>(E, nNodes, eg, Meg, dims, m);
544 return m;
545}
546
550template <CMesh TMesh, class TDerivedeg, class TDerivedwg, class TDerivedrhog, class TDerivedNeg>
552 TMesh const& mesh,
553 Eigen::DenseBase<TDerivedeg> const& eg,
554 Eigen::DenseBase<TDerivedwg> const& wg,
555 Eigen::DenseBase<TDerivedrhog> const& rhog,
556 Eigen::MatrixBase<TDerivedNeg> const& Neg,
557 int dims = 1)
558{
560 mesh.E,
561 static_cast<typename TMesh::IndexType>(mesh.X.cols()),
562 eg.derived(),
563 wg.derived(),
564 rhog.derived(),
565 Neg.derived(),
566 dims);
567}
568
573template <CMesh TMesh, class TDerivedeg, class TDerivedMe>
575 TMesh const& mesh,
576 Eigen::DenseBase<TDerivedeg> const& eg,
577 Eigen::MatrixBase<TDerivedMe> const& Meg,
578 int dims = 1)
579{
581 mesh.E,
582 static_cast<typename TMesh::IndexType>(mesh.X.cols()),
583 eg.derived(),
584 Meg.derived(),
585 dims);
586}
587
588// Implementation of GemmMass
589template <
590 CElement TElement,
591 int Dims,
592 class TDerivedE,
593 class TDerivedeg,
594 class TDerivedMe,
595 class TDerivedIn,
596 class TDerivedOut>
597inline void GemmMass(
598 Eigen::DenseBase<TDerivedE> const& E,
599 Eigen::Index nNodes,
600 Eigen::DenseBase<TDerivedeg> const& eg,
601 Eigen::MatrixBase<TDerivedMe> const& Meg,
602 int dims,
603 Eigen::MatrixBase<TDerivedIn> const& X,
604 Eigen::DenseBase<TDerivedOut>& Y)
605{
606 PBAT_PROFILE_NAMED_SCOPE("pbat.fem.GemmMass");
607
608 // Check inputs
609 auto const numberOfDofs = dims * nNodes;
610 bool const bDimensionsMatch =
611 (X.cols() == Y.cols()) and (X.rows() == numberOfDofs) and (Y.rows() == numberOfDofs);
612 if (not bDimensionsMatch)
613 {
614 std::string const what = fmt::format(
615 "Expected input and output to have {} rows and same number of columns, but got "
616 "dimensions X,Y=({} x {}), ({} x {})",
617 numberOfDofs,
618 X.rows(),
619 X.cols(),
620 Y.rows(),
621 Y.cols());
622 throw std::invalid_argument(what);
623 }
624
625 if (dims < 1)
626 {
627 std::string const what = fmt::format("Expected dims >= 1, got {} instead", dims);
628 throw std::invalid_argument(what);
629 }
630
631 auto constexpr kNodesPerElement = TElement::kNodes;
632 auto const nQuadPts = eg.size();
633
634 // Apply precomputed element mass matrices
635 for (auto c = 0; c < X.cols(); ++c)
636 {
637 for (auto g = 0; g < nQuadPts; ++g)
638 {
639 auto const e = eg(g);
640 auto const nodes = E.col(e);
641
642 // Get precomputed element mass matrix for this quadrature point
643 auto const Me =
644 Meg.template block<kNodesPerElement, kNodesPerElement>(0, g * kNodesPerElement);
645
646 // Apply to each dimension
647 auto ye = Y.col(c).reshaped(dims, nNodes)(Eigen::placeholders::all, nodes);
648 auto const xe = X.col(c).reshaped(dims, nNodes)(Eigen::placeholders::all, nodes);
649 ye += xe * Me.transpose(); // Mass matrix is symmetric, so transpose doesn't matter
650 }
651 }
652}
653
654// Implementation of MassMatrix
655template <
656 CElement TElement,
657 int Dims,
658 Eigen::StorageOptions Options,
659 class TDerivedE,
660 class TDerivedeg,
661 class TDerivedwg,
662 class TDerivedrhog,
663 class TDerivedNeg>
665 Eigen::DenseBase<TDerivedE> const& E,
666 Eigen::Index nNodes,
667 Eigen::DenseBase<TDerivedeg> const& eg,
668 Eigen::DenseBase<TDerivedwg> const& wg,
669 Eigen::DenseBase<TDerivedrhog> const& rhog,
670 Eigen::MatrixBase<TDerivedNeg> const& Neg,
671 int dims)
672 -> Eigen::SparseMatrix<typename TDerivedNeg::Scalar, Options, typename TDerivedE::Scalar>
673{
674 PBAT_PROFILE_NAMED_SCOPE("pbat.fem.MassMatrix");
675 using ScalarType = typename TDerivedNeg::Scalar;
676 using IndexType = typename TDerivedE::Scalar;
677 using SparseMatrixType = Eigen::SparseMatrix<ScalarType, Options, IndexType>;
678 using Triplet = Eigen::Triplet<ScalarType, IndexType>;
679
680 auto constexpr kNodesPerElement = TElement::kNodes;
681 auto const nQuadPts = eg.size();
682
683 std::vector<Triplet> triplets{};
684 triplets.reserve(
685 static_cast<std::size_t>(kNodesPerElement * kNodesPerElement * nQuadPts * dims));
686
687 auto const numberOfDofs = dims * nNodes;
688 SparseMatrixType M(numberOfDofs, numberOfDofs);
689
690 for (auto g = 0; g < nQuadPts; ++g)
691 {
692 auto const e = eg(g);
693 auto const nodes = E.col(e);
694 auto const w = wg(g); // Includes Jacobian determinant
695 auto const rho = rhog(g);
696
697 // Get shape functions at this quadrature point
698 auto const N = Neg.col(g);
699
700 // Compute element mass matrix: w * rho * N * N^T
701 auto const Me = w * rho * N * N.transpose();
702
703 // Add contributions for each dimension
704 for (auto i = 0; i < kNodesPerElement; ++i)
705 {
706 for (auto j = 0; j < kNodesPerElement; ++j)
707 {
708 for (auto d = 0; d < dims; ++d)
709 {
710 auto const ni = static_cast<IndexType>(dims * nodes(i) + d);
711 auto const nj = static_cast<IndexType>(dims * nodes(j) + d);
712 triplets.emplace_back(ni, nj, Me(i, j));
713 }
714 }
715 }
716 }
717
718 M.setFromTriplets(triplets.begin(), triplets.end());
719 return M;
720}
721
722template <
723 CElement TElement,
724 int Dims,
725 Eigen::StorageOptions Options,
726 class TDerivedE,
727 class TDerivedeg,
728 class TDerivedMe>
730 Eigen::DenseBase<TDerivedE> const& E,
731 Eigen::Index nNodes,
732 Eigen::DenseBase<TDerivedeg> const& eg,
733 Eigen::MatrixBase<TDerivedMe> const& Meg,
734 int dims)
735 -> Eigen::SparseMatrix<typename TDerivedMe::Scalar, Options, typename TDerivedE::Scalar>
736{
737 PBAT_PROFILE_NAMED_SCOPE("pbat.fem.MassMatrix");
738 using ScalarType = typename TDerivedMe::Scalar;
739 using IndexType = typename TDerivedE::Scalar;
740 using SparseMatrixType = Eigen::SparseMatrix<ScalarType, Options, IndexType>;
741 using Triplet = Eigen::Triplet<ScalarType, IndexType>;
742
743 auto constexpr kNodesPerElement = TElement::kNodes;
744 auto const nQuadPts = eg.size();
745
746 std::vector<Triplet> triplets{};
747 triplets.reserve(
748 static_cast<std::size_t>(kNodesPerElement * kNodesPerElement * nQuadPts * dims));
749
750 auto const numberOfDofs = dims * nNodes;
751 SparseMatrixType M(numberOfDofs, numberOfDofs);
752
753 for (auto g = 0; g < nQuadPts; ++g)
754 {
755 auto const e = eg(g);
756 auto const nodes = E.col(e);
757
758 // Get precomputed element mass matrix for this quadrature point
759 auto const Me =
760 Meg.template block<kNodesPerElement, kNodesPerElement>(0, g * kNodesPerElement);
761
762 // Add contributions for each dimension
763 for (auto i = 0; i < kNodesPerElement; ++i)
764 {
765 for (auto j = 0; j < kNodesPerElement; ++j)
766 {
767 for (auto d = 0; d < dims; ++d)
768 {
769 auto const ni = static_cast<IndexType>(dims * nodes(i) + d);
770 auto const nj = static_cast<IndexType>(dims * nodes(j) + d);
771 triplets.emplace_back(ni, nj, Me(i, j));
772 }
773 }
774 }
775 }
776
777 M.setFromTriplets(triplets.begin(), triplets.end());
778 return M;
779}
780
781// Implementation of ToLumpedMassMatrix
782template <
783 CElement TElement,
784 int Dims,
785 class TDerivedE,
786 class TDerivedeg,
787 class TDerivedwg,
788 class TDerivedrhog,
789 class TDerivedNeg,
790 class TDerivedOut>
792 Eigen::DenseBase<TDerivedE> const& E,
793 Eigen::Index nNodes,
794 Eigen::DenseBase<TDerivedeg> const& eg,
795 Eigen::DenseBase<TDerivedwg> const& wg,
796 Eigen::DenseBase<TDerivedrhog> const& rhog,
797 Eigen::MatrixBase<TDerivedNeg> const& Neg,
798 int dims,
799 Eigen::PlainObjectBase<TDerivedOut>& m)
800{
801 auto const numberOfDofs = dims * nNodes;
802 m.resize(numberOfDofs, 1);
803 m.setZero();
804 auto constexpr kNodesPerElement = TElement::kNodes;
805 auto const nQuadPts = eg.size();
806 for (auto g = 0; g < nQuadPts; ++g)
807 {
808 auto const e = eg(g);
809 auto const nodes = E.col(e);
810 auto const w = wg(g); // Includes Jacobian determinant
811 auto const rho = rhog(g);
812
813 // Get shape functions at this quadrature point
814 auto const N = Neg.col(g);
815
816 // Compute element mass matrix: w * rho * N * N^T
817 auto const Me = (w * rho * N * N.transpose()).eval();
818
819 // Add lumped contributions (sum rows)
820 for (auto i = 0; i < kNodesPerElement; ++i)
821 {
822 for (auto j = 0; j < kNodesPerElement; ++j)
823 {
824 for (auto d = 0; d < dims; ++d)
825 {
826 auto const ni = dims * nodes(i) + d;
827 m(ni) += Me(i, j);
828 }
829 }
830 }
831 }
832}
833
834template <
835 CElement TElement,
836 int Dims,
837 class TDerivedE,
838 class TDerivedeg,
839 class TDerivedMe,
840 class TDerivedOut>
842 Eigen::DenseBase<TDerivedE> const& E,
843 Eigen::Index nNodes,
844 Eigen::DenseBase<TDerivedeg> const& eg,
845 Eigen::MatrixBase<TDerivedMe> const& Meg,
846 int dims,
847 Eigen::PlainObjectBase<TDerivedOut>& m)
848{
849 auto const numberOfDofs = dims * nNodes;
850 m.resize(numberOfDofs, 1);
851 m.setZero();
852 auto constexpr kNodesPerElement = TElement::kNodes;
853 auto const nQuadPts = eg.size();
854 for (auto g = 0; g < nQuadPts; ++g)
855 {
856 auto const e = eg(g);
857 auto const nodes = E.col(e);
858
859 // Get precomputed element mass matrix for this quadrature point
860 auto const Me =
861 Meg.template block<kNodesPerElement, kNodesPerElement>(0, g * kNodesPerElement);
862
863 // Add lumped contributions (sum rows)
864 for (auto i = 0; i < kNodesPerElement; ++i)
865 {
866 for (auto j = 0; j < kNodesPerElement; ++j)
867 {
868 for (auto d = 0; d < dims; ++d)
869 {
870 auto const ni = dims * nodes(i) + d;
871 m(ni) += Me(i, j);
872 }
873 }
874 }
875 }
876}
877
878} // namespace pbat::fem
879
880#endif // PBAT_FEM_MASS_H
FEM shape functions and gradients.
Eigen adaptors for ranges.
Finite Element Method (FEM)
Definition Concepts.h:19
auto ElementMassMatrices(Eigen::MatrixBase< TN > const &Neg, Eigen::DenseBase< TDerivedwg > const &wg, Eigen::DenseBase< TDerivedrhog > const &rhog)
Compute a matrix of horizontally stacked element mass matrices for all quadrature points.
Definition Mass.h:103
void ToElementMassMatrices(Eigen::MatrixBase< TN > const &Neg, Eigen::DenseBase< TDerivedwg > const &wg, Eigen::DenseBase< TDerivedrhog > const &rhog, Eigen::PlainObjectBase< TDerivedOut > &Meg)
Compute a matrix of horizontally stacked element mass matrices for all quadrature points.
Definition Mass.h:71
auto MassMatrix(Eigen::DenseBase< TDerivedE > const &E, Eigen::Index nNodes, Eigen::DenseBase< TDerivedeg > const &eg, Eigen::DenseBase< TDerivedwg > const &wg, Eigen::DenseBase< TDerivedrhog > const &rhog, Eigen::MatrixBase< TDerivedNeg > const &Neg, int dims=1) -> Eigen::SparseMatrix< typename TDerivedNeg::Scalar, Options, typename TDerivedE::Scalar >
Construct the mass matrix operator's sparse matrix representation.
Definition Mass.h:664
auto LumpedMassMatrix(Eigen::DenseBase< TDerivedE > const &E, Eigen::Index nNodes, Eigen::DenseBase< TDerivedeg > const &eg, Eigen::DenseBase< TDerivedwg > const &wg, Eigen::DenseBase< TDerivedrhog > const &rhog, Eigen::MatrixBase< TDerivedNeg > const &Neg, int dims=1) -> Eigen::Vector< typename TDerivedNeg::Scalar, Eigen::Dynamic >
Compute lumped mass matrix's diagonal vector (allocates output vector)
Definition Mass.h:500
auto ElementMassMatrix(Eigen::MatrixBase< TN > const &N, typename TN::Scalar w, typename TN::Scalar rho)
Compute the element mass matrix at a single quadrature point.
Definition Mass.h:43
void GemmMass(Eigen::DenseBase< TDerivedE > const &E, Eigen::Index nNodes, Eigen::DenseBase< TDerivedeg > const &eg, Eigen::MatrixBase< TDerivedMe > const &Me, int dims, Eigen::MatrixBase< TDerivedIn > const &X, Eigen::DenseBase< TDerivedOut > &Y)
Compute mass matrix-matrix multiply using precomputed element mass matrices.
Definition Mass.h:597
void ToLumpedMassMatrix(Eigen::DenseBase< TDerivedE > const &E, Eigen::Index nNodes, Eigen::DenseBase< TDerivedeg > const &eg, Eigen::DenseBase< TDerivedwg > const &wg, Eigen::DenseBase< TDerivedrhog > const &rhog, Eigen::MatrixBase< TDerivedNeg > const &Neg, int dims, Eigen::PlainObjectBase< TDerivedOut > &m)
Compute lumped mass matrix's diagonal vector into existing output vector.
Definition Mass.h:791
Profiling utilities for the Physics-Based Animation Toolkit (PBAT)