OilfieldToolsNet/IFEM
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

changed: do not resize the newton matrix up front

Browse Source 
this way we can use its non-emptyness to signal that
it has already been assembled.
This commit is contained in:
Arne Morten Kvarving 2022-02-16 08:46:05 +01:00
parent bbc3333380
commit a3555e8942
3 changed files with 25 additions and 20 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

21
src/ASM/BlockElmMats.C
View File

@ -73,20 +73,14 @@ bool BlockElmMats::redimOffDiag (size_t blkIndex, char symmetric)
}
bool BlockElmMats::redimNewtonMat ()
bool BlockElmMats::finalize ()
{
// Calculate the total number of element dofs
size_t neldof = 0;
neldof = 0;
size_t nDiagB = blockInfo.size();
for (size_t i = 0; i < nDiagB; i++)
neldof += blockInfo[i].ncmp*basisInfo[blockInfo[i].basis-1].nen;
// Set the element Newton matrix dimension
if (!A.empty() && A.front().empty())
A.front().resize(neldof,neldof);
if (!b.empty() && b.front().empty())
b.front().resize(neldof);
// Calculate the offset of each block sub-matrix
size_t idof = 1;
for (size_t j = 0; j < basisInfo.size(); j++)
@ -103,13 +97,18 @@ bool BlockElmMats::redimNewtonMat ()
idof += ndof;
}
return A.empty() ? true : idof-1 == A.front().rows();
return idof-1 == neldof;
}
const Matrix& BlockElmMats::getNewtonMatrix () const
{
if (!A.front().empty())
return A.front();
Matrix& N = const_cast<Matrix&>(A.front());
// Set the element Newton matrix dimension
N.resize(neldof,neldof);
size_t ib, jb, kb;
size_t nDiagB = blockInfo.size();
@ -174,7 +173,11 @@ const Matrix& BlockElmMats::getNewtonMatrix () const
const Vector& BlockElmMats::getRHSVector () const
{
if (!b.front().empty())
return b.front();
Vector& R = const_cast<Vector&>(b.front());
R.resize(neldof);
for (size_t ib = 0; ib < blockInfo.size(); ib++)
{

8
src/ASM/BlockElmMats.h
View File

@ -32,7 +32,8 @@ public:
//! \brief The constructor initializes the block information arrays.
//! \param[in] nBlk Number of matrix blocks (in each direction, row & column)
//! \param[in] nBas Number of bases (> 1 for mixed problems)
BlockElmMats(size_t nBlk, size_t nBas) : blockInfo(nBlk), basisInfo(nBas) {}
BlockElmMats(size_t nBlk, size_t nBas) :
neldof(0), blockInfo(nBlk), basisInfo(nBas) {}
//! \brief Empty destructor.
virtual ~BlockElmMats() {}
@ -52,10 +53,10 @@ public:
//! \details This method must not be invoked until the \a redim method has
//! been invoked for all the diagonal blocks.
bool redimOffDiag(size_t blkIndex, char symmetric = 1);
//! \brief Sets the dimension of the Newton matrix.
//! \brief Calculates the dimension of the Newton matrix.
//! \details This method must not be invoked until the \a redim method has
//! been invoked for all the diagonal blocks.
bool redimNewtonMat();
bool finalize();
//! \brief Returns the element-level Newton matrix.
virtual const Matrix& getNewtonMatrix() const;
@ -84,6 +85,7 @@ private:
Basis() : nen(0), ncmp(0) {}
};
size_t neldof; //!< Size of newton matrix
std::vector<Block> blockInfo; //!< Block information
std::vector<Basis> basisInfo; //!< Basis information
std::vector<char> symmFlag; //!< Symmetry flag for off-diagonal blocks

16
src/ASM/Test/TestBlockElmMats.C
View File

@ -22,7 +22,7 @@ TEST(TestBlockElmMats, 1Basis2BlocksDiag)
mats.resize(3, 3);
ASSERT_TRUE(mats.redim(1, 2, 1));
ASSERT_TRUE(mats.redim(2, 2, 1));
mats.redimNewtonMat();
mats.finalize();
mats.A[1].fill(1);
mats.A[2].fill(2);
@ -50,7 +50,7 @@ TEST(TestBlockElmMats, 1Basis2BlocksSymmetric)
ASSERT_TRUE(mats.redim(1, 2, 1));
ASSERT_TRUE(mats.redim(2, 2, 1));
ASSERT_TRUE(mats.redimOffDiag(3, 1));
mats.redimNewtonMat();
mats.finalize();
mats.A[1].fill(1);
mats.A[2].fill(2);
@ -82,7 +82,7 @@ TEST(TestBlockElmMats, 1Basis2BlocksSkewSymmetric)
ASSERT_TRUE(mats.redim(1, 2, 1));
ASSERT_TRUE(mats.redim(2, 2, 1));
ASSERT_TRUE(mats.redimOffDiag(3, -1));
mats.redimNewtonMat();
mats.finalize();
mats.A[1].fill(1);
mats.A[2].fill(2);
@ -115,7 +115,7 @@ TEST(TestBlockElmMats, 1Basis2BlocksFull)
ASSERT_TRUE(mats.redim(2, 2, 1));
ASSERT_TRUE(mats.redimOffDiag(3, 0));
ASSERT_TRUE(mats.redimOffDiag(4, 0));
mats.redimNewtonMat();
mats.finalize();
mats.A[1].fill(1);
mats.A[2].fill(2);
@ -146,7 +146,7 @@ TEST(TestBlockElmMats, 2Basis2BlocksDiag)
mats.resize(3, 3);
ASSERT_TRUE(mats.redim(1, 2, 2, 1));
ASSERT_TRUE(mats.redim(2, 2, 1, 2));
mats.redimNewtonMat();
mats.finalize();
mats.A[1].fill(1);
mats.A[2].fill(2);
@ -173,7 +173,7 @@ TEST(TestBlockElmMats, 2Basis2BlocksSymmetric)
ASSERT_TRUE(mats.redim(1, 2, 2, 1));
ASSERT_TRUE(mats.redim(2, 2, 1, 2));
ASSERT_TRUE(mats.redimOffDiag(3, 1));
mats.redimNewtonMat();
mats.finalize();
mats.A[1].fill(1);
mats.A[2].fill(2);
@ -201,7 +201,7 @@ TEST(TestBlockElmMats, 2Basis2BlocksSkewSymmetric)
ASSERT_TRUE(mats.redim(1, 2, 2, 1));
ASSERT_TRUE(mats.redim(2, 2, 1, 2));
ASSERT_TRUE(mats.redimOffDiag(3, -1));
mats.redimNewtonMat();
mats.finalize();
mats.A[1].fill(1);
mats.A[2].fill(2);
@ -230,7 +230,7 @@ TEST(TestBlockElmMats, 2Basis2BlocksFull)
ASSERT_TRUE(mats.redim(2, 2, 1, 2));
ASSERT_TRUE(mats.redimOffDiag(3, 0));
ASSERT_TRUE(mats.redimOffDiag(4, 0));
mats.redimNewtonMat();
mats.finalize();
mats.A[1].fill(1);
mats.A[2].fill(2);