OGS
MathLib::LinAlg Namespace Reference

Detailed Description

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 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)
 

Function Documentation

◆ axpby() [1/2]

template<typename MatrixOrVector >
void MathLib::LinAlg::axpby ( MatrixOrVector &  y,
double const  a,
double const  b,
MatrixOrVector const &  x 
)

Computes \( y = a \cdot x + b \cdot y \).

Definition at line 69 of file LinAlg.h.

71{
72 y = a*x + b*y;
73}

◆ axpby() [2/2]

void MathLib::LinAlg::axpby ( PETScVector y,
PetscScalar const  a,
PetscScalar const  b,
PETScVector const &  x 
)

Definition at line 64 of file LinAlg.cpp.

66{
67 // TODO check sizes
68 VecAXPBY(y.getRawVector(), a, b, x.getRawVector());
69}
PETSc_Vec & getRawVector()
Exposes the underlying PETSc vector.
Definition: PETScVector.h:171

References MathLib::PETScVector::getRawVector().

◆ axpy() [1/3]

template<typename MatrixOrVector >
void MathLib::LinAlg::axpy ( MatrixOrVector &  y,
double const  a,
MatrixOrVector const &  x 
)

Computes \( y = a \cdot x + y \).

Definition at line 62 of file LinAlg.h.

63{
64 y += a*x;
65}

◆ axpy() [2/3]

void MathLib::LinAlg::axpy ( PETScMatrix Y,
PetscScalar const  a,
PETScMatrix const &  X 
)

Definition at line 132 of file LinAlg.cpp.

133{
134 // TODO check sizes
135 // TODO sparsity pattern, currently they are assumed to be different (slow)
136 MatAXPY(Y.getRawMatrix(), a, X.getRawMatrix(), DIFFERENT_NONZERO_PATTERN);
137}
Mat & getRawMatrix()
Get matrix reference.
Definition: PETScMatrix.h:86

References MathLib::PETScMatrix::getRawMatrix().

◆ axpy() [3/3]

◆ aypx() [1/3]

template<typename MatrixOrVector >
void MathLib::LinAlg::aypx ( MatrixOrVector &  y,
double const  a,
MatrixOrVector const &  x 
)

Computes \( y = a \cdot y + x \).

Definition at line 55 of file LinAlg.h.

56{
57 y = a*y + x;
58}

◆ aypx() [2/3]

void MathLib::LinAlg::aypx ( PETScMatrix Y,
PetscScalar const  a,
PETScMatrix const &  X 
)

Definition at line 124 of file LinAlg.cpp.

125{
126 // TODO check sizes
127 // TODO sparsity pattern, currently they are assumed to be different (slow)
128 MatAYPX(Y.getRawMatrix(), a, X.getRawMatrix(), DIFFERENT_NONZERO_PATTERN);
129}

References MathLib::PETScMatrix::getRawMatrix().

◆ aypx() [3/3]

void MathLib::LinAlg::aypx ( PETScVector y,
PetscScalar const  a,
PETScVector const &  x 
)

Definition at line 50 of file LinAlg.cpp.

51{
52 // TODO check sizes
53 VecAYPX(y.getRawVector(), a, x.getRawVector());
54}

References MathLib::PETScVector::getRawVector().

Referenced by ProcessLib::ComponentTransport::ComponentTransportProcess::solveReactionEquation().

◆ componentwiseDivide() [1/2]

template<typename MatrixOrVector >
void MathLib::LinAlg::componentwiseDivide ( MatrixOrVector &  w,
MatrixOrVector const &  x,
MatrixOrVector const &  y 
)

Computes \(w = x/y\) componentwise.

◆ componentwiseDivide() [2/2]

template<>
void MathLib::LinAlg::componentwiseDivide ( PETScVector w,
PETScVector const &  x,
PETScVector const &  y 
)

Definition at line 74 of file LinAlg.cpp.

76{
77 VecPointwiseDivide(w.getRawVector(), x.getRawVector(), y.getRawVector());
78}

References MathLib::PETScVector::getRawVector().

Referenced by NumLib::LocalLinearLeastSquaresExtrapolator::extrapolate().

◆ copy() [1/3]

template<typename MatrixOrVector >
void MathLib::LinAlg::copy ( MatrixOrVector const &  x,
MatrixOrVector &  y 
)

Copies x to y.

Definition at line 41 of file LinAlg.h.

42{
43 y = x;
44}

◆ copy() [2/3]

void MathLib::LinAlg::copy ( PETScMatrix const &  A,
PETScMatrix B 
)

Definition at line 112 of file LinAlg.cpp.

113{
114 B = A;
115}

◆ copy() [3/3]

◆ finalizeAssembly() [1/3]

template<typename Matrix >
void MathLib::LinAlg::finalizeAssembly ( Matrix &  )

◆ finalizeAssembly() [2/3]

◆ finalizeAssembly() [3/3]

void MathLib::LinAlg::finalizeAssembly ( PETScVector x)

Definition at line 168 of file LinAlg.cpp.

169{
171}
void finalizeAssembly()
Perform MPI collection of assembled entries in buffer.

References MathLib::PETScVector::finalizeAssembly().

◆ matMult() [1/2]

template<typename Matrix , typename Vector >
void MathLib::LinAlg::matMult ( Matrix const &  A,
Vector const &  x,
Vector &  y 
)

Computes \( y = A \cdot x \).

Note
x must not be the same object as y. This restirction has been chosen in order to fulfill the requirements of the respective PETSc function.

Definition at line 119 of file LinAlg.h.

