From cad4b9295e6ab60fa5557de1c655f617099e382a Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Wed, 7 Jun 2017 17:48:56 +0100 Subject: [PATCH] Serial Block GGI SAMG interface --- .../GGIBlockAMGInterfaceField.C | 6 +- .../blockVectorNSAMGInterfaceFields.C | 8 +- .../GGIBlockSAMGInterfaceField.C | 255 ++++++++++++++++++ .../GGIBlockSAMGInterfaceField.H | 224 +++++++++++++++ .../GGIBlockSAMGInterfaceFields.C | 42 +++ .../GGIBlockSAMGInterfaceFields.H | 63 +++++ .../ggiSAMGInterface/ggiSAMGInterface.C | 108 ++++---- 7 files changed, 652 insertions(+), 54 deletions(-) create mode 100644 src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceField.C create mode 100644 src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceField.H create mode 100644 src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceFields.C create mode 100644 src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceFields.H diff --git a/src/foam/matrices/blockLduMatrix/BlockAMG/BlockAMGInterfaceFields/GGIBlockAMGInterfaceField/GGIBlockAMGInterfaceField.C b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockAMGInterfaceFields/GGIBlockAMGInterfaceField/GGIBlockAMGInterfaceField.C index b6353bdca..cedf2639d 100644 --- a/src/foam/matrices/blockLduMatrix/BlockAMG/BlockAMGInterfaceFields/GGIBlockAMGInterfaceField/GGIBlockAMGInterfaceField.C +++ b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockAMGInterfaceFields/GGIBlockAMGInterfaceField/GGIBlockAMGInterfaceField.C @@ -72,7 +72,7 @@ void Foam::GGIBlockAMGInterfaceField::agglomerateBlockType const scalarField& restrictWeights = ggiInterface_.restrictWeights(); // Restrict coefficients - forAll(restrictAddressing, ffi) + forAll (restrictAddressing, ffi) { zoneCoarseCoeffs[restrictAddressing[ffi]] += restrictWeights[ffi]*zoneFineCoeffs[fineAddressing[ffi]]; @@ -237,14 +237,14 @@ void Foam::GGIBlockAMGInterfaceField::updateInterfaceMatrix if (switchToLhs) { - forAll(faceCells, elemI) + forAll (faceCells, elemI) { result[faceCells[elemI]] += pnf[elemI]; } } else { - forAll(faceCells, elemI) + forAll (faceCells, elemI) { result[faceCells[elemI]] -= pnf[elemI]; } diff --git a/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/BlockSAMGInterfaceField/blockVectorNSAMGInterfaceFields.C b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/BlockSAMGInterfaceField/blockVectorNSAMGInterfaceFields.C index a0e031e02..b10691652 100644 --- a/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/BlockSAMGInterfaceField/blockVectorNSAMGInterfaceFields.C +++ b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/BlockSAMGInterfaceField/blockVectorNSAMGInterfaceFields.C @@ -34,6 +34,7 @@ Author #include "BlockSAMGInterfaceField.H" #include "ProcessorBlockSAMGInterfaceField.H" +#include "GGIBlockSAMGInterfaceField.H" #include "VectorNFieldTypes.H" #include "ExpandTensorNField.H" @@ -46,12 +47,15 @@ namespace Foam #define makeTemplateTypeNameAndDebug(type, Type, args...) \ \ -typedef BlockSAMGInterfaceField block##Type##SAMGInterfaceField; \ -defineNamedTemplateTypeNameAndDebug(block##Type##SAMGInterfaceField, 0); \ +typedef BlockSAMGInterfaceField block##Type##SAMGInterfaceField; \ +defineNamedTemplateTypeNameAndDebug(block##Type##SAMGInterfaceField, 0); \ defineTemplateRunTimeSelectionTable(block##Type##SAMGInterfaceField, lduInterface); \ \ typedef ProcessorBlockSAMGInterfaceField block##Type##ProcessorSAMGInterfaceField; \ makeBlockSAMGInterfaceField(block##Type##SAMGInterfaceField, block##Type##ProcessorSAMGInterfaceField); \ + \ +typedef GGIBlockSAMGInterfaceField block##Type##GGISAMGInterfaceField; \ +makeBlockSAMGInterfaceField(block##Type##SAMGInterfaceField, block##Type##GGISAMGInterfaceField); \ forAllVectorNTypes(makeTemplateTypeNameAndDebug); diff --git a/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceField.C b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceField.C new file mode 100644 index 000000000..00ffaec75 --- /dev/null +++ b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceField.C @@ -0,0 +1,255 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | foam-extend: Open Source CFD + \\ / O peration | Version: 4.0 + \\ / A nd | Web: http://www.foam-extend.org + \\/ M anipulation | For copyright notice see file Copyright +------------------------------------------------------------------------------- +License + This file is part of foam-extend. + + foam-extend is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation, either version 3 of the License, or (at your + option) any later version. + + foam-extend is distributed in the hope that it will be useful, but + WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + General Public License for more details. + + You should have received a copy of the GNU General Public License + along with foam-extend. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "GGIBlockSAMGInterfaceField.H" +#include "ggiLduInterfaceField.H" +#include "addToRunTimeSelectionTable.H" +#include "blockLduMatrices.H" + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +template +template +void Foam::GGIBlockSAMGInterfaceField::selectBlockType +( + Field& coarseCoeffs, + const Foam::Field& fineCoeffs +) const +{ + // Note: reconsider better parallel communication here. + // Currently expanding to full zone size + // HJ, 16/Mar/2016 + + // Get fine interface + const ggiLduInterface& fineGgiInterface = ggiInterface_.fineGgiInterface(); + + // Reassemble fine coefficients to full fine zone size + // No need to initialise to zero, as only local coefficients + // are used. HJ, 9/Jun/2016 + Field zoneFineCoeffs(fineGgiInterface.zoneSize()); + + const labelList& fineZa = fineGgiInterface.zoneAddressing(); + + forAll (fineZa, i) + { + zoneFineCoeffs[fineZa[i]] = fineCoeffs[i]; + } + + // Reduce zone data is not required: all coefficients are local + // HJ, 9/Jun/2016 + + Field zoneCoarseCoeffs + ( + ggiInterface_.zoneSize(), + pTraits::zero + ); + + // Get addressing from the fine interface + const labelField& fineAddressing = ggiInterface_.fineAddressing(); + const labelField& restrictAddressing = ggiInterface_.restrictAddressing(); + const scalarField& restrictWeights = ggiInterface_.restrictWeights(); + + // Restrict coefficients + forAll (restrictAddressing, ffi) + { + zoneCoarseCoeffs[restrictAddressing[ffi]] += + restrictWeights[ffi]*zoneFineCoeffs[fineAddressing[ffi]]; + } + + // Filter zone coefficients to local field + const labelList& za = ggiInterface_.zoneAddressing(); + + forAll (za, i) + { + coarseCoeffs[i] = zoneCoarseCoeffs[za[i]]; + } +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template +Foam::GGIBlockSAMGInterfaceField::GGIBlockSAMGInterfaceField +( + const SAMGInterface& SAMGCp, + const BlockLduInterfaceField& fineInterfaceField +) +: + BlockSAMGInterfaceField(SAMGCp, fineInterfaceField), + ggiInterface_(refCast(SAMGCp)), + doTransform_(false), + fieldTransferBuffer_() +{ + // If the interface based on a patch this must be taken care specially of + if (isA >(fineInterfaceField)) + { + const GGIBlockLduInterfaceField& p = + refCast > + ( + fineInterfaceField + ); + + doTransform_ = p.doTransform(); + } + else if (isA(fineInterfaceField)) + { + const ggiLduInterfaceField& p = + refCast(fineInterfaceField); + + doTransform_ = p.doTransform(); + } + else + { + FatalErrorIn("GGIBlockSAMGInterfaceField Constructor") + << "fineInterface must be of ggi type and either" << endl + << " GGIBlockLduInterfaceField or " << endl + << " ggiFvPatchField " << endl + << abort(FatalError); + } +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +template +Foam::GGIBlockSAMGInterfaceField::~GGIBlockSAMGInterfaceField() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template +Foam::tmp > +Foam::GGIBlockSAMGInterfaceField::selectBlockCoeffs +( + const Foam::CoeffField& fineCoeffs +) const +{ + tmp > tcoarseCoeffs(new CoeffField(size())); + CoeffField& coarseCoeffs = tcoarseCoeffs(); + + typedef CoeffField TypeCoeffField; + + typedef typename TypeCoeffField::linearTypeField linearTypeField; + typedef typename TypeCoeffField::squareTypeField squareTypeField; + + // Added weights to account for non-integral matching + if (fineCoeffs.activeType() == blockCoeffBase::SQUARE) + { + squareTypeField& activeCoarseCoeffs = coarseCoeffs.asSquare(); + const squareTypeField& activeFineCoeffs = fineCoeffs.asSquare(); + + this->selectBlockType(activeCoarseCoeffs, activeFineCoeffs); + } + else if (fineCoeffs.activeType() == blockCoeffBase::LINEAR) + { + linearTypeField& activeCoarseCoeffs = coarseCoeffs.asLinear(); + const linearTypeField& activeFineCoeffs = fineCoeffs.asLinear(); + + this->selectBlockType(activeCoarseCoeffs, activeFineCoeffs); + } + else + { + FatalErrorIn + ( + "Foam::tmp >\n" + "Foam::GGIBlockSAMGInterfaceField::selectBlockCoeffs\n" + "(\n" + " const Foam::CoeffField& fineCoeffs\n" + ") const" + ) << "Scalar type agglomeration currently not handled" + << abort(FatalError); + } + + return tcoarseCoeffs; +} + + +template +void Foam::GGIBlockSAMGInterfaceField::initInterfaceMatrixUpdate +( + const Field& psiInternal, + Field&, + const BlockLduMatrix&, + const CoeffField&, + const Pstream::commsTypes, + const bool switchToLhs +) const +{ + // This must have a reduce in it. HJ, 15/May/2009 + Field pif = ggiInterface_.interfaceInternalField(psiInternal); + + fieldTransferBuffer_ = ggiInterface_.fastReduce(pif); +} + + +template +void Foam::GGIBlockSAMGInterfaceField::updateInterfaceMatrix +( + const Field& psiInternal, + Field& result, + const BlockLduMatrix& matrix, + const CoeffField& coeffs, + const Pstream::commsTypes commsType, + const bool switchToLhs +) const +{ + // Get interface from shadow + const GGIBlockSAMGInterfaceField& shadowInterface = + refCast > + ( + matrix.interfaces()[ggiInterface_.shadowIndex()] + ); + + Field pnf = shadowInterface.fieldTransferBuffer(); + + // Complex (VectorN) transformation happens here. + // HJ, 17/Feb/2016 +// transformCoupleField(pnf, cmpt); + + // Multiply neighbour field with coeffs and re-use pnf for result + // of multiplication + multiply(pnf, coeffs, pnf); + + const unallocLabelList& faceCells = ggiInterface_.faceCells(); + + if (switchToLhs) + { + forAll (faceCells, elemI) + { + result[faceCells[elemI]] += pnf[elemI]; + } + } + else + { + forAll (faceCells, elemI) + { + result[faceCells[elemI]] -= pnf[elemI]; + } + } +} + + +// ************************************************************************* // diff --git a/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceField.H b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceField.H new file mode 100644 index 000000000..35df44b2f --- /dev/null +++ b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceField.H @@ -0,0 +1,224 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | foam-extend: Open Source CFD + \\ / O peration | Version: 4.0 + \\ / A nd | Web: http://www.foam-extend.org + \\/ M anipulation | For copyright notice see file Copyright +------------------------------------------------------------------------------- +License + This file is part of foam-extend. + + foam-extend is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation, either version 3 of the License, or (at your + option) any later version. + + foam-extend is distributed in the hope that it will be useful, but + WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + General Public License for more details. + + You should have received a copy of the GNU General Public License + along with foam-extend. If not, see . + +Class + Foam::GGIBlockSAMGInterfaceField + +Description + AMG selected GGI interface field. + +Author + Hrvoje Jasak, Wikki Ltd. All rights reserved + +SourceFiles + GGIBlockSAMGInterfaceField.C + +\*---------------------------------------------------------------------------*/ + +#ifndef GGIBlockSAMGInterfaceField_H +#define GGIBlockSAMGInterfaceField_H + +#include "BlockSAMGInterfaceField.H" +#include "ggiSAMGInterface.H" +#include "GGIBlockLduInterfaceField.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class GGIBlockSAMGInterfaceField Declaration +\*---------------------------------------------------------------------------*/ + +template +class GGIBlockSAMGInterfaceField +: + public BlockSAMGInterfaceField, + public GGIBlockLduInterfaceField +{ + // Private data + + //- Local reference cast into the ggi interface + const ggiSAMGInterface& ggiInterface_; + + //- Is the transform required + bool doTransform_; + + //- Field transfer buffer + mutable Field fieldTransferBuffer_; + + + // Private Member Functions + + //- Disallow default bitwise copy construct + GGIBlockSAMGInterfaceField(const GGIBlockSAMGInterfaceField&); + + //- Disallow default bitwise assignment + void operator=(const GGIBlockSAMGInterfaceField&); + + + //- Select type + template + void selectBlockType + ( + Field& coarseCoeffs, + const Foam::Field& fineCoeffs + ) const; + + +public: + + //- Runtime type information + TypeName("ggi"); + + + // Constructors + + //- Construct from SAMG interface and fine level interface field + GGIBlockSAMGInterfaceField + ( + const SAMGInterface& SAMGCp, + const BlockLduInterfaceField& fineInterfaceField + ); + + + //- Destructor + virtual ~GGIBlockSAMGInterfaceField(); + + + // Member Functions + + // Access + + //- Return size + label size() const + { + return ggiInterface_.size(); + } + + + // Agglomeration + + // Klas Jareteg: 2013-02-06. Moved by HJ, 16/Mar/2016 + //- Agglomerating for the CoeffField fine-level coefficients + virtual tmp > selectBlockCoeffs + ( + const CoeffField& fineCoeffs + ) const; + + + // Block coupled interface matrix update + + //- Transform given patch component field + virtual void transformCoupleField + ( + scalarField& f, + const direction cmpt + ) const + { + GGIBlockLduInterfaceField::transformCoupleField + ( + f, + cmpt + ); + } + + //- Transform neighbour field + virtual void transformCoupleField + ( + Field& f + ) const + { + GGIBlockLduInterfaceField::transformCoupleField(f); + } + + //- Initialise neighbour matrix update + virtual void initInterfaceMatrixUpdate + ( + const Field&, + Field&, + const BlockLduMatrix&, + const CoeffField&, + const Pstream::commsTypes commsType, + const bool switchToLhs + ) const; + + //- Update result field based on interface functionality + virtual void updateInterfaceMatrix + ( + const Field&, + Field&, + const BlockLduMatrix&, + const CoeffField&, + const Pstream::commsTypes commsType, + const bool switchToLhs + ) const; + + + //- Ggi interface functions + + //- Does the interface field perform the transfromation + virtual bool doTransform() const + { + return doTransform_; + } + + //- Return face transformation tensor + virtual const tensorField& forwardT() const + { + return ggiInterface_.forwardT(); + } + + //- Return neighbour-cell transformation tensor + virtual const tensorField& reverseT() const + { + return ggiInterface_.reverseT(); + } + + + // Transfer buffer access + + //- Return contents of the field transfer buffer + const Field& fieldTransferBuffer() const + { + return fieldTransferBuffer_; + } +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository +# include "GGIBlockSAMGInterfaceField.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceFields.C b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceFields.C new file mode 100644 index 000000000..63dc99482 --- /dev/null +++ b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceFields.C @@ -0,0 +1,42 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | foam-extend: Open Source CFD + \\ / O peration | Version: 4.0 + \\ / A nd | Web: http://www.foam-extend.org + \\/ M anipulation | For copyright notice see file Copyright +------------------------------------------------------------------------------- +License + This file is part of foam-extend. + + foam-extend is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation, either version 3 of the License, or (at your + option) any later version. + + foam-extend is distributed in the hope that it will be useful, but + WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + General Public License for more details. + + You should have received a copy of the GNU General Public License + along with foam-extend. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "GGIBlockSAMGInterfaceFields.H" +#include "blockSAMGInterfaceFields.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +makeBlockSAMGInterfaceFields(GGIBlockSAMGInterfaceField); + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceFields.H b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceFields.H new file mode 100644 index 000000000..25de7be74 --- /dev/null +++ b/src/foam/matrices/blockLduMatrix/BlockAMG/BlockSAMGInterfaceFields/GGIBlockSAMGInterfaceField/GGIBlockSAMGInterfaceFields.H @@ -0,0 +1,63 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | foam-extend: Open Source CFD + \\ / O peration | Version: 4.0 + \\ / A nd | Web: http://www.foam-extend.org + \\/ M anipulation | For copyright notice see file Copyright +------------------------------------------------------------------------------- +License + This file is part of foam-extend. + + foam-extend is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation, either version 3 of the License, or (at your + option) any later version. + + foam-extend is distributed in the hope that it will be useful, but + WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + General Public License for more details. + + You should have received a copy of the GNU General Public License + along with foam-extend. If not, see . + +Class + GGIBlockSAMGInterfaceField + +Description + Typedefs for block coefficient ggi SAMG interface fields + +Author + Hrvoje Jasak, Wikki Ltd. All rights reserved. + +SourceFiles + GGIBlockSAMGInterfaceFields.C + +\*---------------------------------------------------------------------------*/ + +#ifndef GGIBlockSAMGInterfaceFields_H +#define GGIBlockSAMGInterfaceFields_H + +#include "GGIBlockSAMGInterfaceField.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +typedef GGIBlockSAMGInterfaceField GGIBlockSAMGInterfaceFieldScalar; +typedef GGIBlockSAMGInterfaceField GGIBlockSAMGInterfaceFieldVector; +typedef GGIBlockSAMGInterfaceField GGIBlockSAMGInterfaceFieldTensor; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/foam/matrices/lduMatrix/solvers/AMG/interfaces/SAMGInterfaces/ggiSAMGInterface/ggiSAMGInterface.C b/src/foam/matrices/lduMatrix/solvers/AMG/interfaces/SAMGInterfaces/ggiSAMGInterface/ggiSAMGInterface.C index 5934c9ccd..67775b4a7 100644 --- a/src/foam/matrices/lduMatrix/solvers/AMG/interfaces/SAMGInterfaces/ggiSAMGInterface/ggiSAMGInterface.C +++ b/src/foam/matrices/lduMatrix/solvers/AMG/interfaces/SAMGInterfaces/ggiSAMGInterface/ggiSAMGInterface.C @@ -282,6 +282,20 @@ Foam::ggiSAMGInterface::ggiSAMGInterface // Note: local addressing contains only local faces const labelList& fineZa = fineGgiInterface_.zoneAddressing(); + // Expand master prolongation to zone + crMatrix masterExpandProlongation(interfaceProlongation); + + //HJ, HERE: expand master prolongation to zone + // Note: master is now the size of local zone + if (!fineGgiInterface_.localParallel()) + { + // Not line this: expand without communication HJ, HERE + // fineGgiInterface_.expandCrMatrixToZone + // ( + // masterExpandProlongation + // ); + } + // Create crMatrix for neighbour faces. Note: expandCrMatrixToZone will // expand the matrix to zone size, including communications. // Faces which are not used locally will be marked by empty rows @@ -298,10 +312,6 @@ Foam::ggiSAMGInterface::ggiSAMGInterface ); } - // - const crAddressing& prolongationCr = interfaceProlongation.crAddr(); - const crAddressing& nbrExpandCr = nbrExpandProlongation.crAddr(); - // Create addressing for neighbour processors. Note: expandAddrToZone will // expand the addressing to zone size. HJ, 13/Jun/2016 labelField neighbourExpandProc @@ -369,6 +379,23 @@ Foam::ggiSAMGInterface::ggiSAMGInterface // Count the number of agglomeration pairs label nAgglomPairs = 0; + // Switching prolongation matrices + const crMatrix* masterP = NULL; + const crMatrix* neighbourP = NULL; + + if (fineGgiInterface_.master()) + { + // Grab prolongation matrix + masterP = &masterExpandProlongation; + neighbourP = &nbrExpandProlongation; + } + else + { + // Grab prolongation matrix + masterP = &nbrExpandProlongation; + neighbourP = &masterExpandProlongation; + } + // On the fine level, addressing is made in a labelListList if (fineGgiInterface_.fineLevel()) { @@ -455,32 +482,31 @@ Foam::ggiSAMGInterface::ggiSAMGInterface // restriction - with included weights from GGI for ( - label indexR = prolongationCr.rowStart()[curSide]; - indexR < prolongationCr.rowStart()[curSide + 1]; + label indexR = masterP->crAddr().rowStart()[curSide]; + indexR < masterP->crAddr().rowStart()[curSide + 1]; indexR++ ) { // Grab weight from restriction - scalar rWeight = interfaceProlongation.coeffs()[indexR]; + scalar rWeight = masterP->coeffs()[indexR]; - // HJ, replace nbrInterfaceProlongation with nbrExpandProlongation for ( - label indexP = nbrExpandCr.rowStart()[nbrSide]; - indexP < nbrExpandCr.rowStart()[nbrSide + 1]; + label indexP = neighbourP->crAddr().rowStart()[nbrSide]; + indexP < neighbourP->crAddr().rowStart()[nbrSide + 1]; indexP++ ) { // Grab weight from prolongation - scalar pWeight = nbrInterfaceProlongation.coeffs()[indexP]; + scalar pWeight = neighbourP->coeffs()[indexP]; // Code in the current master and slave - used for // identifying the face curMaster = - prolongationCr.column()[indexR] + masterP->crAddr().column()[indexR] + procOffset*curMasterProc; curSlave = - nbrExpandCr.column()[indexP] + neighbourP->crAddr().column()[indexP] + procOffset*curSlaveProc; if (neighboursTable.found(curMaster)) @@ -520,7 +546,10 @@ Foam::ggiSAMGInterface::ggiSAMGInterface nbrFound = true; curFaceFaces[curNbrI].append(ffI); curFaceFaceNbrs[curNbrI].append(nbrI); - curFaceWeights[curNbrI].append(curNW*pWeight*rWeight); + curFaceWeights[curNbrI].append + ( + curNW*pWeight*rWeight + ); // New agglomeration pair found in already // existing pair @@ -609,16 +638,6 @@ Foam::ggiSAMGInterface::ggiSAMGInterface } } // end for all current neighbours } // end for all fine faces - if (fineGgiInterface_.master()) - { - Info<< "MASTER: "; - } - else - { - Info<< "SLAVE: "; - } - Info<< "Done fine level, 1. Created " << nAgglomPairs << " pairs and " - << nCoarseFaces << " faces" << endl; } //------------------------------------------------------------------------------ // FINE LEVEL - no GGI weights! @@ -685,30 +704,29 @@ Foam::ggiSAMGInterface::ggiSAMGInterface // restriction - with included weights from GGI for ( - label indexR = prolongationCr.rowStart()[curSide]; - indexR < prolongationCr.rowStart()[curSide + 1]; + label indexR = masterP->crAddr().rowStart()[curSide]; + indexR < masterP->crAddr().rowStart()[curSide + 1]; indexR++ ) { // Grab weight from restriction - scalar rWeight = interfaceProlongation.coeffs()[indexR]; + scalar rWeight = masterP->coeffs()[indexR]; - // HJ, replace nbrInterfaceProlongation with nbrExpandProlongation for ( - label indexP = nbrExpandCr.rowStart()[nbrSide]; - indexP < nbrExpandCr.rowStart()[nbrSide + 1]; + label indexP = neighbourP->crAddr().rowStart()[nbrSide]; + indexP < neighbourP->crAddr().rowStart()[nbrSide + 1]; indexP++ ) { // Grab weight from prolongation - scalar pWeight = nbrInterfaceProlongation.coeffs()[indexP]; + scalar pWeight = neighbourP->coeffs()[indexP]; // Code in the current master and slave - used for // identifying the face - curMaster = prolongationCr.column()[indexR] + curMaster = masterP->crAddr().column()[indexR] + procOffset*curMasterProc; - curSlave = nbrExpandCr.column()[indexP] + curSlave = neighbourP->crAddr().column()[indexP] + procOffset*curSlaveProc; if (neighboursTable.found(curMaster)) @@ -835,8 +853,6 @@ Foam::ggiSAMGInterface::ggiSAMGInterface } } } // end for all fine faces - Info<< "Done coarse level, 1. Created " << nAgglomPairs << " pairs and " - << nCoarseFaces << " faces" << endl; } // end of else in fine level (coarse level) // Since only local faces are analysed, lists can now be resized @@ -894,7 +910,7 @@ Foam::ggiSAMGInterface::ggiSAMGInterface // Sort makes sure the order is identical on both sides. // HJ, 20/Feb/2009 and 6/Jun/2016 sort(contents); - Info<< "START MATRIX ASSEMBLY" << endl; + // Note: Restriction is done on master side only because this is where // the local zone is created. HJ, 1/Aug/2016 if (master()) @@ -983,8 +999,6 @@ Foam::ggiSAMGInterface::ggiSAMGInterface nProcFaces++; } } - Info<< "MASTER ASSEMBLY: Created " << nAgglomPairs << " pairs and " - << nProcFaces << " master faces" << endl; // No need to resize arrays only local faces are used // HJ, 1/Aug/2016 @@ -1137,11 +1151,7 @@ Foam::ggiSAMGInterface::ggiSAMGInterface nProcFaces++; } } - Info<< "SLAVE ASSEMBLY: Created " << nAgglomPairs << " pairs and " - << nProcFaces << " slave faces" << endl; - } - Info<< "END MATRIX ASSEMBLY" << endl; } @@ -1337,7 +1347,6 @@ void Foam::ggiSAMGInterface::expandAddrToZone(labelField& lf) const void Foam::ggiSAMGInterface::expandCrMatrixToZone(crMatrix&) const { - notImplemented("expandCrMatrixToZone"); // Code missing: collapse crMatrices into a zone crMatrix if (!localParallel()) { @@ -1353,9 +1362,8 @@ void Foam::ggiSAMGInterface::initProlongationTransfer const crMatrix& filteredP ) const { - // Send prolongation matrix, using IOstream operators - //OPstream toNbr(Pstream::blocking, neighbProcNo()); - //toNbr<< filteredP; + // crMatrix transfer is local without global reduction + crMatrixTransferBuffer_ = filteredP; } @@ -1366,11 +1374,13 @@ Foam::ggiSAMGInterface::prolongationTransfer const crMatrix& filteredP ) const { - //IPstream fromNbr(Pstream::blocking, neighbProcNo()); - - autoPtr tnbrP(new crMatrix(5,5, labelList(5,1))); + autoPtr tnbrP + ( + new crMatrix(shadowInterface().crMatrixTransferBuffer()) + ); return tnbrP; } + // ************************************************************************* //