preCICE v3.1.2
Loading...
Searching...
No Matches
EigenHelperFunctions.hpp
Go to the documentation of this file.
1#pragma once
2
3#include <Eigen/Core>
4#include <iterator>
5#include <vector>
6#include "utils/algorithm.hpp"
7#include "utils/assertion.hpp"
8
9namespace precice {
10namespace utils {
11
12void shiftSetFirst(Eigen::MatrixXd &A, const Eigen::VectorXd &v);
13
14void appendFront(Eigen::MatrixXd &A, Eigen::VectorXd &v);
15
16void removeColumnFromMatrix(Eigen::MatrixXd &A, int col);
17
18void append(Eigen::VectorXd &v, double value);
19
20template <typename Derived1>
21void append(
22 Eigen::MatrixXd & A,
23 const Eigen::PlainObjectBase<Derived1> &B)
24{
25 int n = A.rows(), m = A.cols();
26 if (n <= 0 && m <= 0) {
27 A = B;
28 } else {
29 PRECICE_ASSERT(B.rows() == n, B.rows(), A.rows());
30 A.conservativeResize(n, m + B.cols());
31 for (auto i = m; i < m + B.cols(); i++)
32 A.col(i) = B.col(i - m);
33 }
34}
35
36template <typename Derived1>
37void append(
38 Eigen::VectorXd & v,
39 const Eigen::PlainObjectBase<Derived1> &app)
40{
41 int n = v.size();
42 if (n <= 0) {
43 v = app;
44 } else {
45 PRECICE_ASSERT(app.cols() == 1, app.cols());
46 v.conservativeResize(n + app.size());
47 for (int i = 0; i < app.size(); i++)
48 v(n + i) = app(i);
49 }
50}
51
57template <typename Derived>
58auto firstN(const Eigen::PlainObjectBase<Derived> &val, unsigned n) -> const Eigen::Map<const Eigen::Matrix<typename Derived::Scalar, 1, Eigen::Dynamic>>
59{
60 return {val.data(), std::min<Eigen::Index>(n, val.size())};
61}
62
63template <typename Derived, typename Iter = const typename Derived::Scalar *, typename Size = typename std::iterator_traits<Iter>::difference_type>
64const RangePreview<Iter> previewRange(Size n, const Eigen::PlainObjectBase<Derived> &eigen)
65{
66 auto begin = eigen.data();
67 auto end = begin;
68 std::advance(end, eigen.size());
69 return {n, begin, end};
70}
71
72template <typename DerivedLHS, typename DerivedRHS>
73bool componentWiseLess(const Eigen::PlainObjectBase<DerivedLHS> &lhs, const Eigen::PlainObjectBase<DerivedRHS> &rhs)
74{
75 if (lhs.size() != rhs.size())
76 return false;
77
78 for (int i = 0; i < rhs.size(); ++i) {
79 if (lhs[i] != rhs[i]) {
80 return lhs[i] < rhs[i];
81 }
82 }
83 return false;
84}
85
87 template <typename DerivedLHS, typename DerivedRHS>
88 bool operator()(const Eigen::PlainObjectBase<DerivedLHS> &lhs, const Eigen::PlainObjectBase<DerivedRHS> &rhs) const
89 {
90 return componentWiseLess(lhs, rhs);
91 }
92};
93
94} // namespace utils
95} // namespace precice
T advance(T... args)
#define PRECICE_ASSERT(...)
Definition assertion.hpp:87
void removeColumnFromMatrix(Eigen::MatrixXd &A, int col)
void appendFront(Eigen::MatrixXd &A, Eigen::VectorXd &v)
void append(Eigen::VectorXd &v, double value)
const RangePreview< Iter > previewRange(Size n, const Range &range)
auto firstN(const Eigen::PlainObjectBase< Derived > &val, unsigned n) -> const Eigen::Map< const Eigen::Matrix< typename Derived::Scalar, 1, Eigen::Dynamic > >
void shiftSetFirst(Eigen::MatrixXd &A, const Eigen::VectorXd &v)
bool componentWiseLess(const Eigen::PlainObjectBase< DerivedLHS > &lhs, const Eigen::PlainObjectBase< DerivedRHS > &rhs)
Main namespace of the precice library.
bool operator()(const Eigen::PlainObjectBase< DerivedLHS > &lhs, const Eigen::PlainObjectBase< DerivedRHS > &rhs) const
The RangePreview object used as a lazy proxy struct for proviewing the content of a Range.