120{
121 assert(&x != &y);
122 y = A*x;
123}

◆ matMult() [2/2]

◆ matMultAdd() [1/2]

template<typename Matrix , typename Vector >
void MathLib::LinAlg::matMultAdd ( Matrix const &  A,
Vector const &  v1,
Vector const &  v2,
Vector &  v3 
)

Computes \( v_3 = A \cdot v_1 + v_2 \).

Note
x must not be the same object as y. This restirction has been chosen in order to fulfill the requirements of the respective PETSc function.

Definition at line 132 of file LinAlg.h.

133{
134 assert(&v1 != &v3);
135 v3 = v2 + A*v1;
136}

◆ matMultAdd() [2/2]

void MathLib::LinAlg::matMultAdd ( PETScMatrix const &  A,
PETScVector const &  v1,
PETScVector const &  v2,
PETScVector v3 
)

Definition at line 152 of file LinAlg.cpp.

154{
155 // TODO check sizes
156 assert(&v1 != &v3);
157 if (!v3.getRawVector())
158 v3.shallowCopy(v1);
159 MatMultAdd(A.getRawMatrix(), v1.getRawVector(), v2.getRawVector(),
160 v3.getRawVector());
161}

References MathLib::PETScMatrix::getRawMatrix(), MathLib::PETScVector::getRawVector(), and MathLib::PETScVector::shallowCopy().

Referenced by NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeResidual(), and NumLib::MatrixTranslatorGeneral< ODESystemTag::FirstOrderImplicitQuasilinear >::computeRhs().

◆ norm()

template<typename MatrixOrVector >
double MathLib::LinAlg::norm ( MatrixOrVector const &  x,
MathLib::VecNormType  type 
)

Definition at line 93 of file LinAlg.h.

94{
95 switch (type) {
97 return norm1(x);
99 return norm2(x);
101 return normMax(x);
102 default:
103 OGS_FATAL("Invalid norm type given.");
104 }
105}
#define OGS_FATAL(...)
Definition: Error.h:26
double norm1(PETScVector const &x)
Definition: LinAlg.cpp:83
double normMax(PETScVector const &x)
Definition: LinAlg.cpp:103
double norm2(PETScVector const &x)
Definition: LinAlg.cpp:93

References MathLib::INFINITY_N, norm1(), MathLib::NORM1, norm2(), MathLib::NORM2, normMax(), and OGS_FATAL.

Referenced by NumLib::ConvergenceCriterionDeltaX::checkDeltaX(), NumLib::ConvergenceCriterionResidual::checkDeltaX(), NumLib::ConvergenceCriterionResidual::checkResidual(), NumLib::computeRelativeNorm(), norm1(), norm2(), and normMax().

◆ norm1() [1/2]

template<typename MatrixOrVector >
double MathLib::LinAlg::norm1 ( MatrixOrVector const &  x)

Computes the Manhattan norm of x.

◆ norm1() [2/2]

template<>
double MathLib::LinAlg::norm1 ( PETScVector const &  x)

Definition at line 83 of file LinAlg.cpp.

84{
85 PetscScalar norm = 0.;
86 VecNorm(x.getRawVector(), NORM_1, &norm);
87 return norm;
88}
double norm(MatrixOrVector const &x, MathLib::VecNormType type)
Definition: LinAlg.h:93

References MathLib::PETScVector::getRawVector(), and norm().

Referenced by norm().

◆ norm2() [1/2]

template<typename MatrixOrVector >
double MathLib::LinAlg::norm2 ( MatrixOrVector const &  x)

Computes the Euclidean norm of x.

◆ norm2() [2/2]

template<>
double MathLib::LinAlg::norm2 ( PETScVector const &  x)

Definition at line 93 of file LinAlg.cpp.

94{
95 PetscScalar norm = 0.;
96 VecNorm(x.getRawVector(), NORM_2, &norm);
97 return norm;
98}

References MathLib::PETScVector::getRawVector(), and norm().

Referenced by norm().

◆ normMax() [1/2]

template<typename MatrixOrVector >
double MathLib::LinAlg::normMax ( MatrixOrVector const &  x)

Computes the maximum norm of x.

◆ normMax() [2/2]

template<>
double MathLib::LinAlg::normMax ( PETScVector const &  x)

Definition at line 103 of file LinAlg.cpp.

104{
105 PetscScalar norm = 0.;
106 VecNorm(x.getRawVector(), NORM_INFINITY, &norm);
107 return norm;
108}

References MathLib::PETScVector::getRawVector(), and norm().

Referenced by norm().

◆ scale() [1/3]

template<typename MatrixOrVector >
void MathLib::LinAlg::scale ( MatrixOrVector &  x,
double const  a 
)

Scales x by a.

Definition at line 48 of file LinAlg.h.

49{
50 x *= a;
51}

◆ scale() [2/3]

void MathLib::LinAlg::scale ( PETScMatrix A,
PetscScalar const  a 
)

Definition at line 118 of file LinAlg.cpp.

119{
120 MatScale(A.getRawMatrix(), a);
121}

References MathLib::PETScMatrix::getRawMatrix().

◆ scale() [3/3]

◆ set()

void MathLib::LinAlg::set ( PETScVector x,
PetscScalar const  a 
)

Definition at line 32 of file LinAlg.cpp.

33{
34 VecSet(x.getRawVector(), a);
35}

References MathLib::PETScVector::getRawVector().

Referenced by ProcessLib::SmallDeformation::writeMaterialForces().

◆ setLocalAccessibleVector()

void MathLib::LinAlg::setLocalAccessibleVector ( PETScVector const &  x)