Minor extendedLdu matrices consistency rewrite compared to ldu matrices

This commit is contained in:
Vuko Vukcevic 2015-09-10 09:33:02 +02:00 committed by Hrvoje Jasak
parent 03347200a9
commit 22a2e59c3b
3 changed files with 71 additions and 166 deletions

View file

@ -45,146 +45,73 @@ void Foam::extendedBlockLduMatrix<Type>::mapOffDiagCoeffs
// Get reference to faceMap in extended addressing
const unallocLabelList& faceMap = extLduAddr_.faceMap();
// Avoid assuming it's upper if the matrix is symmetric
if (blockLdum.thereIsUpper())
// Matrix is considered symmetric if the upper is allocated and lower
// is not allocated. Allocating extended upper only.
extendedUpperPtr_ = new TypeCoeffField
(
extLduAddr_.extendedUpperAddr().size()
);
TypeCoeffField& extUpper = *extendedUpperPtr_;
// Get upper coeffs from underlying lduMatrix
const TypeCoeffField& upper = blockLdum.upper();
if (upper.activeType() == blockCoeffBase::SCALAR)
{
// Allocate extended upper only
extendedUpperPtr_ = new TypeCoeffField
(
extLduAddr_.extendedUpperAddr().size()
);
TypeCoeffField& extUpper = *extendedUpperPtr_;
// Helper type definition
typedef typename CoeffField<Type>::scalarTypeField activeType;
// Get upper coeffs from underlying lduMatrix
const TypeCoeffField& upper = blockLdum.upper();
// Get references to fields
const activeType& activeUpper = upper.asScalar();
activeType& activeExtUpper = extUpper.asScalar();
if (upper.activeType() == blockCoeffBase::SCALAR)
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (upper, faceI)
{
// Helper type definition
typedef typename CoeffField<Type>::scalarTypeField activeType;
// Get references to fields
const activeType& activeUpper = upper.asScalar();
activeType& activeExtUpper = extUpper.asScalar();
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (upper, faceI)
{
activeExtUpper[faceMap[faceI]] = activeUpper[faceI];
}
activeExtUpper[faceMap[faceI]] = activeUpper[faceI];
}
else if (upper.activeType() == blockCoeffBase::LINEAR)
}
else if (upper.activeType() == blockCoeffBase::LINEAR)
{
// Helper type definition
typedef typename CoeffField<Type>::linearTypeField activeType;
// Get references to fields
const activeType& activeUpper = upper.asLinear();
activeType& activeExtUpper = extUpper.asLinear();
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (upper, faceI)
{
// Helper type definition
typedef typename CoeffField<Type>::linearTypeField activeType;
// Get references to fields
const activeType& activeUpper = upper.asLinear();
activeType& activeExtUpper = extUpper.asLinear();
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (upper, faceI)
{
activeExtUpper[faceMap[faceI]] = activeUpper[faceI];
}
activeExtUpper[faceMap[faceI]] = activeUpper[faceI];
}
else if (upper.activeType() == blockCoeffBase::SQUARE)
{
// Helper type definition
typedef typename CoeffField<Type>::squareTypeField activeType;
}
else if (upper.activeType() == blockCoeffBase::SQUARE)
{
// Helper type definition
typedef typename CoeffField<Type>::squareTypeField activeType;
// Get references to fields
const activeType& activeUpper = upper.asSquare();
activeType& activeExtUpper = extUpper.asSquare();
// Get references to fields
const activeType& activeUpper = upper.asSquare();
activeType& activeExtUpper = extUpper.asSquare();
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (upper, faceI)
{
activeExtUpper[faceMap[faceI]] = activeUpper[faceI];
}
}
else
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (upper, faceI)
{
FatalErrorIn
(
"extendedBlockLduMatrix(lduMatrix&, label, polyMesh&)"
) << "Problem between ordinary block matrix and extended"
<< " block matrix upper coeffs type morphing."
<< abort(FatalError);
activeExtUpper[faceMap[faceI]] = activeUpper[faceI];
}
}
else
{
// Allocate extended lower only
extendedLowerPtr_ = new TypeCoeffField
FatalErrorIn
(
extLduAddr_.extendedLowerAddr().size()
);
TypeCoeffField& extLower = *extendedLowerPtr_;
// Get lower coeffs from underlying lduMatrix
const TypeCoeffField& lower = blockLdum.lower();
if (lower.activeType() == blockCoeffBase::SCALAR)
{
// Helper type definition
typedef typename CoeffField<Type>::scalarTypeField activeType;
// Get references to fields
const activeType& activeLower = lower.asScalar();
activeType& activeExtLower = extLower.asScalar();
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (lower, faceI)
{
activeExtLower[faceMap[faceI]] = activeLower[faceI];
}
}
else if (lower.activeType() == blockCoeffBase::LINEAR)
{
// Helper type definition
typedef typename CoeffField<Type>::linearTypeField activeType;
// Get references to fields
const activeType& activeLower = lower.asLinear();
activeType& activeExtLower = extLower.asLinear();
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (lower, faceI)
{
activeExtLower[faceMap[faceI]] = activeLower[faceI];
}
}
else if (lower.activeType() == blockCoeffBase::SQUARE)
{
// Helper type definition
typedef typename CoeffField<Type>::squareTypeField activeType;
// Get references to fields
const activeType& activeLower = lower.asSquare();
activeType& activeExtLower = extLower.asSquare();
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (lower, faceI)
{
activeExtLower[faceMap[faceI]] = activeLower[faceI];
}
}
else
{
FatalErrorIn
(
"extendedBlockLduMatrix(lduMatrix&, label, polyMesh&)"
) << "Problem between ordinary block matrix and extended"
<< " block matrix lower coeffs type morphing."
<< abort(FatalError);
}
"extendedBlockLduMatrix(lduMatrix&, label, polyMesh&)"
) << "Problem between ordinary block matrix and extended"
<< " block matrix upper coeffs type morphing."
<< abort(FatalError);
}
}
else

