exaFOAM/foam-extend4.1-coherent-io
Archived
1
1
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity

Changed crMatrix and crAddressing member naming

Browse source
This commit is contained in:
Hrvoje Jasak 2017-05-10 13:19:15 +01:00
parent 9e9168c389
commit bae31764dc
6 changed files with 249 additions and 231 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

77
src/foam/matrices/crMatrix/crAddressing.C
View file

@ -36,28 +36,30 @@ Author
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::crAddressing::setRowCount(const labelList& count)
void Foam::crAddressing::setRowSizes(const labelList& rowSizes)
{
if (count.size() != nRows_)
if (rowSizes.size() != nRows_)
{
FatalErrorIn("void crAddressing::setRowCount(const labelList& count)")
<< "Incorrect size of count: nRows =" << nRows_
<< " count = " << count.size()
FatalErrorIn
(
"void crAddressing::setRowSizes(const labelList& rowSizes)"
) << "Incorrect size of rowSizes: nRows =" << nRows_
<< " rowSizes = " << rowSizes.size()
<< abort(FatalError);
}
// Accumulate count
row_[0] = 0;
// Accumulate rowSizes
rowStart_[0] = 0;
forAll (count, i)
forAll (rowSizes, i)
{
row_[i + 1] = row_[i] + count[i];
rowStart_[i + 1] = rowStart_[i] + rowSizes[i];
}
// Resize and clear column array
col_.setSize(row_[nRows_]);
col_ = 0;
column_.setSize(rowStart_[nRows_]);
column_ = 0;
}
@ -68,16 +70,16 @@ Foam::crAddressing::crAddressing
(
const label nRows,
const label nCols,
const labelList& count
const labelList& rowSizes
)
:
refCount(),
nRows_(nRows),
nCols_(nCols),
row_(nRows + 1),
col_(0)
rowStart_(nRows + 1),
column_(0)
{
setRowCount(count);
setRowSizes(rowSizes);
}
@ -93,8 +95,8 @@ Foam::crAddressing::crAddressing
refCount(),
nRows_(nRows),
nCols_(nCols),
row_(row),
col_(col)
rowStart_(row),
column_(col)
{}
@ -104,8 +106,8 @@ Foam::crAddressing::crAddressing(const crAddressing& a)
refCount(),
nRows_(a.nRows_),
nCols_(a.nCols_),
row_(a.row_),
col_(a.col_)
rowStart_(a.rowStart_),
column_(a.column_)
{}
@ -115,8 +117,8 @@ Foam::crAddressing::crAddressing(Istream& is)
refCount(),
nRows_(readLabel(is)),
nCols_(readLabel(is)),
row_(is),
col_(is)
rowStart_(is),
column_(is)
{}
@ -124,31 +126,34 @@ Foam::crAddressing::crAddressing(Istream& is)
Foam::tmp<Foam::crAddressing> Foam::crAddressing::T() const
{
const labelList& myRow = row();
const labelList& myCol = col();
const labelList& myRow = rowStart();
const labelList& myCol = column();
labelList trCount(nCols(), 0);
labelList trRowSizes(nCols(), 0);
// Count number of entries in a row
for (label i = 0; i < nRows(); i++)
{
for (label ip = myRow[i]; ip < myRow[i + 1]; ip++)
{
trCount[myCol[ip]]++;
trRowSizes[myCol[ip]]++;
}
}
// Create transpose addressing
tmp<crAddressing> ttranspose(new crAddressing(nCols(), nRows(), trCount));
tmp<crAddressing> ttranspose
(
new crAddressing(nCols(), nRows(), trRowSizes)
);
crAddressing& transpose = ttranspose();
// Set coefficients
const labelList& trRow = transpose.row();
labelList& trCol = transpose.col();
const labelList& trRow = transpose.rowStart();
labelList& trCol = transpose.column();
trCol = 0;
// Reset count to use as counter
trCount = 0;
// Reset rowSizes to use as counter
trRowSizes = 0;
label j;
@ -158,9 +163,9 @@ Foam::tmp<Foam::crAddressing> Foam::crAddressing::T() const
{
j = myCol[ip];
trCol[trRow[j] + trCount[j]] = i;
trCol[trRow[j] + trRowSizes[j]] = i;
trCount[j]++;
trRowSizes[j]++;
}
}
@ -184,8 +189,8 @@ void Foam::crAddressing::operator=(const crAddressing& rhs)
nRows_ = rhs.nRows_;
nCols_ = rhs.nCols_;
row_ = rhs.row_;
col_ = rhs.col_;
rowStart_ = rhs.rowStart_;
column_ = rhs.column_;
}
@ -194,8 +199,8 @@ void Foam::crAddressing::operator=(const crAddressing& rhs)
Foam::Ostream& Foam::operator<<(Ostream& os, const crAddressing& a)
{
os << a.nRows_ << tab << a.nCols_ << nl
<< a.row_ << nl
<< a.col_ << endl;
<< a.rowStart_ << nl
<< a.column_ << endl;
return os;
}

