From caf8776f9b6795c84c6e07121868a59af3a29044 Mon Sep 17 00:00:00 2001
From: Henry Weller <http://cfd.direct>
Date: Tue, 22 Mar 2016 14:13:48 +0000
Subject: [PATCH] SquareMatrix, SymmetricSquareMatrix: Changed the constructor
from size to require only n This avoids the need to check that the m and n
dimensions are the same.
---
applications/test/Matrix/Test-Matrix.C | 14 +++++--
.../viewFactorsGen/viewFactorsGen.C | 3 +-
.../surfaceFeatureExtract.C | 2 +-
src/ODE/ODESolvers/SIBS/SIBS.C | 6 +--
.../matrices/SquareMatrix/SquareMatrix.C | 2 +-
.../matrices/SquareMatrix/SquareMatrix.H | 16 +++++---
.../matrices/SquareMatrix/SquareMatrixI.H | 40 ++++++++++++-------
.../SymmetricSquareMatrix.C | 33 ++++++++++-----
.../SymmetricSquareMatrix.H | 16 +++++---
.../SymmetricSquareMatrixI.H | 35 ++++++++--------
.../lduMatrix/lduMatrix/lduMatrixTemplates.C | 2 +-
.../scalarMatrices/scalarMatricesTemplates.C | 2 +-
.../matrices/simpleMatrix/simpleMatrix.C | 2 +-
.../radiationModels/viewFactor/viewFactor.C | 11 ++---
14 files changed, 110 insertions(+), 74 deletions(-)
diff --git a/applications/test/Matrix/Test-Matrix.C b/applications/test/Matrix/Test-Matrix.C
index ad48481fa15..95cfa4c1af5 100644
--- a/applications/test/Matrix/Test-Matrix.C
+++ b/applications/test/Matrix/Test-Matrix.C
@@ -25,6 +25,7 @@ License
#include "scalarMatrices.H"
#include "vector.H"
+#include "tensor.H"
#include "IFstream.H"
using namespace Foam;
@@ -50,7 +51,7 @@ int main(int argc, char *argv[])
Info<< max(hmm) << endl;
Info<< min(hmm) << endl;
- SquareMatrix<scalar> hmm2(3, 3, 1.0);
+ SquareMatrix<scalar> hmm2(3, I);
hmm = hmm2;
@@ -72,7 +73,14 @@ int main(int argc, char *argv[])
Info<< hmm5 << endl;
{
- scalarSymmetricSquareMatrix symmMatrix(3, 3, 0);
+ RectangularMatrix<scalar> rm1(5, 6, 3.1);
+ rm1(0, 1) = 4.5;
+ RectangularMatrix<scalar> rm1b(rm1.block(2, 2, 0, 0));
+ Info<< "rm1b = " << rm1b << endl;
+ }
+
+ {
+ scalarSymmetricSquareMatrix symmMatrix(3, Zero);
symmMatrix(0, 0) = 4;
symmMatrix(1, 0) = 12;
@@ -104,7 +112,7 @@ int main(int argc, char *argv[])
}
{
- scalarSquareMatrix squareMatrix(3, 3, 0);
+ scalarSquareMatrix squareMatrix(3, Zero);
squareMatrix(0, 0) = 4;
squareMatrix(0, 1) = 12;
diff --git a/applications/utilities/preProcessing/viewFactorsGen/viewFactorsGen.C b/applications/utilities/preProcessing/viewFactorsGen/viewFactorsGen.C
index ed00d1a6cf3..d5f295a62b5 100644
--- a/applications/utilities/preProcessing/viewFactorsGen/viewFactorsGen.C
+++ b/applications/utilities/preProcessing/viewFactorsGen/viewFactorsGen.C
@@ -675,7 +675,6 @@ int main(int argc, char *argv[])
// Matrix sum in j(Fij) for each i (if enclosure sum = 1)
scalarSquareMatrix sumViewFactorPatch
(
- totalPatches,
totalPatches,
0.0
);
@@ -888,7 +887,7 @@ int main(int argc, char *argv[])
if (Pstream::master())
{
- scalarSquareMatrix Fmatrix(totalNCoarseFaces, totalNCoarseFaces, 0.0);
+ scalarSquareMatrix Fmatrix(totalNCoarseFaces, 0.0);
labelListList globalFaceFaces(visibleFaceFaces.size());
diff --git a/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtract.C b/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtract.C
index 510305f28cf..2882e6af2af 100644
--- a/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtract.C
+++ b/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtract.C
@@ -303,7 +303,7 @@ triSurfacePointScalarField calcCurvature
faceCoordSys.normalize();
// Construct the matrix to solve
- scalarSymmetricSquareMatrix T(3, 3, 0);
+ scalarSymmetricSquareMatrix T(3, 0);
scalarDiagonalMatrix Z(3, 0);
// Least Squares
diff --git a/src/ODE/ODESolvers/SIBS/SIBS.C b/src/ODE/ODESolvers/SIBS/SIBS.C
index 36acb8a6b1a..a095a7c5073 100644
--- a/src/ODE/ODESolvers/SIBS/SIBS.C
+++ b/src/ODE/ODESolvers/SIBS/SIBS.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -50,7 +50,7 @@ Foam::SIBS::SIBS(const ODESystem& ode, const dictionary& dict)
:
ODESolver(ode, dict),
a_(iMaxX_, 0.0),
- alpha_(kMaxX_, kMaxX_, 0.0),
+ alpha_(kMaxX_, 0.0),
d_p_(n_, kMaxX_, 0.0),
x_p_(kMaxX_, 0.0),
err_(kMaxX_, 0.0),
@@ -60,7 +60,7 @@ Foam::SIBS::SIBS(const ODESystem& ode, const dictionary& dict)
yErr_(n_, 0.0),
dydx0_(n_),
dfdx_(n_, 0.0),
- dfdy_(n_, n_, 0.0),
+ dfdy_(n_, 0.0),
first_(1),
epsOld_(-1.0)
{}
diff --git a/src/OpenFOAM/matrices/SquareMatrix/SquareMatrix.C b/src/OpenFOAM/matrices/SquareMatrix/SquareMatrix.C
index 74439c15def..cf34d78a212 100644
--- a/src/OpenFOAM/matrices/SquareMatrix/SquareMatrix.C
+++ b/src/OpenFOAM/matrices/SquareMatrix/SquareMatrix.C
@@ -35,7 +35,7 @@ Foam::scalar Foam::detDecomposed
const label sign
)
{
- scalar diagProduct = 1.0;
+ Type diagProduct = pTraits<Type>::one;
for (label i = 0; i < matrix.m(); ++i)
{
diff --git a/src/OpenFOAM/matrices/SquareMatrix/SquareMatrix.H b/src/OpenFOAM/matrices/SquareMatrix/SquareMatrix.H
index 270891a9551..39403194d55 100644
--- a/src/OpenFOAM/matrices/SquareMatrix/SquareMatrix.H
+++ b/src/OpenFOAM/matrices/SquareMatrix/SquareMatrix.H
@@ -38,6 +38,7 @@ SourceFiles
#define SquareMatrix_H
#include "Matrix.H"
+#include "Identity.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -64,14 +65,17 @@ public:
//- Construct given number of rows/columns.
inline SquareMatrix(const label n);
- //- Construct given number of rows and columns,
- // It checks that m == n.
- inline SquareMatrix(const label m, const label n);
+ //- Construct given number of rows/columns
+ // initializing all elements to zero
+ inline SquareMatrix(const label n, const zero);
+
+ //- Construct given number of rows/columns
+ // Initializing to the identity matrix
+ inline SquareMatrix(const label n, const Identity<Type>);
//- Construct with given number of rows and rows
- // and value for all elements.
- // It checks that m == n.
- inline SquareMatrix(const label m, const label n, const Type&);
+ // initializing all elements to the given value
+ inline SquareMatrix(const label n, const Type&);
//- Construct from Istream.
inline SquareMatrix(Istream&);
diff --git a/src/OpenFOAM/matrices/SquareMatrix/SquareMatrixI.H b/src/OpenFOAM/matrices/SquareMatrix/SquareMatrixI.H
index 74d7bdb95c2..ffa595d51da 100644
--- a/src/OpenFOAM/matrices/SquareMatrix/SquareMatrixI.H
+++ b/src/OpenFOAM/matrices/SquareMatrix/SquareMatrixI.H
@@ -31,40 +31,51 @@ inline Foam::SquareMatrix<Type>::SquareMatrix()
Matrix<SquareMatrix<Type>, Type>()
{}
+
template<class Type>
inline Foam::SquareMatrix<Type>::SquareMatrix(const label n)
:
Matrix<SquareMatrix<Type>, Type>(n, n)
{}
+
template<class Type>
-inline Foam::SquareMatrix<Type>::SquareMatrix(const label m, const label n)
+inline Foam::SquareMatrix<Type>::SquareMatrix
+(
+ const label n,
+ const zero
+)
:
- Matrix<SquareMatrix<Type>, Type>(m, n)
+ Matrix<SquareMatrix<Type>, Type>(n, n, Zero)
+{}
+
+
+template<class Type>
+inline Foam::SquareMatrix<Type>::SquareMatrix
+(
+ const label n,
+ const Identity<Type>
+)
+:
+ Matrix<SquareMatrix<Type>, Type>(n, n, Zero)
{
- if (m != n)
+ for (label i = 0; i < n; ++i)
{
- FatalErrorInFunction
- << "m != n for constructing a square matrix" << exit(FatalError);
+ this->operator()(i, i) = I;
}
}
+
template<class Type>
inline Foam::SquareMatrix<Type>::SquareMatrix
(
- const label m,
const label n,
const Type& t
)
:
- Matrix<SquareMatrix<Type>, Type>(m, n, t)
-{
- if (m != n)
- {
- FatalErrorInFunction
- << "m != n for constructing a square matrix" << exit(FatalError);
- }
-}
+ Matrix<SquareMatrix<Type>, Type>(n, n, t)
+{}
+
template<class Type>
inline Foam::SquareMatrix<Type>::SquareMatrix(Istream& is)
@@ -72,6 +83,7 @@ inline Foam::SquareMatrix<Type>::SquareMatrix(Istream& is)
Matrix<SquareMatrix<Type>, Type>(is)
{}
+
template<class Type>
inline Foam::autoPtr<Foam::SquareMatrix<Type>>
Foam::SquareMatrix<Type>::clone() const
diff --git a/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrix.C b/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrix.C
index 28ad1d2745d..dcb31ef2413 100644
--- a/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrix.C
+++ b/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrix.C
@@ -33,15 +33,17 @@ Foam::SymmetricSquareMatrix<Type> Foam::invDecomposed
const SymmetricSquareMatrix<Type>& matrix
)
{
- SymmetricSquareMatrix<Type> inv(matrix.m(), matrix.m(), 0.0);
+ const label n = matrix.n();
- for (label i = 0; i < matrix.m(); ++i)
+ SymmetricSquareMatrix<Type> inv(n, Zero);
+
+ for (label i = 0; i < n; ++i)
{
inv(i, i) = 1.0/matrix(i, i);
for (label j = 0; j < i; ++j)
{
- scalar sum = 0.0;
+ Type sum = Zero;
for (label k = j; k < i; k++)
{
@@ -52,7 +54,20 @@ Foam::SymmetricSquareMatrix<Type> Foam::invDecomposed
}
}
- return inv.T()*inv;
+ SymmetricSquareMatrix<Type> result(n, Zero);
+
+ for (label k = 0; k < n; k++)
+ {
+ for (label i = 0; i <= k; i++)
+ {
+ for (label j = 0; j <= k; j++)
+ {
+ result(i, j) += inv(k, i)*inv(k, j);
+ }
+ }
+ }
+
+ return result;
}
@@ -63,7 +78,6 @@ Foam::SymmetricSquareMatrix<Type> Foam::inv
)
{
SymmetricSquareMatrix<Type> matrixTmp(matrix);
-
LUDecompose(matrixTmp);
return invDecomposed(matrixTmp);
@@ -71,9 +85,9 @@ Foam::SymmetricSquareMatrix<Type> Foam::inv
template<class Type>
-Foam::scalar Foam::detDecomposed(const SymmetricSquareMatrix<Type>& matrix)
+Type Foam::detDecomposed(const SymmetricSquareMatrix<Type>& matrix)
{
- scalar diagProduct = 1.0;
+ Type diagProduct = pTraits<Type>::one;
for (label i = 0; i < matrix.m(); ++i)
{
@@ -85,10 +99,9 @@ Foam::scalar Foam::detDecomposed(const SymmetricSquareMatrix<Type>& matrix)
template<class Type>
-Foam::scalar Foam::det(const SymmetricSquareMatrix<Type>& matrix)
+Type Foam::det(const SymmetricSquareMatrix<Type>& matrix)
{
- SymmetricSquareMatrix<Type> matrixTmp = matrix;
-
+ SymmetricSquareMatrix<Type> matrixTmp(matrix);
LUDecompose(matrixTmp);
return detDecomposed(matrixTmp);
diff --git a/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrix.H b/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrix.H
index 193ee1602f2..543d16b8482 100644
--- a/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrix.H
+++ b/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrix.H
@@ -38,6 +38,7 @@ SourceFiles
#define SymmetricSquareMatrix_H
#include "SquareMatrix.H"
+#include "Identity.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -64,12 +65,15 @@ public:
//- Construct given number of rows/columns.
inline SymmetricSquareMatrix(const label n);
- //- Construct with given number of rows/columns
- inline SymmetricSquareMatrix(const label m, const label n);
+ //- Construct given number of rows/columns, initializing to zero
+ inline SymmetricSquareMatrix(const label n, const zero);
+
+ //- Construct given number of rows/columns,
+ inline SymmetricSquareMatrix(const label n, const Identity<Type>);
//- Construct with given number of rows/columns
- // and value for all elements.
- inline SymmetricSquareMatrix(const label m, const label n, const Type&);
+ // initializing all elements to the given value
+ inline SymmetricSquareMatrix(const label n, const Type&);
//- Construct from Istream.
inline SymmetricSquareMatrix(Istream&);
@@ -91,11 +95,11 @@ SymmetricSquareMatrix<Type> inv(const SymmetricSquareMatrix<Type>&);
//- Return the LU decomposed SymmetricSquareMatrix det
template<class Type>
-scalar detDecomposed(const SymmetricSquareMatrix<Type>&);
+Type detDecomposed(const SymmetricSquareMatrix<Type>&);
//- Return the SymmetricSquareMatrix det
template<class Type>
-scalar det(const SymmetricSquareMatrix<Type>&);
+Type det(const SymmetricSquareMatrix<Type>&);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrixI.H b/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrixI.H
index a363026aac5..821daf199bd 100644
--- a/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrixI.H
+++ b/src/OpenFOAM/matrices/SymmetricSquareMatrix/SymmetricSquareMatrixI.H
@@ -42,17 +42,26 @@ inline Foam::SymmetricSquareMatrix<Type>::SymmetricSquareMatrix(const label n)
template<class Type>
inline Foam::SymmetricSquareMatrix<Type>::SymmetricSquareMatrix
(
- const label m,
- const label n
+ const label n,
+ const zero
)
:
- Matrix<SymmetricSquareMatrix<Type>, Type>(m, n)
+ Matrix<SymmetricSquareMatrix<Type>, Type>(n, n, Zero)
+{}
+
+
+template<class Type>
+inline Foam::SymmetricSquareMatrix<Type>::SymmetricSquareMatrix
+(
+ const label n,
+ const Identity<Type>
+)
+:
+ Matrix<SymmetricSquareMatrix<Type>, Type>(n, n, Zero)
{
- if (m != n)
+ for (label i = 0; i < n; i++)
{
- FatalErrorInFunction
- << "m != n for constructing a symmetric square matrix"
- << exit(FatalError);
+ this->operator()(i, i) = pTraits<Type>::one;
}
}
@@ -60,20 +69,12 @@ inline Foam::SymmetricSquareMatrix<Type>::SymmetricSquareMatrix
template<class Type>
inline Foam::SymmetricSquareMatrix<Type>::SymmetricSquareMatrix
(
- const label m,
const label n,
const Type& t
)
:
- Matrix<SymmetricSquareMatrix<Type>, Type>(m, n, t)
-{
- if (m != n)
- {
- FatalErrorInFunction
- << "m != n for constructing a symmetric square matrix"
- << exit(FatalError);
- }
-}
+ Matrix<SymmetricSquareMatrix<Type>, Type>(n, n, t)
+{}
template<class Type>
diff --git a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixTemplates.C b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixTemplates.C
index 739251c3587..7ac02d189f3 100644
--- a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixTemplates.C
+++ b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixTemplates.C
@@ -35,7 +35,7 @@ Foam::tmp<Foam::Field<Type>> Foam::lduMatrix::H(const Field<Type>& psi) const
{
tmp<Field<Type>> tHpsi
(
- new Field<Type>(lduAddr().size(), pTraits<Type>::zero)
+ new Field<Type>(lduAddr().size(), Zero)
);
if (lowerPtr_ || upperPtr_)
diff --git a/src/OpenFOAM/matrices/scalarMatrices/scalarMatricesTemplates.C b/src/OpenFOAM/matrices/scalarMatrices/scalarMatricesTemplates.C
index c6acc739e27..52413a38dda 100644
--- a/src/OpenFOAM/matrices/scalarMatrices/scalarMatricesTemplates.C
+++ b/src/OpenFOAM/matrices/scalarMatrices/scalarMatricesTemplates.C
@@ -87,7 +87,7 @@ void Foam::solve
// Back-substitution
for (label j=m-1; j>=0; j--)
{
- Type ntempvec = pTraits<Type>::zero;
+ Type ntempvec = Zero;
for (label k=j+1; k<m; k++)
{
diff --git a/src/OpenFOAM/matrices/simpleMatrix/simpleMatrix.C b/src/OpenFOAM/matrices/simpleMatrix/simpleMatrix.C
index 960ba0d563a..0f1c6995b98 100644
--- a/src/OpenFOAM/matrices/simpleMatrix/simpleMatrix.C
+++ b/src/OpenFOAM/matrices/simpleMatrix/simpleMatrix.C
@@ -43,7 +43,7 @@ Foam::simpleMatrix<Type>::simpleMatrix
const Type& sourceVal
)
:
- scalarSquareMatrix(mSize, mSize, coeffVal),
+ scalarSquareMatrix(mSize, coeffVal),
source_(mSize, sourceVal)
{}
diff --git a/src/thermophysicalModels/radiation/radiationModels/viewFactor/viewFactor.C b/src/thermophysicalModels/radiation/radiationModels/viewFactor/viewFactor.C
index 2c8a7e2ba10..d70eccfa9c3 100644
--- a/src/thermophysicalModels/radiation/radiationModels/viewFactor/viewFactor.C
+++ b/src/thermophysicalModels/radiation/radiationModels/viewFactor/viewFactor.C
@@ -173,7 +173,7 @@ void Foam::radiation::viewFactor::initialise()
{
Fmatrix_.reset
(
- new scalarSquareMatrix(totalNCoarseFaces_, totalNCoarseFaces_, 0.0)
+ new scalarSquareMatrix(totalNCoarseFaces_, 0.0)
);
if (debug)
@@ -224,12 +224,7 @@ void Foam::radiation::viewFactor::initialise()
{
CLU_.reset
(
- new scalarSquareMatrix
- (
- totalNCoarseFaces_,
- totalNCoarseFaces_,
- 0.0
- )
+ new scalarSquareMatrix(totalNCoarseFaces_, 0.0)
);
pivotIndices_.setSize(CLU_().m());
@@ -529,7 +524,7 @@ void Foam::radiation::viewFactor::calculate()
// Variable emissivity
if (!constEmissivity_)
{
- scalarSquareMatrix C(totalNCoarseFaces_, totalNCoarseFaces_, 0.0);
+ scalarSquareMatrix C(totalNCoarseFaces_, 0.0);
for (label i=0; i<totalNCoarseFaces_; i++)
{
--
GitLab