9#include <fei_ParameterSet.hpp>
10#include "fei_MatrixReducer.hpp"
11#include "fei_EqnComm.hpp"
12#include "fei_Matrix_core.hpp"
13#include "fei_sstream.hpp"
14#include "fei_fstream.hpp"
15#include <fei_CommUtils.hpp>
23 globalAssembleCalled_(false),
24 changedSinceMark_(false)
27 target->getMatrixGraph()->getRowSpace();
28 MPI_Comm comm = vspace->getCommunicator();
29 int numLocal = reducer_->getLocalReducedEqns().size();
31 fei::Matrix_core* target_core =
32 dynamic_cast<fei::Matrix_core*
>(target_.get());
33 if (target_core == NULL) {
34 throw std::runtime_error(
"fei::MatrixReducer ERROR, target matrix not dynamic_cast-able to fei::Matrix_core.");
37 target_core->setEqnComm(eqnComm);
40MatrixReducer::~MatrixReducer()
47 return(target_->parameters(paramset));
53 target_->setMatrixGraph(matrixGraph);
57MatrixReducer::getGlobalNumRows()
const
59 return(target_->getMatrixGraph()->getRowSpace()->getGlobalNumIndices());
63MatrixReducer::getLocalNumRows()
const
65 return(target_->getMatrixGraph()->getRowSpace()->getNumIndices_Owned());
68int MatrixReducer::putScalar(
double scalar)
69{
return(target_->putScalar(scalar)); }
72MatrixReducer::getRowLength(
int row,
int& length)
const
74 if (reducer_->isSlaveEqn(row)) {
75 FEI_OSTRINGSTREAM osstr;
76 osstr <<
"fei::MatrixReducer::getRowLength ERROR, row="<<row<<
" is a slave eqn. You can't get a slave row from the reduced matrix.";
77 throw std::runtime_error(osstr.str());
80 int reducedrow = reducer_->translateToReducedEqn(row);
81 return(target_->getRowLength(reducedrow, length));
85MatrixReducer::copyOutRow(
int row,
int len,
double* coefs,
int* indices)
const
87 if (reducer_->isSlaveEqn(row)) {
88 FEI_OSTRINGSTREAM osstr;
89 osstr <<
"fei::MatrixReducer::copyOutRow ERROR, requested row ("<<row
90 <<
") is a slave eqn. You can't get a slave row from the reduced matrix.";
91 throw std::runtime_error(osstr.str());
94 int reducedrow = reducer_->translateToReducedEqn(row);
95 int err = target_->copyOutRow(reducedrow, len, coefs, indices);
96 for(
int i=0; i<len; ++i) {
97 indices[i] = reducer_->translateFromReducedEqn(indices[i]);
103MatrixReducer::sumIn(
int numRows,
const int* rows,
104 int numCols,
const int* cols,
105 const double*
const* values,
108 int err = reducer_->addMatrixValues(numRows, rows, numCols, cols,
109 values,
true, *target_, format);
114MatrixReducer::copyIn(
int numRows,
const int* rows,
115 int numCols,
const int* cols,
116 const double*
const* values,
119 int err = reducer_->addMatrixValues(numRows, rows, numCols, cols,
120 values,
false, *target_, format);
125MatrixReducer::sumInFieldData(
int fieldID,
129 const double*
const* data,
133 target_->getMatrixGraph()->getRowSpace();
135 target_->getMatrixGraph()->getColSpace();
137 unsigned fieldSize = rowSpace->getFieldSize(fieldID);
138 std::vector<int> indices(fieldSize*2);
139 int* rowIndices = &indices[0];
140 int* colIndices = rowIndices+fieldSize;
142 rowSpace->getGlobalIndices(1, &rowID, idType, fieldID, rowIndices);
143 colSpace->getGlobalIndices(1, &colID, idType, fieldID, colIndices);
145 if (format != FEI_DENSE_ROW) {
146 throw std::runtime_error(
"MatrixReducer: bad format");
149 int err = reducer_->addMatrixValues(fieldSize, rowIndices,
150 fieldSize, colIndices,
151 data,
true, *target_, format);
156MatrixReducer::sumInFieldData(
int fieldID,
164 target_->getMatrixGraph()->getRowSpace();
166 unsigned fieldSize = rowSpace->getFieldSize(fieldID);
168 std::vector<const double*> data_2d(fieldSize);
169 for(
unsigned i=0; i<fieldSize; ++i) {
170 data_2d[i] = &data[i*fieldSize];
173 return(sumInFieldData(fieldID, idType, rowID, colID, &data_2d[0], format));
177MatrixReducer::sumIn(
int blockID,
int connectivityID,
178 const double*
const* values,
182 int numRowIndices, numColIndices, dummy;
183 matGraph->getConnectivityNumIndices(blockID, numRowIndices, numColIndices);
185 std::vector<int> indices(numRowIndices+numColIndices);
186 int* rowIndices = &indices[0];
187 int* colIndices = rowIndices+numRowIndices;
189 matGraph->getConnectivityIndices(blockID, connectivityID,
190 numRowIndices, rowIndices, dummy,
191 numColIndices, colIndices, dummy);
193 return(sumIn(numRowIndices, rowIndices, numColIndices, colIndices,
198MatrixReducer::globalAssemble()
200 reducer_->assembleReducedMatrix(*target_);
201 return(target_->globalAssemble());
205{
return(target_->multiply(x, y)); }
207int MatrixReducer::gatherFromOverlap(
bool accumulate)
209 reducer_->assembleReducedMatrix(*target_);
210 target_->setCommSizes();
211 return(target_->gatherFromOverlap(accumulate));
214int MatrixReducer::writeToFile(
const char* filename,
215 bool matrixMarketFormat)
217 static char mmbanner[] =
"%%MatrixMarket matrix coordinate real general";
218 std::vector<int>& localrows = reducer_->getLocalReducedEqns();
219 int localNumRows = localrows.size();
224 for(
int i=0; i<localNumRows; ++i) {
226 CHK_ERR( target_->getRowLength(localrows[i], len) );
230 MPI_Comm comm = getMatrixGraph()->getRowSpace()->getCommunicator();
233 int globalNumRows = 0;
236 int globalNumCols = globalNumRows;
242 FEI_OFSTREAM* outFile = NULL;
244 outFile =
new FEI_OFSTREAM(filename, IOS_OUT);
245 FEI_OFSTREAM& ofs = *outFile;
246 if (matrixMarketFormat) {
247 ofs << mmbanner << FEI_ENDL;
248 ofs <<globalNumRows<<
" " <<globalNumCols<<
" " <<globalNNZ<<FEI_ENDL;
251 ofs <<globalNumRows<<
" " <<globalNumCols<<FEI_ENDL;
254 else outFile =
new FEI_OFSTREAM(filename, IOS_APP);
256 outFile->setf(IOS_SCIENTIFIC, IOS_FLOATFIELD);
257 outFile->precision(13);
258 FEI_OFSTREAM& ofs = *outFile;
261 std::vector<int> work_indices;
262 std::vector<double> work_data1D;
264 for(
int i=0; i<localNumRows; ++i) {
265 int row = localrows[i];
266 CHK_ERR( target_->getRowLength(row, rowLength) );
268 work_indices.resize(rowLength);
269 work_data1D.resize(rowLength);
271 int* indPtr = &work_indices[0];
272 double* coefPtr = &work_data1D[0];
274 CHK_ERR( target_->copyOutRow(row, rowLength, coefPtr, indPtr) );
276 for(
int j=0; j<rowLength; ++j) {
277 if (matrixMarketFormat) {
278 ofs << row+1 <<
" "<<indPtr[j]+1<<
" "<<coefPtr[j]<<FEI_ENDL;
281 ofs << row <<
" "<<indPtr[j]<<
" "<<coefPtr[j]<<FEI_ENDL;
292int MatrixReducer::writeToStream(FEI_OSTREAM& ostrm,
293 bool matrixMarketFormat)
295 return(target_->writeToStream(ostrm, matrixMarketFormat));
298void MatrixReducer::markState()
299{ target_->markState(); }
301bool MatrixReducer::changedSinceMark()
302{
return(target_->changedSinceMark()); }
int localProc(MPI_Comm comm)
int numProcs(MPI_Comm comm)
int GlobalSum(MPI_Comm comm, std::vector< T > &local, std::vector< T > &global)