OGS
|
Some general linear algebra functionality.
By using the provided functions linear algebra capabilities can be used for different matrix and vector types in a way that is agnostic towards the specific type used.
For documentation, refer to that of the templated method. All specializations or overload must behave in the same way.
Functions | |
void | setLocalAccessibleVector (PETScVector const &x) |
void | set (PETScVector &x, PetscScalar const a) |
void | copy (PETScVector const &x, PETScVector &y) |
void | scale (PETScVector &x, PetscScalar const a) |
void | aypx (PETScVector &y, PetscScalar const a, PETScVector const &x) |
void | axpy (PETScVector &y, PetscScalar const a, PETScVector const &x) |
void | axpby (PETScVector &y, PetscScalar const a, PetscScalar const b, PETScVector const &x) |
template<> | |
void | componentwiseDivide (PETScVector &w, PETScVector const &x, PETScVector const &y) |
template<> | |
double | norm1 (PETScVector const &x) |
template<> | |
double | norm2 (PETScVector const &x) |
template<> | |
double | normMax (PETScVector const &x) |
void | copy (PETScMatrix const &A, PETScMatrix &B) |
void | scale (PETScMatrix &A, PetscScalar const a) |
void | aypx (PETScMatrix &Y, PetscScalar const a, PETScMatrix const &X) |
void | axpy (PETScMatrix &Y, PetscScalar const a, PETScMatrix const &X) |
void | matMult (PETScMatrix const &A, PETScVector const &x, PETScVector &y) |
void | matMultAdd (PETScMatrix const &A, PETScVector const &v1, PETScVector const &v2, PETScVector &v3) |
void | linearSysNormalize (PETScMatrix const &, PETScMatrix &, PETScVector const &, PETScVector &) |
void | finalizeAssembly (PETScMatrix &A) |
void | finalizeAssembly (PETScVector &x) |
template<typename MatrixOrVector > | |
void | copy (MatrixOrVector const &x, MatrixOrVector &y) |
Copies x to y . | |
template<typename MatrixOrVector > | |
void | scale (MatrixOrVector &x, double const a) |
Scales x by a . | |
template<typename MatrixOrVector > | |
void | aypx (MatrixOrVector &y, double const a, MatrixOrVector const &x) |
Computes \( y = a \cdot y + x \). | |
template<typename MatrixOrVector > | |
void | axpy (MatrixOrVector &y, double const a, MatrixOrVector const &x) |
Computes \( y = a \cdot x + y \). | |
template<typename MatrixOrVector > | |
void | axpby (MatrixOrVector &y, double const a, double const b, MatrixOrVector const &x) |
Computes \( y = a \cdot x + b \cdot y \). | |
template<typename MatrixOrVector > | |
void | componentwiseDivide (MatrixOrVector &w, MatrixOrVector const &x, MatrixOrVector const &y) |
Computes \(w = x/y\) componentwise. | |
template<typename MatrixOrVector > | |
double | norm1 (MatrixOrVector const &x) |
Computes the Manhattan norm of x . | |
template<typename MatrixOrVector > | |
double | norm2 (MatrixOrVector const &x) |
Computes the Euclidean norm of x . | |
template<typename MatrixOrVector > | |
double | normMax (MatrixOrVector const &x) |
Computes the maximum norm of x . | |
template<typename MatrixOrVector > | |
double | norm (MatrixOrVector const &x, MathLib::VecNormType type) |
template<typename Matrix > | |
void | finalizeAssembly (Matrix &) |
template<typename Matrix , typename Vector > | |
void | matMult (Matrix const &A, Vector const &x, Vector &y) |
template<typename Matrix , typename Vector > | |
void | matMultAdd (Matrix const &A, Vector const &v1, Vector const &v2, Vector &v3) |
template<typename VectorType > | |
double | computeRelativeNorm (VectorType const &x, VectorType const &y, MathLib::VecNormType norm_type) |
void MathLib::LinAlg::axpby | ( | MatrixOrVector & | y, |
double const | a, | ||
double const | b, | ||
MatrixOrVector const & | x ) |
void MathLib::LinAlg::axpby | ( | PETScVector & | y, |
PetscScalar const | a, | ||
PetscScalar const | b, | ||
PETScVector const & | x ) |
Definition at line 64 of file LinAlg.cpp.
References MathLib::PETScVector::getRawVector().
void MathLib::LinAlg::axpy | ( | MatrixOrVector & | y, |
double const | a, | ||
MatrixOrVector const & | x ) |
void MathLib::LinAlg::axpy | ( | PETScMatrix & | Y, |
PetscScalar const | a, | ||
PETScMatrix const & | X ) |
Definition at line 139 of file LinAlg.cpp.
References MathLib::PETScMatrix::getRawMatrix().
void MathLib::LinAlg::axpy | ( | PETScVector & | y, |
PetscScalar const | a, | ||
PETScVector const & | x ) |
Definition at line 57 of file LinAlg.cpp.
References MathLib::PETScVector::getRawVector().
Referenced by ProcessLib::AssemblyMixin< Process >::assembleWithJacobian(), NumLib::NonlinearSolver< NonlinearSolverTag::Picard >::calculateNonEquilibriumInitialResiduum(), NumLib::StaggeredCoupling::checkCouplingConvergence(), NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeA(), computeRelativeNorm(), NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeResidual(), NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::solve(), and NumLib::NonlinearSolver< NonlinearSolverTag::Picard >::solve().
void MathLib::LinAlg::aypx | ( | MatrixOrVector & | y, |
double const | a, | ||
MatrixOrVector const & | x ) |
void MathLib::LinAlg::aypx | ( | PETScMatrix & | Y, |
PetscScalar const | a, | ||
PETScMatrix const & | X ) |
Definition at line 131 of file LinAlg.cpp.
References MathLib::PETScMatrix::getRawMatrix().
void MathLib::LinAlg::aypx | ( | PETScVector & | y, |
PetscScalar const | a, | ||
PETScVector const & | x ) |
Definition at line 50 of file LinAlg.cpp.
References MathLib::PETScVector::getRawVector().
Referenced by ProcessLib::ComponentTransport::ComponentTransportProcess::solveReactionEquation().
void MathLib::LinAlg::componentwiseDivide | ( | MatrixOrVector & | w, |
MatrixOrVector const & | x, | ||
MatrixOrVector const & | y ) |
Computes \(w = x/y\) componentwise.
void MathLib::LinAlg::componentwiseDivide | ( | PETScVector & | w, |
PETScVector const & | x, | ||
PETScVector const & | y ) |
Definition at line 81 of file LinAlg.cpp.
References MathLib::PETScVector::getRawVector().
Referenced by NumLib::LocalLinearLeastSquaresExtrapolator::extrapolate().
double MathLib::LinAlg::computeRelativeNorm | ( | VectorType const & | x, |
VectorType const & | y, | ||
MathLib::VecNormType | norm_type ) |
Compute the relative norm between two vectors by \( \|x-y\|/\|x\| \).
x | Vector x |
y | Vector y |
norm_type | The norm type of global vector |
Definition at line 297 of file LinAlg.h.
References axpy(), copy(), MathLib::INVALID, norm(), and OGS_FATAL.
Referenced by ProcessLib::TimeLoop::computeTimeStepping().
void MathLib::LinAlg::copy | ( | MatrixOrVector const & | x, |
MatrixOrVector & | y ) |
void MathLib::LinAlg::copy | ( | PETScMatrix const & | A, |
PETScMatrix & | B ) |
Definition at line 119 of file LinAlg.cpp.
void MathLib::LinAlg::copy | ( | PETScVector const & | x, |
PETScVector & | y ) |
Definition at line 37 of file LinAlg.cpp.
References MathLib::PETScVector::getRawVector(), and MathLib::PETScVector::shallowCopy().
Referenced by ProcessLib::AssembledMatrixCache::assemble(), NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeA(), NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeJacobian(), computeRelativeNorm(), NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeResidual(), ProcessLib::TimeLoop::computeTimeStepping(), NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >::getResidual(), NumLib::SimpleMatrixVectorProvider::getVector(), NumLib::SimpleMatrixVectorProvider::getVector(), NumLib::BackwardEuler::getWeightedOldX(), NumLib::StaggeredCoupling::initializeCoupledSolutions(), NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::normalizeAandRhs(), ProcessLib::Process::setInitialConditions(), NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::solve(), NumLib::NonlinearSolver< NonlinearSolverTag::Picard >::solve(), NumLib::PETScNonlinearSolver::solve(), ProcessLib::ComponentTransport::ComponentTransportProcess::solveReactionEquation(), ProcessLib::PhaseField::PhaseFieldProcess< DisplacementDim >::updateConstraints(), and NumLib::StaggeredCoupling::updatePreviousSolution().
void MathLib::LinAlg::finalizeAssembly | ( | Matrix & | ) |
void MathLib::LinAlg::finalizeAssembly | ( | PETScMatrix & | A | ) |
Definition at line 198 of file LinAlg.cpp.
References MathLib::PETScMatrix::finalizeAssembly().
Referenced by detail::applyKnownSolutions(), NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >::assemble(), NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Picard >::assemble(), NumLib::NonlinearSolver< NonlinearSolverTag::Newton >::calculateNonEquilibriumInitialResiduum(), NumLib::NonlinearSolver< NonlinearSolverTag::Picard >::calculateNonEquilibriumInitialResiduum(), NumLib::LocalLinearLeastSquaresExtrapolator::calculateResiduals(), NumLib::LocalLinearLeastSquaresExtrapolator::extrapolate(), ProcessLib::Process::setInitialConditions(), and ProcessLib::SmallDeformation::writeMaterialForces().
void MathLib::LinAlg::finalizeAssembly | ( | PETScVector & | x | ) |
Definition at line 203 of file LinAlg.cpp.
References MathLib::PETScVector::finalizeAssembly().
void MathLib::LinAlg::linearSysNormalize | ( | PETScMatrix const & | , |
PETScMatrix & | , | ||
PETScVector const & | , | ||
PETScVector & | ) |
Definition at line 170 of file LinAlg.cpp.
References OGS_FATAL.
Referenced by NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::normalizeAandRhs().
void MathLib::LinAlg::matMult | ( | Matrix const & | A, |
Vector const & | x, | ||
Vector & | y ) |
void MathLib::LinAlg::matMult | ( | PETScMatrix const & | A, |
PETScVector const & | x, | ||
PETScVector & | y ) |
Definition at line 149 of file LinAlg.cpp.
References MathLib::PETScMatrix::getRawMatrix(), MathLib::PETScVector::getRawVector(), and MathLib::PETScVector::shallowCopy().
Referenced by NumLib::NonlinearSolver< NonlinearSolverTag::Picard >::calculateNonEquilibriumInitialResiduum(), NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeResidual(), NumLib::NonlinearSolver< NonlinearSolverTag::Picard >::solve(), and ProcessLib::ComponentTransport::ComponentTransportProcess::solveReactionEquation().
void MathLib::LinAlg::matMultAdd | ( | Matrix const & | A, |
Vector const & | v1, | ||
Vector const & | v2, | ||
Vector & | v3 ) |
void MathLib::LinAlg::matMultAdd | ( | PETScMatrix const & | A, |
PETScVector const & | v1, | ||
PETScVector const & | v2, | ||
PETScVector & | v3 ) |
Definition at line 159 of file LinAlg.cpp.
References MathLib::PETScMatrix::getRawMatrix(), MathLib::PETScVector::getRawVector(), and MathLib::PETScVector::shallowCopy().
Referenced by NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeResidual(), and NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeRhs().
double MathLib::LinAlg::norm | ( | MatrixOrVector const & | x, |
MathLib::VecNormType | type ) |
Definition at line 94 of file LinAlg.h.
References MathLib::INFINITY_N, MathLib::NORM1, norm1(), MathLib::NORM2, norm2(), normMax(), and OGS_FATAL.
Referenced by NumLib::ConvergenceCriterionDeltaX::checkDeltaX(), NumLib::ConvergenceCriterionResidual::checkDeltaX(), NumLib::ConvergenceCriterionResidual::checkResidual(), computeRelativeNorm(), norm1(), norm2(), and normMax().
double MathLib::LinAlg::norm1 | ( | MatrixOrVector const & | x | ) |
Computes the Manhattan norm of x
.
double MathLib::LinAlg::norm1 | ( | PETScVector const & | x | ) |
Definition at line 90 of file LinAlg.cpp.
References MathLib::PETScVector::getRawVector(), and norm().
Referenced by norm().
double MathLib::LinAlg::norm2 | ( | MatrixOrVector const & | x | ) |
Computes the Euclidean norm of x
.
double MathLib::LinAlg::norm2 | ( | PETScVector const & | x | ) |
Definition at line 100 of file LinAlg.cpp.
References MathLib::PETScVector::getRawVector(), and norm().
Referenced by norm().
double MathLib::LinAlg::normMax | ( | MatrixOrVector const & | x | ) |
Computes the maximum norm of x
.
double MathLib::LinAlg::normMax | ( | PETScVector const & | x | ) |
Definition at line 110 of file LinAlg.cpp.
References MathLib::PETScVector::getRawVector(), and norm().
Referenced by norm().
void MathLib::LinAlg::scale | ( | MatrixOrVector & | x, |
double const | a ) |
void MathLib::LinAlg::scale | ( | PETScMatrix & | A, |
PetscScalar const | a ) |
Definition at line 125 of file LinAlg.cpp.
References MathLib::PETScMatrix::getRawMatrix().
void MathLib::LinAlg::scale | ( | PETScVector & | x, |
PetscScalar const | a ) |
Definition at line 44 of file LinAlg.cpp.
References MathLib::PETScVector::getRawVector().
Referenced by NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeResidual(), NumLib::TimeDiscretizedODESystem< ODESystemTag::FirstOrderImplicitQuasilinear, NonlinearSolverTag::Newton >::getResidual(), NumLib::BackwardEuler::getWeightedOldX(), and ProcessLib::PhaseField::PhaseFieldProcess< DisplacementDim >::postNonLinearSolverConcreteProcess().
void MathLib::LinAlg::set | ( | PETScVector & | x, |
PetscScalar const | a ) |
Definition at line 32 of file LinAlg.cpp.
References MathLib::PETScVector::getRawVector().
Referenced by ProcessLib::SmallDeformation::writeMaterialForces().
void MathLib::LinAlg::setLocalAccessibleVector | ( | PETScVector const & | x | ) |
Set local accessible vector in order to get entries. Call this before call operator[] or get(...) of x.
Definition at line 27 of file LinAlg.cpp.
References MathLib::PETScVector::setLocalAccessibleVector().
Referenced by NumLib::LocalLinearLeastSquaresExtrapolator::calculateResidualElement(), ProcessLib::Process::computeSecondaryVariable(), ChemistryLib::PhreeqcIOData::PhreeqcIO::getConcentration(), ProcessLib::Process::postIteration(), ProcessLib::TES::TESProcess::postIterationConcreteProcess(), ProcessLib::LIE::HydroMechanics::HydroMechanicsProcess< GlobalDim >::postTimestepConcreteProcess(), ProcessLib::Process::preIteration(), ProcessLib::Process::preOutput(), ProcessLib::Process::preTimestep(), ChemistryLib::PhreeqcIOData::AqueousSolution::print(), ChemistryLib::PhreeqcIOData::anonymous_namespace{PhreeqcIO.cpp}::setAqueousSolution(), ProcessLib::Process::setInitialConditions(), anonymous_namespace{Process.cpp}::setLocalAccessibleVectors(), NumLib::NonlinearSolver< NonlinearSolverTag::Picard >::solve(), ProcessLib::ComponentTransport::ComponentTransportProcess::solveReactionEquation(), NumLib::transformVariableFromGlobalVector(), NumLib::transformVariableFromGlobalVector(), ProcessLib::PhaseField::PhaseFieldProcess< DisplacementDim >::updateConstraints(), and ProcessLib::SmallDeformation::writeMaterialForces().