View file

@ -74,46 +74,23 @@ Foam::extendedLduMatrix::extendedLduMatrix
// Get reference to faceMap in extended addressing
const unallocLabelList& faceMap = extLduAddr_.faceMap();
// Avoid assuming it's upper if the matrix is symmetric
if (ldum.hasUpper())
// Matrix is considered symmetric if the upper is allocated and lower
// is not allocated. Allocating extended upper only.
extendedUpperPtr_ = new scalarField
(
extLduAddr_.extendedUpperAddr().size(),
0.0
);
scalarField& extUpper = *extendedUpperPtr_;
// Get upper coeffs from underlying lduMatrix
const scalarField& upper = ldum.upper();
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (upper, faceI)
{
// Allocate extended upper only
extendedUpperPtr_ = new scalarField
(
extLduAddr_.extendedUpperAddr().size(),
0.0
);
scalarField& extUpper = *extendedUpperPtr_;
// Get upper coeffs from underlying lduMatrix
const scalarField& upper = ldum.upper();
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (upper, faceI)
{
extUpper[faceMap[faceI]] = upper[faceI];
}
}
else
{
// Allocate extended lower only
extendedLowerPtr_ = new scalarField
(
extLduAddr_.extendedLowerAddr().size(),
0.0
);
scalarField& extLower = *extendedLowerPtr_;
// Get lower coeffs from underlying lduMatrix
const scalarField& lower = ldum.lower();
// Copy non-zero coeffs from basic lduMatrix into corresponding
// positions
forAll (lower, faceI)
{
extLower[faceMap[faceI]] = lower[faceI];
}
extUpper[faceMap[faceI]] = upper[faceI];
}
}
else

View file

@ -49,8 +49,9 @@ namespace Foam
addasymMatrixConstructorToTable<ILUCp>
addILUCpPreconditionerAsymMatrixConstructorToTable_;
// Add to symmetric constructor table as well until we implement Choleskyp
// preconditioner. VV, 27/Jun/2015.
// Add to symmetric constructor table as well. Cholesky with fill in would
// yield the same sparseness pattern as the original matrix, hence it is not
// implemented. VV, 10/Sep/2015.
lduPreconditioner::
addsymMatrixConstructorToTable<ILUCp>
addILUCpPreconditionerSymMatrixConstructorToTable_;