31
src/foam/matrices/crMatrix/crAddressing.H
View file

@ -64,15 +64,15 @@ class crAddressing
//- Row array.
// Provides start index for each row, dimensioned to nRows + 1
labelList row_;
labelList rowStart_;
//- Column array
labelList col_;
//- Column index array
labelList column_;
// Private Member Functions
//- Set row count
void setRowCount(const labelList& count);
//- Set row sizes
void setRowSizes(const labelList& count);
public:
@ -84,7 +84,7 @@ public:
(
const label nRows,
const label nCols,
const labelList& count
const labelList& nEntries
);
//- Construct from components
@ -92,8 +92,8 @@ public:
(
const label nRows,
const label nCols,
const labelList& row,
const labelList& col
const labelList& rowStart,
const labelList& column
);
//- Construct as copy
@ -125,28 +125,28 @@ public:
//- Return number of coefficients
label nEntries() const
{
return col_.size();
return column_.size();
}
//- Return row array
const labelList& row() const
const labelList& rowStart() const
{
return row_;
return rowStart_;
}
//- Return column array
const labelList& col() const
const labelList& column() const
{
return col_;
return column_;
}
// Edit
//- Return column array
labelList& col()
labelList& column()
{
return col_;
return column_;
}
@ -160,6 +160,7 @@ public:
void operator=(const crAddressing&);
// IOstream Operators
friend Ostream& operator<<(Ostream&, const crAddressing&);

12
src/foam/matrices/crMatrix/crMatrix.C
View file

