20 #ifndef OPM_AUTODIFFMATRIX_HEADER_INCLUDED
21 #define OPM_AUTODIFFMATRIX_HEADER_INCLUDED
23 #include <opm/common/utility/platform_dependent/disable_warnings.h>
25 #include <Eigen/Eigen>
26 #include <Eigen/Sparse>
28 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
30 #include <opm/common/ErrorMacros.hpp>
91 explicit AutoDiffMatrix(
const Eigen::DiagonalMatrix<double, Eigen::Dynamic>& d)
95 diag_(d.diagonal().array().data(), d.diagonal().array().data() + d.
rows()),
98 assert(rows_ == cols_);
144 std::swap(type_, other.type_);
145 std::swap(rows_, other.rows_);
146 std::swap(cols_, other.cols_);
147 diag_.swap(other.diag_);
148 sparse_.swap(other.sparse_);
161 assert(rows_ == rhs.rows_);
162 assert(cols_ == rhs.cols_);
171 return addII(*
this, rhs);
173 return addDI(rhs, *
this);
175 return addSI(rhs, *
this);
177 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << rhs.type_);
184 return addDI(*
this, rhs);
186 return addDD(*
this, rhs);
188 return addSD(rhs, *
this);
190 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << rhs.type_);
197 return addSI(*
this, rhs);
199 return addSD(*
this, rhs);
201 return addSS(*
this, rhs);
203 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << rhs.type_);
206 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << rhs.type_);
222 assert(cols_ == rhs.rows_);
235 return mulDD(*
this, rhs);
237 return mulDS(*
this, rhs);
239 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << rhs.type_);
248 return mulSD(*
this, rhs);
250 return mulSS(*
this, rhs);
252 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << rhs.type_);
255 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << rhs.type_);
278 *
this = *
this + (rhs * -1.0);
300 retval.type_ = Diagonal;
301 retval.diag_.assign(rows_, rhs);
307 for (
double& elem : retval.diag_) {
315 retval.sparse_ *= rhs;
319 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << type_);
341 retval.type_ = Diagonal;
342 retval.diag_.assign(rows_, 1.0/rhs);
348 for (
double& elem : retval.diag_) {
356 retval.sparse_ /= rhs;
360 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << type_);
374 Eigen::VectorXd
operator*(
const Eigen::VectorXd& rhs)
const
376 assert(cols_ == rhs.size());
379 return Eigen::VectorXd::Zero(rows_);
384 const Eigen::VectorXd d = Eigen::Map<const Eigen::VectorXd>(diag_.data(), rows_);
385 return d.cwiseProduct(rhs);
388 return sparse_ * rhs;
390 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << type_);
401 assert(lhs.type_ == Identity);
402 assert(rhs.type_ == Identity);
404 retval.type_ = Diagonal;
405 retval.rows_ = lhs.rows_;
406 retval.cols_ = rhs.cols_;
407 retval.diag_.assign(lhs.rows_, 2.0);
414 static_cast<void>(rhs);
415 assert(lhs.type_ == Diagonal);
416 assert(rhs.type_ == Identity);
418 for (
int r = 0; r < lhs.rows_; ++r) {
419 retval.diag_[r] += 1.0;
427 assert(lhs.type_ == Diagonal);
428 assert(rhs.type_ == Diagonal);
430 for (
int r = 0; r < lhs.rows_; ++r) {
431 retval.diag_[r] += rhs.diag_[r];
439 static_cast<void>(rhs);
440 assert(lhs.type_ == Sparse);
441 assert(rhs.type_ == Identity);
443 retval.sparse_ += spdiag(Eigen::VectorXd::Ones(lhs.rows_));;
450 assert(lhs.type_ == Sparse);
451 assert(rhs.type_ == Diagonal);
453 retval.sparse_ += spdiag(rhs.diag_);
460 assert(lhs.type_ == Sparse);
461 assert(rhs.type_ == Sparse);
463 retval.sparse_ += rhs.sparse_;
473 assert(lhs.type_ == Diagonal);
474 assert(rhs.type_ == Diagonal);
476 for (
int r = 0; r < lhs.rows_; ++r) {
477 retval.diag_[r] *= rhs.diag_[r];
485 assert(lhs.type_ == Diagonal);
486 assert(rhs.type_ == Sparse);
488 retval.type_ = Sparse;
489 retval.rows_ = lhs.rows_;
490 retval.cols_ = rhs.cols_;
498 assert(lhs.type_ == Sparse);
499 assert(rhs.type_ == Diagonal);
501 retval.type_ = Sparse;
502 retval.rows_ = lhs.rows_;
503 retval.cols_ = rhs.cols_;
511 assert(lhs.type_ == Sparse);
512 assert(rhs.type_ == Sparse);
514 retval.type_ = Sparse;
515 retval.rows_ = lhs.rows_;
516 retval.cols_ = rhs.cols_;
529 template<
class Scalar,
int Options,
class Index>
530 void toSparse(Eigen::SparseMatrix<Scalar, Options, Index>& s)
const
534 s = Eigen::SparseMatrix<Scalar, Options, Index>(rows_, cols_);
537 s = spdiag(Eigen::VectorXd::Ones(rows_));
586 return sparse_.nonZeros();
588 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << type_);
598 double coeff(
const int row,
const int col)
const
604 return (row == col) ? 1.0 : 0.0;
606 return (row == col) ? diag_[row] : 0.0;
608 return sparse_.coeff(row, col);
610 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << type_);
624 if (type_ != Sparse) {
632 SparseRep& mutable_sparse =
const_cast<SparseRep&
>(sparse_);
640 enum AudoDiffMatrixType { Zero, Identity, Diagonal, Sparse };
642 AudoDiffMatrixType type_;
677 const int n = d.size();
679 mat.reserve(Eigen::ArrayXi::Ones(n, 1));
680 for (SparseRep::Index i = 0; i < n; ++i) {
682 mat.insert(i, i) = d[i];
724 #endif // OPM_AUTODIFFMATRIX_HEADER_INCLUDED
AutoDiffMatrix operator*(const AutoDiffMatrix &rhs) const
Definition: AutoDiffMatrix.hpp:220
AutoDiff< Scalar > operator*(const AutoDiff< Scalar > &lhs, const AutoDiff< Scalar > &rhs)
Definition: AutoDiff.hpp:211
static AutoDiffMatrix mulDD(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:471
static AutoDiffMatrix createIdentity(const int num_rows_cols)
Definition: AutoDiffMatrix.hpp:81
Eigen::VectorXd operator*(const Eigen::VectorXd &rhs) const
Definition: AutoDiffMatrix.hpp:374
Definition: AdditionalObjectDeleter.hpp:22
static AutoDiffMatrix addDI(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:412
AutoDiffMatrix & operator+=(const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:265
int cols() const
Definition: AutoDiffMatrix.hpp:563
const SparseRep & getSparse() const
Definition: AutoDiffMatrix.hpp:623
void fastSparseDiagProduct(const Eigen::SparseMatrix< double > &lhs, const std::vector< double > &rhs, Eigen::SparseMatrix< double > &res)
Definition: fastSparseProduct.hpp:164
void toSparse(Eigen::SparseMatrix< Scalar, Options, Index > &s) const
Definition: AutoDiffMatrix.hpp:530
static AutoDiffMatrix mulSS(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:509
double coeff(const int row, const int col) const
Definition: AutoDiffMatrix.hpp:598
int nonZeros() const
Definition: AutoDiffMatrix.hpp:576
void fastSparseProduct(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs, AutoDiffMatrix &res)
Definition: AutoDiffMatrix.hpp:696
AutoDiffMatrix operator*(const double rhs) const
Definition: AutoDiffMatrix.hpp:292
Eigen::SparseMatrix< double > SparseRep
Definition: AutoDiffMatrix.hpp:47
AutoDiffMatrix & operator-=(const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:276
static AutoDiffMatrix addSD(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:448
void swap(AutoDiffMatrix &other)
Definition: AutoDiffMatrix.hpp:142
Definition: AutoDiffMatrix.hpp:43
AutoDiffMatrix operator+(const AutoDiffMatrix &rhs) const
Definition: AutoDiffMatrix.hpp:159
AutoDiffMatrix & operator=(const AutoDiffMatrix &other)=default
ADB::V V
Definition: BlackoilModelBase_impl.hpp:84
static AutoDiffMatrix mulSD(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:496
AutoDiffMatrix & operator=(AutoDiffMatrix &&other)
Definition: AutoDiffMatrix.hpp:134
std::vector< double > DiagRep
Definition: AutoDiffMatrix.hpp:46
AutoDiffMatrix(AutoDiffMatrix &&other)
Definition: AutoDiffMatrix.hpp:122
AutoDiffMatrix operator/(const double rhs) const
Definition: AutoDiffMatrix.hpp:333
static AutoDiffMatrix addSS(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:458
static AutoDiffMatrix addSI(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:437
AutoDiffMatrix()
Definition: AutoDiffMatrix.hpp:53
static AutoDiffMatrix mulDS(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:483
static AutoDiffMatrix addDD(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:425
static AutoDiffMatrix addII(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs)
Definition: AutoDiffMatrix.hpp:399
AutoDiffMatrix(const Eigen::SparseMatrix< double > &s)
Definition: AutoDiffMatrix.hpp:106
AutoDiffMatrix(const int num_rows, const int num_cols)
Definition: AutoDiffMatrix.hpp:66
AutoDiffMatrix(const Eigen::DiagonalMatrix< double, Eigen::Dynamic > &d)
Definition: AutoDiffMatrix.hpp:91
void fastDiagSparseProduct(const std::vector< double > &lhs, const Eigen::SparseMatrix< double > &rhs, Eigen::SparseMatrix< double > &res)
Definition: fastSparseProduct.hpp:145
int rows() const
Definition: AutoDiffMatrix.hpp:553