4 #ifndef PANZER_L2_PROJECTION_IMPL_HPP
5 #define PANZER_L2_PROJECTION_IMPL_HPP
8 #include "Tpetra_CrsGraph.hpp"
9 #include "Tpetra_MultiVector.hpp"
10 #include "Tpetra_CrsMatrix.hpp"
14 #include "Panzer_TpetraLinearObjFactory.hpp"
30 const std::vector<std::string>& elementBlockNames,
33 targetBasisDescriptor_ = targetBasis;
34 integrationDescriptor_ = integrationDescriptor;
36 connManager_ = connManager;
37 elementBlockNames_ = elementBlockNames;
38 worksetContainer_ = worksetContainer;
40 useUserSuppliedBasisValues_ =
false;
43 targetGlobalIndexer_ =
47 for (
const auto& eBlock : elementBlockNames_) {
48 std::vector<shards::CellTopology> topologies;
49 connManager_->getElementBlockTopologies(topologies);
50 std::vector<std::string> ebNames;
51 connManager_->getElementBlockIds(ebNames);
52 const auto search = std::find(ebNames.cbegin(),ebNames.cend(),eBlock);
54 const int index = std::distance(ebNames.cbegin(),search);
55 const auto& cellTopology = topologies[index];
57 auto intrepidBasis = panzer::createIntrepid2Basis<PHX::Device,double,double>(targetBasisDescriptor_.getType(),
58 targetBasisDescriptor_.getOrder(),
62 targetGlobalIndexer_->addField(eBlock,targetBasisDescriptor_.getType(),fieldPattern);
65 targetGlobalIndexer_->buildGlobalUnknowns();
72 useUserSuppliedBasisValues_ =
true;
77 for (
const auto& eblock : elementBlockNames_) {
78 TEUCHOS_ASSERT(basisValues_[eblock]->getBasisDescriptor()==targetBasisDescriptor_);
84 {
return targetGlobalIndexer_;}
88 const std::unordered_map<std::string,double>* elementBlockMultipliers)
90 PANZER_FUNC_TIME_MONITOR_DIFF(
"L2Projection::Build Mass Matrix",TopBuildMassMatrix);
93 if (elementBlockMultipliers !=
nullptr) {
94 TEUCHOS_ASSERT(elementBlockMultipliers->size() == elementBlockNames_.size());
98 std::vector<Teuchos::RCP<const panzer::GlobalIndexer>> indexers;
99 indexers.push_back(targetGlobalIndexer_);
105 ownedMatrix->resumeFill();
106 ownedMatrix->setAllToScalar(0.0);
107 ghostedMatrix->resumeFill();
108 ghostedMatrix->setAllToScalar(0.0);
110 const int fieldIndex = targetGlobalIndexer_->getFieldNum(targetBasisDescriptor_.getType());
112 const bool is_scalar = targetBasisDescriptor_.getType()==
"HGrad" || targetBasisDescriptor_.getType()==
"Const" || targetBasisDescriptor_.getType()==
"HVol";
116 auto M = ghostedMatrix->getLocalMatrixDevice();
117 for (
const auto& block : elementBlockNames_) {
119 double ebMultiplier = 1.0;
120 if (elementBlockMultipliers !=
nullptr)
121 ebMultiplier = elementBlockMultipliers->find(block)->second;
126 int num_cells_owned_ghosted_virtual = 0;
127 int num_cells_owned = 0;
129 if (useUserSuppliedBasisValues_) {
131 auto tmp = connManager_->getElementBlock(block);
137 Kokkos::deep_copy(cell_local_ids_nonconst,cell_local_ids_host);
138 cell_local_ids = cell_local_ids_nonconst;
140 bv_ptr = basisValues_[block].get();
141 num_cells_owned_ghosted_virtual = cell_local_ids.extent(0);
142 num_cells_owned = cell_local_ids.extent(0);
149 const auto worksets = worksetContainer_->getWorksets(wd);
151 if (worksets->size() == 0)
156 const auto& workset = (*worksets)[0];
157 bv_ptr = &(workset.getBasisValues(targetBasisDescriptor_,integrationDescriptor_));
158 num_cells_owned_ghosted_virtual = workset.numOwnedCells()+workset.numGhostCells()+workset.numVirtualCells();
159 num_cells_owned = workset.numOwnedCells();
160 cell_local_ids = workset.getLocalCellIDs();
162 const auto& basisValues = *bv_ptr;
165 const auto weightedBasis = basisValues.getBasisValues(
true).get_static_view();
167 const std::vector<panzer::LocalOrdinal>&
offsets = targetGlobalIndexer_->getGIDFieldOffsets(block,fieldIndex);
168 PHX::View<panzer::LocalOrdinal*> kOffsets(
"MassMatrix: Offsets",offsets.size());
169 auto kOffsets_h = Kokkos::create_mirror_view(kOffsets);
171 for(
const auto& i : offsets)
172 kOffsets_h(i) = offsets[i];
174 Kokkos::deep_copy(kOffsets, kOffsets_h);
177 PHX::View<panzer::LocalOrdinal**> localIds(
"MassMatrix: LocalIds", num_cells_owned_ghosted_virtual,
178 targetGlobalIndexer_->getElementBlockGIDCount(block));
181 const auto cellLocalIdsNoGhost = Kokkos::subview(cell_local_ids,std::make_pair(0,num_cells_owned));
183 targetGlobalIndexer_->getElementLIDs(cellLocalIdsNoGhost,localIds);
185 const int numBasisPoints =
static_cast<int>(weightedBasis.extent(1));
187 Kokkos::parallel_for(num_cells_owned,KOKKOS_LAMBDA (
const int& cell) {
188 double total_mass = 0.0, trace = 0.0;
190 panzer::LocalOrdinal cLIDs[256];
191 const int numIds =
static_cast<int>(localIds.extent(1));
192 for(
int i=0;i<numIds;++i)
193 cLIDs[i] = localIds(cell,i);
195 double vals[256]={0.0};
196 const int numQP =
static_cast<int>(unweightedBasis.extent(2));
198 for (
int row=0; row < numBasisPoints; ++row) {
199 for (
int col=0; col < numIds; ++col) {
200 for (
int qp=0; qp < numQP; ++qp) {
201 auto tmp = unweightedBasis(cell,row,qp) * weightedBasis(cell,col,qp) * ebMultiplier;
209 for (
int row=0; row < numBasisPoints; ++row) {
210 for (
int col=0; col < numBasisPoints; ++col)
213 int offset = kOffsets(row);
214 panzer::LocalOrdinal lid = localIds(cell,offset);
217 for (
int qp=0; qp < numQP; ++qp)
218 vals[col] += unweightedBasis(cell,row,qp) * weightedBasis(cell,col,qp) * ebMultiplier * total_mass / trace;
220 M.sumIntoValues(lid,cLIDs,numIds,vals,
true,
true);
225 Kokkos::parallel_for(num_cells_owned,KOKKOS_LAMBDA (
const int& cell) {
226 panzer::LocalOrdinal cLIDs[256];
227 const int numIds =
static_cast<int>(localIds.extent(1));
228 for(
int i=0;i<numIds;++i)
229 cLIDs[i] = localIds(cell,i);
232 const int numQP =
static_cast<int>(unweightedBasis.extent(2));
234 for (
int row=0; row < numBasisPoints; ++row) {
235 int offset = kOffsets(row);
236 panzer::LocalOrdinal lid = localIds(cell,offset);
238 for (
int col=0; col < numIds; ++col) {
240 for (
int qp=0; qp < numQP; ++qp)
241 vals[col] += unweightedBasis(cell,row,qp) * weightedBasis(cell,col,qp) * ebMultiplier;
243 M.sumIntoValues(lid,cLIDs,numIds,vals,
true,
true);
250 auto M = ghostedMatrix->getLocalMatrixDevice();
251 for (
const auto& block : elementBlockNames_) {
253 double ebMultiplier = 1.0;
254 if (elementBlockMultipliers !=
nullptr)
255 ebMultiplier = elementBlockMultipliers->find(block)->second;
260 int num_cells_owned_ghosted_virtual = 0;
261 int num_cells_owned = 0;
263 if (useUserSuppliedBasisValues_) {
265 auto tmp = connManager_->getElementBlock(block);
271 Kokkos::deep_copy(cell_local_ids_nonconst,cell_local_ids_host);
272 cell_local_ids = cell_local_ids_nonconst;
274 bv_ptr = basisValues_[block].get();
275 num_cells_owned_ghosted_virtual = cell_local_ids.extent(0);
276 num_cells_owned = cell_local_ids.extent(0);
283 const auto worksets = worksetContainer_->getWorksets(wd);
285 if (worksets->size() == 0)
290 const auto& workset = (*worksets)[0];
291 bv_ptr = &(workset.getBasisValues(targetBasisDescriptor_,integrationDescriptor_));
292 num_cells_owned_ghosted_virtual = workset.numOwnedCells()+workset.numGhostCells()+workset.numVirtualCells();
293 num_cells_owned = workset.numOwnedCells();
294 cell_local_ids = workset.getLocalCellIDs();
296 const auto& basisValues = *bv_ptr;
299 const auto weightedBasis = basisValues.getVectorBasisValues(
true).get_static_view();
301 const std::vector<panzer::LocalOrdinal>&
offsets = targetGlobalIndexer_->getGIDFieldOffsets(block,fieldIndex);
302 PHX::View<panzer::LocalOrdinal*> kOffsets(
"MassMatrix: Offsets",offsets.size());
303 auto kOffsets_h = Kokkos::create_mirror_view(kOffsets);
305 for(
const auto& i : offsets)
306 kOffsets_h(i) = offsets[i];
308 Kokkos::deep_copy(kOffsets, kOffsets_h);
311 PHX::View<panzer::LocalOrdinal**> localIds(
"MassMatrix: LocalIds", num_cells_owned_ghosted_virtual,
312 targetGlobalIndexer_->getElementBlockGIDCount(block));
315 const PHX::View<const int*> cellLocalIdsNoGhost = Kokkos::subview(cell_local_ids,std::make_pair(0,num_cells_owned));
317 targetGlobalIndexer_->getElementLIDs(cellLocalIdsNoGhost,localIds);
319 const int numBasisPoints =
static_cast<int>(weightedBasis.extent(1));
320 Kokkos::parallel_for(num_cells_owned,KOKKOS_LAMBDA (
const int& cell) {
322 panzer::LocalOrdinal cLIDs[256];
323 const int numIds =
static_cast<int>(localIds.extent(1));
324 for(
int i=0;i<numIds;++i)
325 cLIDs[i] = localIds(cell,i);
328 const int numQP =
static_cast<int>(unweightedBasis.extent(2));
330 for (
int qp=0; qp < numQP; ++qp) {
331 for (
int row=0; row < numBasisPoints; ++row) {
332 int offset = kOffsets(row);
333 panzer::LocalOrdinal lid = localIds(cell,offset);
335 for (
int col=0; col < numIds; ++col){
337 for(
int dim=0; dim < static_cast<int>(weightedBasis.extent(3)); ++dim)
338 vals[col] += unweightedBasis(cell,row,qp,dim) * weightedBasis(cell,col,qp,dim) * ebMultiplier;
341 M.sumIntoValues(lid,cLIDs,numIds,vals,
true,
true);
349 PANZER_FUNC_TIME_MONITOR_DIFF(
"Exporting of mass matrix",ExportMM);
350 auto map = factory.
getMap(0);
351 ghostedMatrix->fillComplete(map,map);
353 ownedMatrix->doExport(*ghostedMatrix, *exporter, Tpetra::ADD);
354 ownedMatrix->fillComplete();
362 PANZER_FUNC_TIME_MONITOR_DIFF(
"L2Projection<panzer::LocalOrdinal,panzer::GlobalOrdinal>::buildInverseLumpedMassMatrix",buildInvLMM);
364 const auto massMatrix = this->buildMassMatrix(
true);
365 const auto lumpedMassMatrix =
rcp(
new Tpetra::MultiVector<double,panzer::LocalOrdinal,panzer::GlobalOrdinal,panzer::TpetraNodeType>(massMatrix->getDomainMap(),1,
true));
366 const auto tmp =
rcp(
new Tpetra::MultiVector<double,panzer::LocalOrdinal,panzer::GlobalOrdinal,panzer::TpetraNodeType>(massMatrix->getRangeMap(),1,
false));
369 PANZER_FUNC_TIME_MONITOR_DIFF(
"Apply",Apply);
370 massMatrix->apply(*tmp,*lumpedMassMatrix);
373 PANZER_FUNC_TIME_MONITOR_DIFF(
"Reciprocal",Reciprocal);
374 lumpedMassMatrix->reciprocal(*lumpedMassMatrix);
376 return lumpedMassMatrix;
381 const Teuchos::RCP<
const Tpetra::Map<panzer::LocalOrdinal,panzer::GlobalOrdinal,panzer::TpetraNodeType>>& inputOwnedSourceMap,
382 const std::string& sourceFieldName,
384 const int directionIndex)
391 using MapType = Tpetra::Map<panzer::LocalOrdinal,panzer::GlobalOrdinal,panzer::TpetraNodeType>;
392 using GraphType = Tpetra::CrsGraph<panzer::LocalOrdinal,panzer::GlobalOrdinal,panzer::TpetraNodeType>;
393 using ExportType = Tpetra::Export<panzer::LocalOrdinal,panzer::GlobalOrdinal,panzer::TpetraNodeType>;
394 using MatrixType = Tpetra::CrsMatrix<double,panzer::LocalOrdinal,panzer::GlobalOrdinal,panzer::TpetraNodeType>;
402 std::vector<panzer::GlobalOrdinal> indices;
403 targetGlobalIndexer_->getOwnedAndGhostedIndices(indices);
409 std::vector<panzer::GlobalOrdinal> indices;
417 std::vector<size_t> nEntriesPerRow(ghostedTargetMap->getLocalNumElements(),0);
418 std::vector<std::string> elementBlockIds;
419 targetGlobalIndexer_->getElementBlockIds(elementBlockIds);
420 std::vector<std::string>::const_iterator blockItr;
421 for (blockItr=elementBlockIds.begin();blockItr!=elementBlockIds.end();++blockItr) {
422 std::string blockId = *blockItr;
423 const std::vector<panzer::LocalOrdinal> & elements = targetGlobalIndexer_->getElementBlock(blockId);
425 std::vector<panzer::GlobalOrdinal> row_gids;
426 std::vector<panzer::GlobalOrdinal> col_gids;
428 for(std::size_t elmt=0;elmt<elements.size();elmt++) {
429 targetGlobalIndexer_->getElementGIDs(elements[elmt],row_gids);
431 for(std::size_t row=0;row<row_gids.size();row++) {
432 panzer::LocalOrdinal lid =
433 ghostedTargetMap->getLocalElement(row_gids[row]);
434 nEntriesPerRow[lid] += col_gids.size();
440 RCP<GraphType> ghostedGraph =
rcp(
new GraphType(ghostedTargetMap,ghostedSourceMap,nEntriesPerRowView));
442 for (blockItr=elementBlockIds.begin();blockItr!=elementBlockIds.end();++blockItr) {
443 std::string blockId = *blockItr;
444 const std::vector<panzer::LocalOrdinal> & elements = targetGlobalIndexer_->getElementBlock(blockId);
446 std::vector<panzer::GlobalOrdinal> row_gids;
447 std::vector<panzer::GlobalOrdinal> col_gids;
449 for(std::size_t elmt=0;elmt<elements.size();elmt++) {
450 targetGlobalIndexer_->getElementGIDs(elements[elmt],row_gids);
452 for(std::size_t row=0;row<row_gids.size();row++)
453 ghostedGraph->insertGlobalIndices(row_gids[row],col_gids);
459 std::vector<panzer::GlobalOrdinal> indices;
460 targetGlobalIndexer_->getOwnedIndices(indices);
466 std::vector<panzer::GlobalOrdinal> indices;
472 ghostedGraph->fillComplete(ownedSourceMap,ownedTargetMap);
480 ownedGraph->doExport(*ghostedGraph, *exporter, Tpetra::INSERT);
481 ownedGraph->fillComplete(ownedSourceMap,ownedTargetMap);
488 ghostedMatrix->setAllToScalar(0.0);
489 ownedMatrix->setAllToScalar(0.0);
494 for (
const auto& block : elementBlockNames_) {
497 const auto& worksets = worksetContainer_->getWorksets(wd);
498 for (
const auto& workset : *worksets) {
501 const auto& targetBasisValues = workset.getBasisValues(targetBasisDescriptor_,integrationDescriptor_);
502 const auto& targetWeightedBasis = targetBasisValues.getBasisValues(
true).get_static_view();
505 const auto& sourceBasisValues = workset.getBasisValues(sourceBasisDescriptor,integrationDescriptor_);
506 PHX::View<const double***> sourceUnweightedScalarBasis;
507 PHX::View<const double****> sourceUnweightedVectorBasis;
508 bool useRankThreeBasis =
false;
509 if ( (sourceBasisDescriptor.
getType() ==
"HGrad") || (sourceBasisDescriptor.
getType() ==
"Const") || (sourceBasisDescriptor.
getType() ==
"HVol") ) {
510 if (directionIndex == -1) {
511 sourceUnweightedScalarBasis = sourceBasisValues.getBasisValues(
false).get_static_view();
512 useRankThreeBasis =
true;
515 sourceUnweightedVectorBasis = sourceBasisValues.getGradBasisValues(
false).get_static_view();
519 sourceUnweightedVectorBasis = sourceBasisValues.getVectorBasisValues(
false).get_static_view();
523 PHX::View<panzer::LocalOrdinal**> targetLocalIds(
"buildRHSMatrix: targetLocalIds", workset.numOwnedCells(),
524 targetGlobalIndexer_->getElementBlockGIDCount(block));
525 PHX::View<panzer::LocalOrdinal**> sourceLocalIds(
"buildRHSMatrix: sourceLocalIds", workset.numOwnedCells(),
529 const auto cellLocalIdsNoGhost = Kokkos::subview(workset.cell_local_ids_k,std::make_pair(0,workset.numOwnedCells()));
530 targetGlobalIndexer_->getElementLIDs(cellLocalIdsNoGhost,targetLocalIds);
531 sourceGlobalIndexer.
getElementLIDs(cellLocalIdsNoGhost,sourceLocalIds);
535 PHX::View<panzer::LocalOrdinal*> targetFieldOffsets;
537 const auto fieldIndex = targetGlobalIndexer_->getFieldNum(targetBasisDescriptor_.getType());
538 const std::vector<panzer::LocalOrdinal>&
offsets = targetGlobalIndexer_->getGIDFieldOffsets(block,fieldIndex);
539 targetFieldOffsets = PHX::View<panzer::LocalOrdinal*>(
"L2Projection:buildRHS:targetFieldOffsets",offsets.size());
540 const auto hostOffsets = Kokkos::create_mirror_view(targetFieldOffsets);
541 for(
size_t i=0; i < offsets.size(); ++i)
542 hostOffsets(i) = offsets[i];
543 Kokkos::deep_copy(targetFieldOffsets,hostOffsets);
546 PHX::View<panzer::LocalOrdinal*> sourceFieldOffsets;
548 const auto fieldIndex = sourceGlobalIndexer.
getFieldNum(sourceFieldName);
551 "ERROR: panzer::L2Projection::buildRHSMatrix() - The source field, \""
552 << sourceFieldName <<
"\", does not exist in element block \""
554 sourceFieldOffsets = PHX::View<panzer::LocalOrdinal*>(
"L2Projection:buildRHS:sourceFieldOffsets",offsets.size());
555 const auto hostOffsets = Kokkos::create_mirror_view(sourceFieldOffsets);
556 for(
size_t i=0; i <offsets.size(); ++i)
557 hostOffsets(i) = offsets[i];
558 Kokkos::deep_copy(sourceFieldOffsets,hostOffsets);
561 const auto localMatrix = ghostedMatrix->getLocalMatrixDevice();
562 const int numRows =
static_cast<int>(targetWeightedBasis.extent(1));
565 if (useRankThreeBasis) {
566 tmpNumCols =
static_cast<int>(sourceUnweightedScalarBasis.extent(1));
567 tmpNumQP =
static_cast<int>(sourceUnweightedScalarBasis.extent(2));
570 tmpNumCols =
static_cast<int>(sourceUnweightedVectorBasis.extent(1));
571 tmpNumQP =
static_cast<int>(sourceUnweightedVectorBasis.extent(2));
573 const int numCols = tmpNumCols;
574 const int numQP = tmpNumQP;
576 Kokkos::parallel_for(Kokkos::RangePolicy<PHX::Device>(0,workset.numOwnedCells()),KOKKOS_LAMBDA (
const int& cell) {
577 panzer::LocalOrdinal cLIDs[256];
579 for (
int row = 0; row < numRows; ++row) {
580 const int rowOffset = targetFieldOffsets(row);
581 const int rowLID = targetLocalIds(cell,rowOffset);
582 for (
int col = 0; col < numCols; ++col)
585 for (
int col = 0; col < numCols; ++col) {
586 for (
int qp = 0; qp < numQP; ++qp) {
587 const int colOffset = sourceFieldOffsets(col);
588 const int colLID = sourceLocalIds(cell,colOffset);
590 if (useRankThreeBasis)
591 vals[col] += sourceUnweightedScalarBasis(cell,col,qp) * targetWeightedBasis(cell,row,qp);
593 vals[col] += sourceUnweightedVectorBasis(cell,col,qp,directionIndex) * targetWeightedBasis(cell,row,qp);
596 localMatrix.sumIntoValues(rowLID,cLIDs,numCols,vals,
true,
true);
602 ghostedMatrix->fillComplete(ownedSourceMap,ownedTargetMap);
603 ownedMatrix->resumeFill();
604 ownedMatrix->doExport(*ghostedMatrix,*exporter,Tpetra::ADD);
605 ownedMatrix->fillComplete(ownedSourceMap,ownedTargetMap);
Teuchos::RCP< CrsMatrixType > getTpetraMatrix(int i, int j) const
bool is_null(const boost::shared_ptr< T > &p)
Teuchos::RCP< Tpetra::CrsMatrix< double, panzer::LocalOrdinal, panzer::GlobalOrdinal, Tpetra::KokkosCompat::KokkosDeviceWrapperNode< PHX::Device > > > buildRHSMatrix(const panzer::GlobalIndexer &sourceDOFManager, const Teuchos::RCP< const Tpetra::Map< panzer::LocalOrdinal, panzer::GlobalOrdinal, Tpetra::KokkosCompat::KokkosDeviceWrapperNode< PHX::Device >>> &ownedSourceMap, const std::string &sourceFieldName, const panzer::BasisDescriptor &sourceBasisDescriptor, const int vectorOrGradientDirectionIndex=-1)
Allocates, fills and returns a rectangular matrix for L2 projection of a scalar field, one dimension of gradient (for hgrad basis), or one dimension of a vector field onto the target scalar basis. If you wish to project all values of a vector field or all the gradients of a scalar field, then you will need three separate RHS matrices to form the RHS for each independently. The vectors must be independent Tpetra vectors to solve multiple right hand sides with the linear solver.
virtual Teuchos::RCP< const ExportType > getGhostedExport(int j) const
get exporter for converting an overalapped object to a "normal" object
void setup(const panzer::BasisDescriptor &targetBasis, const panzer::IntegrationDescriptor &integrationDescriptor, const Teuchos::RCP< const Teuchos::MpiComm< int >> &comm, const Teuchos::RCP< const panzer::ConnManager > &connManager, const std::vector< std::string > &elementBlockNames, const Teuchos::RCP< panzer::WorksetContainer > worksetContainer=Teuchos::null)
Setup base objects for L2 Projections - requires target scalar basis and creates worksets if not supp...
const Kokkos::View< const panzer::LocalOrdinal *, Kokkos::LayoutRight, PHX::Device > getElementLIDs(panzer::LocalOrdinal localElmtId) const
virtual int getElementBlockGIDCount(const std::size_t &blockIndex) const =0
How any GIDs are associate with each element in a particular element block.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Teuchos::RCP< panzer::GlobalIndexer > getTargetGlobalIndexer() const
Returns the target global indexer. Will be null if setup() has not been called.
PHX::View< const int * > offsets
virtual void getOwnedIndices(std::vector< panzer::GlobalOrdinal > &indices) const =0
Get the set of indices owned by this processor.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
const std::string & getType() const
Get type of basis.
ConstArray_CellBasisIP getBasisValues(const bool weighted, const bool cache=true, const bool force=false) const
Get the basis values evaluated at mesh points.
Kokkos::View< typename Sacado::ScalarType< ScalarT >::type **, Kokkos::LayoutRight, PHX::Device > vals
virtual int getFieldNum(const std::string &str) const =0
Get the number used for access to this field.
virtual void getOwnedAndGhostedIndices(std::vector< panzer::GlobalOrdinal > &indices) const =0
Get the set of owned and ghosted indices for this processor.
virtual Teuchos::RCP< const MapType > getMap(int i) const
get the map from the matrix
void useBasisValues(const std::map< std::string, Teuchos::RCP< panzer::BasisValues2< double >>> &map_eblock_to_bv)
Override using the panzer::WorksetContainer and instead use the registered BasisValues object...
ConstArray_CellBasisIPDim getVectorBasisValues(const bool weighted, const bool cache=true, const bool force=false) const
Get the vector basis values evaluated at mesh points.
KOKKOS_FORCEINLINE_FUNCTION array_type get_static_view()
Teuchos::RCP< CrsMatrixType > getGhostedTpetraMatrix(int i, int j) const
virtual const std::vector< int > & getGIDFieldOffsets(const std::string &blockId, int fieldNum) const =0
Use the field pattern so that you can find a particular field in the GIDs array.
virtual void getElementGIDs(panzer::LocalOrdinal localElmtId, std::vector< panzer::GlobalOrdinal > &gids, const std::string &blockIdHint="") const =0
Get the global IDs for a particular element. This function overwrites the gids variable.
Teuchos::RCP< Tpetra::MultiVector< double, panzer::LocalOrdinal, panzer::GlobalOrdinal, Tpetra::KokkosCompat::KokkosDeviceWrapperNode< PHX::Device > > > buildInverseLumpedMassMatrix()
Allocates, fills and returns a Tpetra::MultiVector containing the inverse lumped mass matrix values a...
#define TEUCHOS_ASSERT(assertion_test)
Workset size is set to the total number of local elements in the MPI process.
Teuchos::RCP< Tpetra::CrsMatrix< double, panzer::LocalOrdinal, panzer::GlobalOrdinal, Tpetra::KokkosCompat::KokkosDeviceWrapperNode< PHX::Device > > > buildMassMatrix(bool use_lumping=false, const std::unordered_map< std::string, double > *elementBlockMultipliers=nullptr)
Allocates, fills and returns a mass matrix for L2 projection onto a target basis. ...