@ -109,8 +109,8 @@ Foam::tmp<Foam::crMatrix> Foam::crMatrix::T() const
{
// My addressing
const label myNRows = crAddr().nRows();
const labelList& myRow = crAddr().row();
const labelList& myCol = crAddr().col();
const labelList& myRow = crAddr().rowStart();
const labelList& myCol = crAddr().column();
const scalarField& myCoeffs = coeffs();
// Create transpose
@ -123,8 +123,8 @@ Foam::tmp<Foam::crMatrix> Foam::crMatrix::T() const
tCoeffs = 0;
// Transpose addressing
const labelList& tRow = transpose.crAddr().row();
const labelList& tCol = transpose.crAddr().col();
const labelList& tRow = transpose.crAddr().rowStart();
const labelList& tCol = transpose.crAddr().column();
// Set transpose coefficients
@ -153,8 +153,8 @@ Foam::tmp<Foam::crMatrix> Foam::crMatrix::T() const
// Calculate b += A*x
void Foam::crMatrix::dotPlus(scalarField& b, const scalarField& x) const
{
const labelList& row = crAddr_.row();
const labelList& col = crAddr_.col();
const labelList& row = crAddr_.rowStart();
const labelList& col = crAddr_.column();
forAll (b, i)
{

8
src/foam/matrices/crMatrix/crMatrix.H
View file

@ -83,8 +83,8 @@ public:
(
const label nRows,
const label nCols,
const labelList& row,
const labelList& col
const labelList& rowStart,
const labelList& column
);
//- Construct as copy
@ -123,9 +123,9 @@ public:
}
//- Return column array to be set
labelList& col()
labelList& column()
{
return crAddr_.col();
return crAddr_.column();
}

339
src/foam/matrices/lduMatrix/solvers/AMG/interfaces/SAMGInterfaces/processorSAMGInterface/processorSAMGInterface.C
View file

@ -45,210 +45,215 @@ namespace Foam
Foam::processorSAMGInterface::processorSAMGInterface
(
const lduPrimitiveMesh& lduMesh,
const crMatrix& prolongation,
const lduInterfacePtrsList& coarseInterfaces,
const lduInterface& fineInterface,
const labelField& localRestrictAddressing,
const labelField& neighbourRestrictAddressing
const labelField& localRowLabel,
const labelField& neighbourRowLabel
)
:
SAMGInterface(lduMesh),
SAMGInterface(lduMesh, prolongation),
fineProcInterface_(refCast<const processorLduInterface>(fineInterface)),
comm_(fineProcInterface_.comm()),
tag_(fineProcInterface_.tag())
{
Pout<< "Creating processor SAMG interface" << endl;
/* HJ, Code missing here
// Make a lookup table of entries for owner/neighbour
HashTable<SLList<label>, label, Hash<label> > neighboursTable
// Analyse the local and neighbour row label:
// local coarse, remote coarse = regular coarse face
// local coarse, remote fine = local expanded face: receive prolonged data
// local fine, remote coarse = neighbour expanded face: send prolonged data
// Algorithm
// Go through local row labels and examine cases
// 1) coarse local and coarse remote:
// create coarse processor face and set faceCells. Set weight to 1
// 2) coarse local and fine remote:
// will receive coarse neighbours and weights from opposite side
// 3) fine local and coarse remote:
// assemble local coarse neighbours and weights from the prolongation
// 4) fine local and fine remote - ignore
//
// On completion of selection:
// send and receive local coarse neighbours for fine local/coarse remote
// sort coarse equations by increasing coarse index from master side
// create faceCells, fineAddressing and fineWeights
// First collect and communicate internal neighbours from the local fine
// side (ie neighbour processor cell is coarse)
HashTable<labelList, label, Hash<label> > neighboursFromLocalFine
(
localRestrictAddressing.size()
Foam::max(128, fineProcInterface_.interfaceSize()/4)
);
// Table of face-sets to be agglomerated
HashTable<SLList<SLList<label> >, label, Hash<label> > faceFaceTable
HashTable<scalarField, label, Hash<label> > weightsFromLocalFine
(
localRestrictAddressing.size()
Foam::max(128, fineProcInterface_.interfaceSize()/4)
);
// Get access to the prolongation addressing and coefficients
const labelList& pRowStart = prolongation.crAddr().rowStart();
const labelList& pColumn = prolongation.crAddr().column();
const scalarField& pCoeffs = prolongation.coeffs();
// Get fine faceCells
const labelList& fineFaceCells = fineInterface.faceCells();
Pout<< "fineFaceCells: " << fineFaceCells << endl;
// Collect local fine to neighbour coarse connections for communication
forAll (localRowLabel, faceI)
{
if
(
localRowLabel[faceI] < 0
&& neighbourRowLabel[faceI] >= 0
)
{
// Found local fine to neighbour coarse interface
// Collect local coarse neighbours and weights from the
// prolongation matrix to send to other processor
const label curStart = pRowStart[fineFaceCells[faceI]];
const label curEnd = pRowStart[fineFaceCells[faceI] + 1];
const label nCoarse = curEnd - curStart;
Pout<< "Eqn: " << fineFaceCells[faceI] << " Span: " << curStart << " " << curEnd << " = " << nCoarse << endl;
labelList nbrs(nCoarse);
scalarField weights(nCoarse);
forAll (nbrs, i)
{
nbrs[i] = pColumn[curStart + i];
weights[i] = pCoeffs[curStart + i];
}
Pout<< "weights: " << weights << endl;
// Insert neighbours under remote coarse index
neighboursFromLocalFine.insert(neighbourRowLabel[faceI], nbrs);
weightsFromLocalFine.insert(neighbourRowLabel[faceI], weights);
}
}
// Receive remote prolongation data
HashTable<labelList, label, Hash<label> > neighboursFromRemoteFine;
HashTable<scalarField, label, Hash<label> > weightsFromRemoteFine;
// Send and receive the addressing from the other side
{
OPstream toNbr(Pstream::blocking, neighbProcNo());
toNbr<< neighboursFromLocalFine << weightsFromLocalFine;
}
{
IPstream fromNbr(Pstream::blocking, neighbProcNo());
neighboursFromRemoteFine =
HashTable<labelList, label, Hash<label> >(fromNbr);
weightsFromRemoteFine =
HashTable<scalarField, label, Hash<label> >(fromNbr);
}
// Assemble connectivity
// Resize arrays to size of fine interface
// Note: it is technically possible to have MORE coarse processor faces
// so a check will be made in the end
faceCells_.setSize(5*fineProcInterface_.interfaceSize());
fineAddressing_.setSize(5*fineProcInterface_.interfaceSize());
fineWeights_.setSize(5*fineProcInterface_.interfaceSize());
// Count coarse faces
label nCoarseFaces = 0;
forAll (localRestrictAddressing, ffi)
// Collect coarse-to-fine connections
forAll (localRowLabel, faceI)
{
label curMaster = -1;
label curSlave = -1;
// Do switching on master/slave indexes based on the owner/neighbour of
// the processor index such that both sides get the same answer.
if (myProcNo() < neighbProcNo())
if
(
localRowLabel[faceI] >= 0
&& neighbourRowLabel[faceI] >= 0
)
{
// Master side
curMaster = localRestrictAddressing[ffi];
curSlave = neighbourRestrictAddressing[ffi];
}
else
{
// Slave side
curMaster = neighbourRestrictAddressing[ffi];
curSlave = localRestrictAddressing[ffi];
}
// Look for the master cell. If it has already got a face,
// add the coefficient to the face. If not, create a new face.
if (neighboursTable.found(curMaster))
{
// Check all current neighbours to see if the current slave already
// exists and if so, add the fine face to the agglomeration.
SLList<label>& curNbrs = neighboursTable.find(curMaster)();
SLList<SLList<label> >& curFaceFaces =
faceFaceTable.find(curMaster)();
bool nbrFound = false;
SLList<label>::iterator nbrsIter = curNbrs.begin();
SLList<SLList<label> >::iterator faceFacesIter =
curFaceFaces.begin();
for
(
;
nbrsIter != curNbrs.end(), faceFacesIter != curFaceFaces.end();
++nbrsIter, ++faceFacesIter
)
{
if (nbrsIter() == curSlave)
{
nbrFound = true;
faceFacesIter().append(ffi);
break;
}
}
if (!nbrFound)
{
curNbrs.append(curSlave);
curFaceFaces.append(SLList<label>(ffi));
// New coarse face created
nCoarseFaces++;
}
}
else
{
// This master has got no neighbours yet. Add a neighbour
// and a coefficient, thus creating a new face
neighboursTable.insert(curMaster, SLList<label>(curSlave));
faceFaceTable.insert
(
curMaster,
SLList<SLList<label> >(SLList<label>(ffi))
);
// New coarse face created
// Found local coarse to neighbour coarse face
Pout<< "face " << faceI << " CC" << endl;
// Create new coarse face
faceCells_[nCoarseFaces] = localRowLabel[faceI];
fineAddressing_[nCoarseFaces] = faceI;
fineWeights_[nCoarseFaces] = 1;
nCoarseFaces++;
}
} // end for all fine faces
faceCells_.setSize(nCoarseFaces, -1);
fineAddressing_.setSize(localRestrictAddressing.size(), -1);
restrictAddressing_.setSize(localRestrictAddressing.size(), -1);
// All weights are equal to 1: integral matching
restrictWeights_.setSize(localRestrictAddressing.size(), 1.0);
labelList contents = neighboursTable.toc();
// Sort makes sure the order is identical on both sides.
// HJ, 20/Feb.2009
sort(contents);
// Reset face counter for re-use
nCoarseFaces = 0;
if (myProcNo() < neighbProcNo())
{
// On master side, the owner addressing is stored in table of contents
forAll (contents, masterI)
else if
(
localRowLabel[faceI] < 0
&& neighbourRowLabel[faceI] >= 0
)
{
SLList<label>& curNbrs = neighboursTable.find(contents[masterI])();
// Found local fine to neighbour coarse face
SLList<SLList<label> >& curFaceFaces =
faceFaceTable.find(contents[masterI])();
// Pick up local prolongation coarse entries and for all
// add a new face with appropriate weight
// Note: faceCells changes due to (internal) local coarse cells
const labelList& curLocalCoarseNbrs =
neighboursFromLocalFine[neighbourRowLabel[faceI]];
SLList<label>::iterator nbrsIter = curNbrs.begin();
SLList<SLList<label> >::iterator faceFacesIter =
curFaceFaces.begin();
for
(
;
nbrsIter != curNbrs.end(), faceFacesIter != curFaceFaces.end();
++nbrsIter, ++faceFacesIter
)
const scalarField& curLocalCoarseWeights =
weightsFromLocalFine[neighbourRowLabel[faceI]];
Pout<< "face " << faceI << " FC, size: " << curLocalCoarseNbrs.size()
<< " W: " << curLocalCoarseWeights << endl;
forAll (curLocalCoarseNbrs, curNbrI)
{
faceCells_[nCoarseFaces] = contents[masterI];
for
(
SLList<label>::iterator facesIter =
faceFacesIter().begin();
facesIter != faceFacesIter().end();
++facesIter
)
{
fineAddressing_[facesIter()] = facesIter();
restrictAddressing_[facesIter()] = nCoarseFaces;
}
// Create new coarse face
faceCells_[nCoarseFaces] = curLocalCoarseNbrs[curNbrI];
fineAddressing_[nCoarseFaces] = faceI;
fineWeights_[nCoarseFaces] = curLocalCoarseWeights[curNbrI];
nCoarseFaces++;
}
}
}
else
{
// On slave side, the owner addressing is stored in linked lists
forAll (contents, masterI)
else if
(
localRowLabel[faceI] >= 0
&& neighbourRowLabel[faceI] < 0
)
{
SLList<label>& curNbrs = neighboursTable.find(contents[masterI])();
// Found local coarse to neighbour fine face
SLList<SLList<label> >& curFaceFaces =
faceFaceTable.find(contents[masterI])();
// Pick up neighbour prolongation coarse entries and for all
// add a new face with appropriate weight
// Note: faceCells remains the same: single local coarse cell
const labelList& curNbrCoarseNbrs =
neighboursFromRemoteFine[localRowLabel[faceI]];
SLList<label>::iterator nbrsIter = curNbrs.begin();
SLList<SLList<label> >::iterator faceFacesIter =
curFaceFaces.begin();
for
(
;
nbrsIter != curNbrs.end(), faceFacesIter != curFaceFaces.end();
++nbrsIter, ++faceFacesIter
)
const scalarField& curNbrCoarseWeights =
weightsFromRemoteFine[localRowLabel[faceI]];
Pout<< "face " << faceI << " CF, size: " << curNbrCoarseNbrs.size()
<< " W: " << curNbrCoarseWeights << endl;
forAll (curNbrCoarseNbrs, curNbrI)
{
faceCells_[nCoarseFaces] = nbrsIter();
for
(
SLList<label>::iterator facesIter = faceFacesIter().begin();
facesIter != faceFacesIter().end();
++facesIter
)
{
fineAddressing_[facesIter()] = facesIter();
restrictAddressing_[facesIter()] = nCoarseFaces;
}
// Create new coarse face
faceCells_[nCoarseFaces] = localRowLabel[faceI];
fineAddressing_[nCoarseFaces] = faceI;
fineWeights_[nCoarseFaces] = curNbrCoarseWeights[curNbrI];
nCoarseFaces++;
}
}
else
{
// Fine to fine. Ignore
}
}
*/
// Check for fixed lists
if (nCoarseFaces > 5*fineProcInterface_.interfaceSize())
{
FatalErrorIn("processorSAMGInterface::processorSAMGInterface(...)")
<< "Coarse SAMG processor siginificantly bigger than fine: "
<< "nCoarseFaces = " << nCoarseFaces
<< " nFineFaces = " << fineProcInterface_.interfaceSize()
<< abort(FatalError);
}
Pout<< "nCoarseFaces: " << nCoarseFaces << endl;
// Resize arrays to final size
faceCells_.setSize(nCoarseFaces);
fineAddressing_.setSize(nCoarseFaces);
fineWeights_.setSize(nCoarseFaces);
}

13
src/foam/matrices/lduMatrix/solvers/AMG/interfaces/SAMGInterfaces/processorSAMGInterface/processorSAMGInterface.H
View file

@ -83,14 +83,15 @@ public:
// Constructors
//- Construct from fine-level interface,
// local and neighbour restrict addressing
// local and neighbour row label
processorSAMGInterface
(
const lduPrimitiveMesh& lduMesh,
const crMatrix& prolongation,
const lduInterfacePtrsList& coarseInterfaces,
const lduInterface& fineInterface,
const labelField& localRestrictAddressing,
const labelField& neighbourRestrictAddressing
const labelField& localRowLabel,
const labelField& neighbourRowLabel
);
@ -102,6 +103,12 @@ public:
// Access
//- Return interface size
virtual label interfaceSize() const
{
return SAMGInterface::size();
}
//- Return true if interface is coupled
virtual bool coupled() const
{