54 std::vector<int> & myRowGidOffsets)
58 std::vector<int> myRowOffsets(per_dof_row_basis.size());
59 for(std::size_t i=0;i<per_dof_row_basis.size();i++) {
60 myRowOffsets[i] = totalSz;
61 totalSz += per_dof_row_basis[i]->size();
65 Teuchos::RCP<Epetra_Map> rowMap = Teuchos::rcp(
new Epetra_Map(-1,totalSz,0,Comm));
67 myRowGidOffsets.resize(per_dof_row_basis.size());
68 for(std::size_t i=0;i<myRowGidOffsets.size();i++)
69 myRowGidOffsets[i] = rowMap->GID(myRowOffsets[i]);
89 const std::vector<int> & myRowGidOffsets,
90 std::vector<int> & myColGidOffsets)
95 for(std::size_t i=0;i<myRowGidOffsets.size();i++)
96 adaptRowGids[i] = myRowGidOffsets[i];
99 Epetra_Import importer(determGraph.ColMap(),determGraph.RowMap());
101 adaptColGids.Import(adaptRowGids,importer,
Insert);
104 myColGidOffsets.resize(adaptColGids.MyLength());
105 for(std::size_t i=0;i<myColGidOffsets.size();i++)
106 myColGidOffsets[i] = adaptColGids[i];
121 int numStochDim = masterBasis->dimension();
124 Epetra_BlockMap blkRowMap(-1,rowMap.NumMyElements(),rowMap.MyGlobalElements(),numStochDim,0,rowMap.Comm());
125 Epetra_BlockMap blkColMap(-1,colMap.NumMyElements(),colMap.MyGlobalElements(),numStochDim,0,colMap.Comm());
128 stochColOrders(blkColMap);
131 for(std::size_t dof=0;dof<per_dof_row_basis.size();dof++) {
132 RCP<const Stokhos::ProductBasis<int,double> > rowBasis = per_dof_row_basis[dof];
133 TEUCHOS_TEST_FOR_EXCEPTION(rowBasis->dimension()!=masterBasis->dimension(),std::invalid_argument,
134 "Stokhos::adapt_utils::buildColBasisFunctions: Row basis must match dimension of master basis!");
136 Teuchos::Array<RCP<const OneDOrthogPolyBasis<int,double> > > onedBasis
137 = rowBasis->getCoordinateBases();
139 TEUCHOS_TEST_FOR_EXCEPTION(onedBasis.size()!=numStochDim,std::logic_error,
140 "Stokhos::adapt_utils::buildColBasisFunctions: Wrong number of dimensions from row basis!");
143 for(
int i=0;i<numStochDim;i++)
144 stochRowOrders[i+dof*numStochDim] = onedBasis[i]->order();
149 stochColOrders.Import(stochRowOrders,importer,
Insert);
151 Teuchos::Array<RCP<const OneDOrthogPolyBasis<int,double> > > oneDMasterBasis
152 = masterBasis->getCoordinateBases();
155 std::vector<int> polyOrder(numStochDim,0);
156 per_dof_col_basis.resize(blkColMap.NumMyElements());
157 for(std::size_t col=0;col<per_dof_col_basis.size();col++) {
158 int colOffset = numStochDim*col;
161 for(
int o=0;o<numStochDim;o++)
162 polyOrder[o] = stochColOrders[colOffset+o];
164 Teuchos::Array<Teuchos::RCP<const Stokhos::OneDOrthogPolyBasis<int,double> > > newBasisArray(numStochDim);
165 for(Teuchos::Ordinal dim=0;dim<oneDMasterBasis.size();dim++) {
166 RCP<const OneDOrthogPolyBasis<int,double> > oneDBasis = oneDMasterBasis[dim];
168 newBasisArray[dim] = oneDBasis->cloneWithOrder(polyOrder[dim]);
171 per_dof_col_basis[col] = Teuchos::rcp(
192 std::vector<int> & myRowGidOffsets,std::vector<int> & myColGidOffsets,
196 TEUCHOS_TEST_FOR_EXCEPTION(
int(per_dof_row_basis.size())!=determGraph.NumMyRows(),std::logic_error,
197 "Stokhos::adapt_utils::buildAdaptedGraph: per_dof_row_basis.size()!=determGraph.NumMyRows()");
199 myRowGidOffsets.clear();
200 myColGidOffsets.clear();
202 std::vector<Teuchos::RCP<const Stokhos::ProductBasis<int,double> > > per_dof_col_basis;
203 Teuchos::RCP<Epetra_Map> rowMap;
214 Teuchos::RCP<const Stokhos::Sparse3Tensor<int,double> > Cijk;
216 Cijk = masterBasis->computeTripleProductTensor();
218 Cijk = masterBasis->computeLinearTripleProductTensor();
221 int maxNNZ = determGraph.MaxNumNonzeros();
222 std::vector<int> determGraphCols(maxNNZ);
223 std::vector<int> graphCols;
224 for(
int lRID=0;lRID<determGraph.NumMyRows();lRID++) {
225 int gRID = determGraph.GRID(lRID);
227 int rowOffsetIndex = myRowGidOffsets[lRID];
230 determGraph.ExtractGlobalRowCopy(gRID,maxNNZ,numIndices,&determGraphCols[0]);
232 for(
int i=0;i<numIndices;i++) {
233 int gCID = determGraphCols[i];
234 int lCID = determGraph.LCID(gCID);
235 int colOffsetIndex = myColGidOffsets[lCID];
238 *per_dof_col_basis[lCID],
239 onlyUseLinear,kExpOrder);
240 for(std::size_t basisRow=0;basisRow<interactGraph.rowCount();basisRow++) {
241 const std::vector<std::size_t> & basisCols = interactGraph.activeIndices(basisRow);
242 graphCols.resize(basisCols.size());
243 for(std::size_t
g=0;
g<basisCols.size();
g++)
244 graphCols[
g] = basisCols[
g] + colOffsetIndex;
246 graph->InsertGlobalIndices(basisRow+rowOffsetIndex,graphCols.size(),&graphCols[0]);
251 graph->FillComplete();
void buildColBasisFunctions(const Epetra_CrsGraph &determGraph, const Teuchos::RCP< const Stokhos::ProductBasis< int, double > > &masterBasis, const std::vector< Teuchos::RCP< const Stokhos::ProductBasis< int, double > > > &per_dof_row_basis, std::vector< Teuchos::RCP< const Stokhos::ProductBasis< int, double > > > &per_dof_col_basis)
Teuchos::RCP< Epetra_CrsGraph > buildAdaptedGraph(const Epetra_CrsGraph &determGraph, const Teuchos::RCP< const Stokhos::ProductBasis< int, double > > &masterBasis, const std::vector< Teuchos::RCP< const Stokhos::ProductBasis< int, double > > > &per_dof_row_basis, bool onlyUseLinear=false, int kExpOrder=-1)