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> 31 #include <opm/autodiff/fastSparseOperations.hpp> 46 typedef std::vector<double> DiagRep;
47 typedef Eigen::SparseMatrix<double> SparseRep;
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_);
267 if( type_ == Sparse && rhs.type_ == Sparse )
269 fastSparseAdd( sparse_, rhs.sparse_ );
284 if( type_ == Sparse && rhs.type_ == Sparse )
286 fastSparseSubstract( sparse_, rhs.sparse_ );
289 *
this = *
this + (rhs * -1.0);
312 retval.type_ = Diagonal;
313 retval.diag_.assign(rows_, rhs);
319 for (
double& elem : retval.diag_) {
327 retval.sparse_ *= rhs;
331 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << type_);
353 retval.type_ = Diagonal;
354 retval.diag_.assign(rows_, 1.0/rhs);
360 for (
double& elem : retval.diag_) {
368 retval.sparse_ /= rhs;
372 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << type_);
386 Eigen::VectorXd
operator*(
const Eigen::VectorXd& rhs)
const 388 assert(cols_ == rhs.size());
391 return Eigen::VectorXd::Zero(rows_);
396 const Eigen::VectorXd d = Eigen::Map<const Eigen::VectorXd>(diag_.data(), rows_);
397 return d.cwiseProduct(rhs);
400 return sparse_ * rhs;
402 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << type_);
413 assert(lhs.type_ == Identity);
414 assert(rhs.type_ == Identity);
416 retval.type_ = Diagonal;
417 retval.rows_ = lhs.rows_;
418 retval.cols_ = rhs.cols_;
419 retval.diag_.assign(lhs.rows_, 2.0);
426 static_cast<void>(rhs);
427 assert(lhs.type_ == Diagonal);
428 assert(rhs.type_ == Identity);
430 for (
int r = 0; r < lhs.rows_; ++r) {
431 retval.diag_[r] += 1.0;
439 assert(lhs.type_ == Diagonal);
440 assert(rhs.type_ == Diagonal);
442 for (
int r = 0; r < lhs.rows_; ++r) {
443 retval.diag_[r] += rhs.diag_[r];
451 static_cast<void>(rhs);
452 assert(lhs.type_ == Sparse);
453 assert(rhs.type_ == Identity);
455 retval.sparse_ += spdiag(Eigen::VectorXd::Ones(lhs.rows_));;
462 assert(lhs.type_ == Sparse);
463 assert(rhs.type_ == Diagonal);
465 retval.sparse_ += spdiag(rhs.diag_);
472 assert(lhs.type_ == Sparse);
473 assert(rhs.type_ == Sparse);
475 retval.sparse_ += rhs.sparse_;
485 assert(lhs.type_ == Diagonal);
486 assert(rhs.type_ == Diagonal);
488 for (
int r = 0; r < lhs.rows_; ++r) {
489 retval.diag_[r] *= rhs.diag_[r];
497 assert(lhs.type_ == Diagonal);
498 assert(rhs.type_ == Sparse);
500 retval.type_ = Sparse;
501 retval.rows_ = lhs.rows_;
502 retval.cols_ = rhs.cols_;
503 fastDiagSparseProduct(lhs.diag_, rhs.sparse_, retval.sparse_);
510 assert(lhs.type_ == Sparse);
511 assert(rhs.type_ == Diagonal);
513 retval.type_ = Sparse;
514 retval.rows_ = lhs.rows_;
515 retval.cols_ = rhs.cols_;
516 fastSparseDiagProduct(lhs.sparse_, rhs.diag_, retval.sparse_);
523 assert(lhs.type_ == Sparse);
524 assert(rhs.type_ == Sparse);
526 retval.type_ = Sparse;
527 retval.rows_ = lhs.rows_;
528 retval.cols_ = rhs.cols_;
541 template<
class Scalar,
int Options,
class Index>
542 void toSparse(Eigen::SparseMatrix<Scalar, Options, Index>& s)
const 546 s = Eigen::SparseMatrix<Scalar, Options, Index>(rows_, cols_);
549 s = spdiag(Eigen::VectorXd::Ones(rows_));
598 return sparse_.nonZeros();
600 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << type_);
610 double coeff(
const int row,
const int col)
const 616 return (row == col) ? 1.0 : 0.0;
618 return (row == col) ? diag_[row] : 0.0;
620 return sparse_.coeff(row, col);
622 OPM_THROW(std::logic_error,
"Invalid AutoDiffMatrix type encountered: " << type_);
636 if (type_ != Sparse) {
644 SparseRep& mutable_sparse =
const_cast<SparseRep&
>(sparse_);
652 enum AudoDiffMatrixType { Zero, Identity, Diagonal, Sparse };
654 AudoDiffMatrixType type_;
665 AutoDiffMatrix(AudoDiffMatrixType type,
int rows_arg,
int cols_arg,
666 DiagRep diag=DiagRep(), SparseRep sparse=SparseRep())
689 const int n = d.size();
691 mat.reserve(Eigen::ArrayXi::Ones(n, 1));
692 for (SparseRep::Index i = 0; i < n; ++i) {
694 mat.insert(i, i) = d[i];
736 #endif // OPM_AUTODIFFMATRIX_HEADER_INCLUDED Eigen::VectorXd operator*(const Eigen::VectorXd &rhs) const
Multiplies an AutoDiffMatrix with a vector.
Definition: AutoDiffMatrix.hpp:386
AutoDiffMatrix operator*(const double rhs) const
Multiplies an AutoDiffMatrix with a scalar.
Definition: AutoDiffMatrix.hpp:304
void toSparse(Eigen::SparseMatrix< Scalar, Options, Index > &s) const
Converts the AutoDiffMatrix to an Eigen SparseMatrix.This might be an expensive operation to perform ...
Definition: AutoDiffMatrix.hpp:542
AutoDiffMatrix(const Eigen::SparseMatrix< double > &s)
Creates a sparse matrix from an Eigen sparse matrix.
Definition: AutoDiffMatrix.hpp:106
static AutoDiffMatrix createIdentity(const int num_rows_cols)
Creates an identity matrix with num_rows_cols x num_rows_cols entries.
Definition: AutoDiffMatrix.hpp:81
AutoDiffMatrix is a wrapper class that optimizes matrix operations.
Definition: AutoDiffMatrix.hpp:43
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
AutoDiffMatrix(const Eigen::DiagonalMatrix< double, Eigen::Dynamic > &d)
Creates a diagonal matrix from an Eigen diagonal matrix.
Definition: AutoDiffMatrix.hpp:91
AutoDiffMatrix operator*(const AutoDiffMatrix &rhs) const
Multiplies two AutoDiffMatrices.
Definition: AutoDiffMatrix.hpp:220
const SparseRep & getSparse() const
Returns the sparse representation of this matrix.
Definition: AutoDiffMatrix.hpp:635
int cols() const
Returns number of columns in the matrix.
Definition: AutoDiffMatrix.hpp:575
AutoDiffMatrix operator/(const double rhs) const
Divides an AutoDiffMatrix by a scalar.
Definition: AutoDiffMatrix.hpp:345
int rows() const
Returns number of rows in the matrix.
Definition: AutoDiffMatrix.hpp:565
int nonZeros() const
Returns number of non-zero elements in the matrix.
Definition: AutoDiffMatrix.hpp:588
AutoDiffMatrix()
Creates an empty zero matrix.
Definition: AutoDiffMatrix.hpp:53
AutoDiffMatrix operator+(const AutoDiffMatrix &rhs) const
Adds two AutoDiffMatrices.
Definition: AutoDiffMatrix.hpp:159
double coeff(const int row, const int col) const
Returns element (row, col) in the matrix.
Definition: AutoDiffMatrix.hpp:610
void fastSparseProduct(const AutoDiffMatrix &lhs, const AutoDiffMatrix &rhs, AutoDiffMatrix &res)
Utility function to lessen code changes required elsewhere.
Definition: AutoDiffMatrix.hpp:708
AutoDiffMatrix(const int num_rows, const int num_cols)
Creates a zero matrix with num_rows x num_cols entries.
Definition: AutoDiffMatrix.hpp:66