46 #ifndef MUELU_GEOMETRICINTERPOLATIONPFACTORY_DEF_HPP 47 #define MUELU_GEOMETRICINTERPOLATIONPFACTORY_DEF_HPP 49 #include "Xpetra_CrsGraph.hpp" 50 #include "Xpetra_CrsMatrixUtils.hpp" 54 #include "MueLu_Aggregates.hpp" 61 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
63 RCP<ParameterList> validParamList = rcp(
new ParameterList());
65 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name)) 68 #undef SET_VALID_ENTRY 71 validParamList->set<RCP<const FactoryBase> >(
"A", Teuchos::null,
72 "Generating factory of the matrix A");
73 validParamList->set<RCP<const FactoryBase> >(
"Aggregates", Teuchos::null,
74 "Aggregates generated by StructuredAggregationFactory used to construct a piece-constant prolongator.");
75 validParamList->set<RCP<const FactoryBase> >(
"prolongatorGraph", Teuchos::null,
76 "Graph generated by StructuredAggregationFactory used to construct a piece-linear prolongator.");
77 validParamList->set<RCP<const FactoryBase> >(
"Coordinates", Teuchos::null,
78 "Fine level coordinates used to construct piece-wise linear prolongator and coarse level coordinates.");
79 validParamList->set<RCP<const FactoryBase> >(
"coarseCoordinatesFineMap", Teuchos::null,
80 "map of the coarse coordinates' GIDs as indexed on the fine mesh.");
81 validParamList->set<RCP<const FactoryBase> >(
"coarseCoordinatesMap", Teuchos::null,
82 "map of the coarse coordinates' GIDs as indexed on the coarse mesh.");
83 validParamList->set<RCP<const FactoryBase> >(
"Nullspace", Teuchos::null,
84 "Fine level nullspace used to construct the coarse level nullspace.");
85 validParamList->set<RCP<const FactoryBase> >(
"lCoarseNodesPerDim", Teuchos::null,
86 "Number of nodes per spatial dimension on the coarse grid.");
87 validParamList->set<RCP<const FactoryBase> >(
"numDimensions", Teuchos::null,
88 "Number of spatial dimensions of the mesh.");
89 validParamList->set<RCP<const FactoryBase> >(
"interpolationOrder", Teuchos::null,
90 "Interpolation order use to construct the prolongator.");
92 return validParamList;
95 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
98 const ParameterList& pL = GetParameterList();
100 Input(fineLevel,
"A");
101 Input(fineLevel,
"Nullspace");
102 Input(fineLevel,
"prolongatorGraph");
103 Input(fineLevel,
"lCoarseNodesPerDim");
104 Input(fineLevel,
"numDimensions");
105 Input(fineLevel,
"interpolationOrder");
107 if( pL.get<
bool>(
"gmg: build coarse coordinates") ||
108 (pL.get<
int>(
"gmg: interpolation order") == 1) ) {
109 Input(fineLevel,
"Coordinates");
110 Input(fineLevel,
"coarseCoordinatesFineMap");
111 Input(fineLevel,
"coarseCoordinatesMap");
116 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
119 return BuildP(fineLevel, coarseLevel);
122 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
128 RCP<Teuchos::FancyOStream> out;
129 if(
const char* dbg = std::getenv(
"MUELU_GEOMETRICINTERPOLATIONPFACTORY_DEBUG")) {
130 out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
131 out->setShowAllFrontMatter(
false).setShowProcRank(
true);
133 out = Teuchos::getFancyOStream(rcp(
new Teuchos::oblackholestream()));
136 *out <<
"Starting GeometricInterpolationPFactory::BuildP." << std::endl;
139 const ParameterList& pL = GetParameterList();
140 const bool buildCoarseCoordinates = pL.get<
bool>(
"gmg: build coarse coordinates");
143 RCP<Matrix> A = Get<RCP<Matrix> >(fineLevel,
"A");
144 RCP<CrsGraph> prolongatorGraph = Get<RCP<CrsGraph> >(fineLevel,
"prolongatorGraph");
145 RCP<realvaluedmultivector_type> fineCoordinates, coarseCoordinates;
147 const int numDimensions = Get<int>(fineLevel,
"numDimensions");
148 const int interpolationOrder = Get<int>(fineLevel,
"interpolationOrder");
152 if(buildCoarseCoordinates || (interpolationOrder == 1)) {
153 RCP<const Map> coarseCoordsFineMap = Get< RCP<const Map> >(fineLevel,
"coarseCoordinatesFineMap");
154 RCP<const Map> coarseCoordsMap = Get< RCP<const Map> >(fineLevel,
"coarseCoordinatesMap");
155 fineCoordinates = Get< RCP<realvaluedmultivector_type> >(fineLevel,
"Coordinates");
156 coarseCoordinates = Xpetra::MultiVectorFactory<real_type,LO,GO,Node>::Build(coarseCoordsFineMap,
157 fineCoordinates->getNumVectors());
158 RCP<const Import> coordsImporter = ImportFactory::Build(fineCoordinates->getMap(),
159 coarseCoordsFineMap);
160 coarseCoordinates->doImport(*fineCoordinates, *coordsImporter, Xpetra::INSERT);
161 coarseCoordinates->replaceMap(coarseCoordsMap);
163 Set(coarseLevel,
"Coordinates", coarseCoordinates);
166 *out <<
"Fine and coarse coordinates have been loaded from the fine level and set on the coarse level." << std::endl;
168 if(interpolationOrder == 0) {
170 BuildConstantP(P, prolongatorGraph, A);
171 }
else if(interpolationOrder == 1) {
174 RCP<realvaluedmultivector_type> ghostCoordinates
175 = Xpetra::MultiVectorFactory<real_type,LO,GO,NO>::Build(prolongatorGraph->getColMap(),
176 fineCoordinates->getNumVectors());
177 RCP<const Import> ghostImporter = ImportFactory::Build(coarseCoordinates->getMap(),
178 prolongatorGraph->getColMap());
179 ghostCoordinates->doImport(*coarseCoordinates, *ghostImporter, Xpetra::INSERT);
181 BuildLinearP(A, prolongatorGraph, fineCoordinates, ghostCoordinates, numDimensions, P);
184 *out <<
"The prolongator matrix has been build." << std::endl;
187 RCP<MultiVector> fineNullspace = Get< RCP<MultiVector> > (fineLevel,
"Nullspace");
188 RCP<MultiVector> coarseNullspace = MultiVectorFactory::Build(P->getDomainMap(),
189 fineNullspace->getNumVectors());
190 P->apply(*fineNullspace, *coarseNullspace, Teuchos::TRANS, Teuchos::ScalarTraits<SC>::one(),
191 Teuchos::ScalarTraits<SC>::zero());
192 Set(coarseLevel,
"Nullspace", coarseNullspace);
194 *out <<
"The coarse nullspace is constructed and set on the coarse level." << std::endl;
196 Array<LO> lNodesPerDir = Get<Array<LO> >(fineLevel,
"lCoarseNodesPerDim");
197 Set(coarseLevel,
"lNodesPerDim", lNodesPerDir);
198 Set(coarseLevel,
"P", P);
200 *out <<
"GeometricInterpolationPFactory::BuildP has completed." << std::endl;
204 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
206 BuildConstantP(RCP<Matrix>& P, RCP<CrsGraph>& prolongatorGraph, RCP<Matrix>& A)
const {
209 RCP<Teuchos::FancyOStream> out;
210 if(
const char* dbg = std::getenv(
"MUELU_GEOMETRICINTERPOLATIONPFACTORY_DEBUG")) {
211 out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
212 out->setShowAllFrontMatter(
false).setShowProcRank(
true);
214 out = Teuchos::getFancyOStream(rcp(
new Teuchos::oblackholestream()));
219 int dofsPerNode = A->GetFixedBlockSize();
220 ArrayView<const GO> initialDomainMapLIDs =
221 prolongatorGraph->getDomainMap()->getNodeElementList();
222 Array<GO> domainMapLIDs(initialDomainMapLIDs.size()*dofsPerNode);
223 for(LO elementIdx = 0; elementIdx < as<LO>(domainMapLIDs.size()); ++elementIdx) {
224 for(
int dof = 0; dof < dofsPerNode; ++dof) {
225 domainMapLIDs[elementIdx*dofsPerNode + dof] =
226 initialDomainMapLIDs[elementIdx]*dofsPerNode + dof;
229 RCP<Map> domainMap = MapFactory::Build(prolongatorGraph->getRowMap()->lib(),
230 prolongatorGraph->getGlobalNumCols()*dofsPerNode,
232 prolongatorGraph->getIndexBase(),
233 prolongatorGraph->getComm(),
234 prolongatorGraph->getRowMap()->getNode());
236 ArrayView<const GO> initialColMapLIDs =
237 prolongatorGraph->getColMap()->getNodeElementList();
238 Array<GO> colMapLIDs(initialColMapLIDs.size()*dofsPerNode);
239 for(LO elementIdx = 0; elementIdx < as<LO>(colMapLIDs.size()); ++elementIdx) {
240 for(
int dof = 0; dof < dofsPerNode; ++dof) {
241 colMapLIDs[elementIdx*dofsPerNode + dof] =
242 initialColMapLIDs[elementIdx]*dofsPerNode + dof;
245 RCP<Map> colMap = MapFactory::Build(prolongatorGraph->getColMap()->lib(),
246 prolongatorGraph->getGlobalNumCols()*dofsPerNode,
248 prolongatorGraph->getIndexBase(),
249 prolongatorGraph->getComm(),
250 prolongatorGraph->getColMap()->getNode());
252 std::vector<size_t> strideInfo(1);
253 strideInfo[0] = dofsPerNode;
254 RCP<const StridedMap> stridedDomainMap = StridedMapFactory::Build(domainMap, strideInfo);
257 P = rcp(
new CrsMatrixWrap(A->getDomainMap(), colMap, 0, Xpetra::StaticProfile));
259 RCP<CrsMatrix> PCrs = rcp_dynamic_cast<CrsMatrixWrap>(P)->getCrsMatrix();
261 ArrayRCP<size_t> iaP;
265 PCrs->allocateAllValues(A->getDomainMap()->getNodeNumElements(), iaP, jaP, valP);
267 ArrayView<size_t> ia = iaP();
268 ArrayView<LO> ja = jaP();
269 ArrayView<SC> val = valP();
274 ArrayView<const LO> colIdx;
275 for(LO rowIdx = 0; rowIdx < prolongatorGraph->getNodeNumRows(); ++rowIdx) {
276 prolongatorGraph->getLocalRowView(rowIdx, colIdx);
277 for(
int dof = 0; dof < dofsPerNode; ++dof) {
278 dofIdx = rowIdx*dofsPerNode + dof;
279 ia[dofIdx + 1] = dofIdx + 1;
280 ja[dofIdx] = colIdx[0]*dofsPerNode + dof;
286 PCrs->setAllValues(iaP, jaP, valP);
287 PCrs->expertStaticFillComplete(domainMap, A->getDomainMap());
290 if (A->IsView(
"stridedMaps") ==
true) {
291 P->CreateView(
"stridedMaps", A->getRowMap(
"stridedMaps"), stridedDomainMap);
293 P->CreateView(
"stridedMaps", P->getRangeMap(), stridedDomainMap);
298 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
301 RCP<realvaluedmultivector_type>& fineCoordinates,
302 RCP<realvaluedmultivector_type>& ghostCoordinates,
303 const int numDimensions, RCP<Matrix>& P)
const {
306 RCP<Teuchos::FancyOStream> out;
307 if(
const char* dbg = std::getenv(
"MUELU_GEOMETRICINTERPOLATIONPFACTORY_DEBUG")) {
308 out = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout));
309 out->setShowAllFrontMatter(
false).setShowProcRank(
true);
311 out = Teuchos::getFancyOStream(rcp(
new Teuchos::oblackholestream()));
314 *out <<
"Entering BuildLinearP" << std::endl;
317 Array<ArrayView<const real_type> > fineCoords(3);
318 Array<ArrayView<const real_type> > ghostCoords(3);
319 real_type realZero = Teuchos::as<real_type>(0.0);
320 Array<real_type> fineZero(fineCoordinates->getLocalLength(), realZero);
321 Array<real_type> ghostZero(ghostCoordinates->getLocalLength(), realZero);
322 for(
int dim = 0; dim < 3; ++dim) {
323 if(dim < numDimensions) {
324 fineCoords[dim] = fineCoordinates->getData(dim)();
325 ghostCoords[dim] = ghostCoordinates->getData(dim)();
327 fineCoords[dim] = fineZero();
328 ghostCoords[dim] = ghostZero();
332 *out <<
"Coordinates extracted from the multivectors!" << std::endl;
335 const int numInterpolationPoints = 1 << numDimensions;
336 const int dofsPerNode = A->GetFixedBlockSize();
337 ArrayView<const GO> initialDomainMapLIDs =
338 prolongatorGraph->getDomainMap()->getNodeElementList();
339 Array<GO> domainMapLIDs(initialDomainMapLIDs.size()*dofsPerNode);
340 for(LO elementIdx = 0; elementIdx < as<LO>(domainMapLIDs.size()); ++elementIdx) {
341 for(
int dof = 0; dof < dofsPerNode; ++dof) {
342 domainMapLIDs[elementIdx*dofsPerNode + dof] =
343 initialDomainMapLIDs[elementIdx]*dofsPerNode + dof;
346 RCP<Map> domainMap = MapFactory::Build(prolongatorGraph->getRowMap()->lib(),
347 prolongatorGraph->getGlobalNumCols()*dofsPerNode,
349 prolongatorGraph->getIndexBase(),
350 prolongatorGraph->getComm(),
351 prolongatorGraph->getRowMap()->getNode());
353 ArrayView<const GO> initialColMapLIDs =
354 prolongatorGraph->getColMap()->getNodeElementList();
355 Array<GO> colMapLIDs(initialColMapLIDs.size()*dofsPerNode);
356 for(LO elementIdx = 0; elementIdx < as<LO>(colMapLIDs.size()); ++elementIdx) {
357 for(
int dof = 0; dof < dofsPerNode; ++dof) {
358 colMapLIDs[elementIdx*dofsPerNode + dof] =
359 initialColMapLIDs[elementIdx]*dofsPerNode + dof;
362 RCP<Map> colMap = MapFactory::Build(prolongatorGraph->getColMap()->lib(),
363 prolongatorGraph->getGlobalNumCols()*dofsPerNode,
365 prolongatorGraph->getIndexBase(),
366 prolongatorGraph->getComm(),
367 prolongatorGraph->getColMap()->getNode());
369 std::vector<size_t> strideInfo(1);
370 strideInfo[0] = dofsPerNode;
371 RCP<const StridedMap> stridedDomainMap = StridedMapFactory::Build(domainMap, strideInfo);
373 *out <<
"The maps of P have been computed" << std::endl;
375 P = rcp(
new CrsMatrixWrap(A->getDomainMap(), colMap, 0, Xpetra::StaticProfile));
377 RCP<CrsMatrix> PCrs = rcp_dynamic_cast<CrsMatrixWrap>(P)->getCrsMatrix();
379 ArrayRCP<size_t> iaP;
383 LO numNonZeroP = prolongatorGraph->getNodeNumEntries()*dofsPerNode;
384 PCrs->allocateAllValues(numNonZeroP, iaP, jaP, valP);
386 *out <<
"dofsPerNode=" << dofsPerNode << std::endl;
387 *out <<
"number of non-zeroes in P: " << numNonZeroP << std::endl;
389 ArrayView<size_t> ia = iaP();
390 ArrayView<LO> ja = jaP();
391 ArrayView<SC> val = valP();
394 LO rowDofIdx, colDofIdx;
395 ArrayView<const LO> colIndices;
396 Array<Array<real_type> > coords(numInterpolationPoints + 1);
397 Array<real_type> stencil(numInterpolationPoints);
398 for(LO rowIdx = 0; rowIdx < prolongatorGraph->getNodeNumRows(); ++rowIdx) {
399 prolongatorGraph->getLocalRowView(rowIdx, colIndices);
403 if(colIndices.size() == 1) {
406 for(
int dof = 0; dof < dofsPerNode; ++dof) {
407 rowDofIdx = rowIdx*dofsPerNode + dof;
408 ia[rowDofIdx + 1] = ia[rowDofIdx] + 1;
410 colDofIdx = ia[rowDofIdx];
411 ja[colDofIdx] = colIndices[0]*dofsPerNode + dof;
412 val[colDofIdx] = 1.0;
416 for(
int dof = 0; dof < dofsPerNode; ++dof) {
417 rowDofIdx = rowIdx*dofsPerNode + dof;
418 ia[rowDofIdx + 1] = ia[rowDofIdx] + colIndices.size();
423 for(
int dim = 0; dim < 3; ++dim) {
424 coords[0][dim] = fineCoords[dim][rowIdx];
426 for(
int interpolationIdx = 0; interpolationIdx < numInterpolationPoints; ++interpolationIdx) {
427 coords[interpolationIdx + 1].resize(3);
428 for(
int dim = 0; dim < 3; ++dim) {
429 coords[interpolationIdx + 1][dim] = ghostCoords[dim][colIndices[interpolationIdx]];
432 ComputeLinearInterpolationStencil(numDimensions, numInterpolationPoints, coords, stencil);
434 for(
int colIdx = 0; colIdx < colIndices.size(); ++colIdx) {
435 colDofIdx = ia[rowDofIdx] + colIdx;
436 ja[colDofIdx] = colIndices[colIdx]*dofsPerNode + dof;
437 val[colDofIdx] = stencil[colIdx];
443 *out <<
"The calculation of the interpolation stencils has completed." << std::endl;
445 Xpetra::CrsMatrixUtils<SC,LO,GO,NO>::sortCrsEntries(ia, ja, val, A->getDomainMap()->lib());
446 PCrs->setAllValues(iaP, jaP, valP);
447 PCrs->expertStaticFillComplete(domainMap, A->getDomainMap());
449 *out <<
"All values in P have been set and expertStaticFillComplete has been performed." << std::endl;
452 if (A->IsView(
"stridedMaps") ==
true) {
453 P->CreateView(
"stridedMaps", A->getRowMap(
"stridedMaps"), stridedDomainMap);
455 P->CreateView(
"stridedMaps", P->getRangeMap(), stridedDomainMap);
461 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
464 const Array<Array<real_type> > coord,
465 Array<real_type>& stencil)
const {
482 Teuchos::SerialDenseMatrix<LO,real_type> Jacobian(numDimensions, numDimensions);
483 Teuchos::SerialDenseVector<LO,real_type> residual(numDimensions);
484 Teuchos::SerialDenseVector<LO,real_type> solutionDirection(numDimensions);
485 Teuchos::SerialDenseVector<LO,real_type> paramCoords(numDimensions);
486 Teuchos::SerialDenseSolver<LO,real_type> problem;
487 int iter = 0, max_iter = 5;
488 real_type functions[4][8], norm_ref = 1.0, norm2 = 1.0, tol = 1.0e-5;
489 paramCoords.size(numDimensions);
491 while( (iter < max_iter) && (norm2 > tol*norm_ref) ) {
494 solutionDirection.size(numDimensions);
495 residual.size(numDimensions);
499 GetInterpolationFunctions(numDimensions, paramCoords, functions);
500 for(LO i = 0; i < numDimensions; ++i) {
501 residual(i) = coord[0][i];
502 for(LO k = 0; k < numInterpolationPoints; ++k) {
503 residual(i) -= functions[0][k]*coord[k+1][i];
506 norm_ref += residual(i)*residual(i);
507 if(i == numDimensions - 1) {
508 norm_ref = std::sqrt(norm_ref);
512 for(LO j = 0; j < numDimensions; ++j) {
513 for(LO k = 0; k < numInterpolationPoints; ++k) {
514 Jacobian(i,j) += functions[j+1][k]*coord[k+1][i];
520 problem.setMatrix(Teuchos::rcp(&Jacobian,
false));
521 problem.setVectors(Teuchos::rcp(&solutionDirection,
false), Teuchos::rcp(&residual,
false));
522 if(problem.shouldEquilibrate()) {problem.factorWithEquilibration(
true);}
525 for(LO i = 0; i < numDimensions; ++i) {
526 paramCoords(i) = paramCoords(i) + solutionDirection(i);
530 GetInterpolationFunctions(numDimensions, paramCoords, functions);
531 for(LO i = 0; i < numDimensions; ++i) {
533 for(LO k = 0; k < numInterpolationPoints; ++k) {
534 tmp -= functions[0][k]*coord[k+1][i];
539 norm2 = std::sqrt(norm2);
543 for(LO i = 0; i < 8; ++i) {
544 stencil[i] = functions[0][i];
549 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
552 const Teuchos::SerialDenseVector<LO, real_type> parametricCoordinates,
554 real_type xi = 0.0, eta = 0.0, zeta = 0.0, denominator = 0.0;
555 if(numDimensions == 1) {
556 xi = parametricCoordinates[0];
558 }
else if(numDimensions == 2) {
559 xi = parametricCoordinates[0];
560 eta = parametricCoordinates[1];
562 }
else if(numDimensions == 3) {
563 xi = parametricCoordinates[0];
564 eta = parametricCoordinates[1];
565 zeta = parametricCoordinates[2];
569 functions[0][0] = (1.0 - xi)*(1.0 - eta)*(1.0 - zeta) / denominator;
570 functions[0][1] = (1.0 + xi)*(1.0 - eta)*(1.0 - zeta) / denominator;
571 functions[0][2] = (1.0 - xi)*(1.0 + eta)*(1.0 - zeta) / denominator;
572 functions[0][3] = (1.0 + xi)*(1.0 + eta)*(1.0 - zeta) / denominator;
573 functions[0][4] = (1.0 - xi)*(1.0 - eta)*(1.0 + zeta) / denominator;
574 functions[0][5] = (1.0 + xi)*(1.0 - eta)*(1.0 + zeta) / denominator;
575 functions[0][6] = (1.0 - xi)*(1.0 + eta)*(1.0 + zeta) / denominator;
576 functions[0][7] = (1.0 + xi)*(1.0 + eta)*(1.0 + zeta) / denominator;
578 functions[1][0] = -(1.0 - eta)*(1.0 - zeta) / denominator;
579 functions[1][1] = (1.0 - eta)*(1.0 - zeta) / denominator;
580 functions[1][2] = -(1.0 + eta)*(1.0 - zeta) / denominator;
581 functions[1][3] = (1.0 + eta)*(1.0 - zeta) / denominator;
582 functions[1][4] = -(1.0 - eta)*(1.0 + zeta) / denominator;
583 functions[1][5] = (1.0 - eta)*(1.0 + zeta) / denominator;
584 functions[1][6] = -(1.0 + eta)*(1.0 + zeta) / denominator;
585 functions[1][7] = (1.0 + eta)*(1.0 + zeta) / denominator;
587 functions[2][0] = -(1.0 - xi)*(1.0 - zeta) / denominator;
588 functions[2][1] = -(1.0 + xi)*(1.0 - zeta) / denominator;
589 functions[2][2] = (1.0 - xi)*(1.0 - zeta) / denominator;
590 functions[2][3] = (1.0 + xi)*(1.0 - zeta) / denominator;
591 functions[2][4] = -(1.0 - xi)*(1.0 + zeta) / denominator;
592 functions[2][5] = -(1.0 + xi)*(1.0 + zeta) / denominator;
593 functions[2][6] = (1.0 - xi)*(1.0 + zeta) / denominator;
594 functions[2][7] = (1.0 + xi)*(1.0 + zeta) / denominator;
596 functions[3][0] = -(1.0 - xi)*(1.0 - eta) / denominator;
597 functions[3][1] = -(1.0 + xi)*(1.0 - eta) / denominator;
598 functions[3][2] = -(1.0 - xi)*(1.0 + eta) / denominator;
599 functions[3][3] = -(1.0 + xi)*(1.0 + eta) / denominator;
600 functions[3][4] = (1.0 - xi)*(1.0 - eta) / denominator;
601 functions[3][5] = (1.0 + xi)*(1.0 - eta) / denominator;
602 functions[3][6] = (1.0 - xi)*(1.0 + eta) / denominator;
603 functions[3][7] = (1.0 + xi)*(1.0 + eta) / denominator;
609 #endif // MUELU_GEOMETRICINTERPOLATIONPFACTORY_DEF_HPP void BuildP(Level &fineLevel, Level &coarseLevel) const
Abstract Build method.
typename Teuchos::ScalarTraits< SC >::magnitudeType real_type
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
void GetInterpolationFunctions(const LO numDimensions, const Teuchos::SerialDenseVector< LO, real_type > parametricCoordinates, real_type functions[4][8]) const
Timer to be used in factories. Similar to Monitor but with additional timers.
Namespace for MueLu classes and methods.
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
Class that holds all level-specific information.
#define SET_VALID_ENTRY(name)
void BuildLinearP(RCP< Matrix > &A, RCP< CrsGraph > &prolongatorGraph, RCP< realvaluedmultivector_type > &fineCoordinates, RCP< realvaluedmultivector_type > &ghostCoordinates, const int numDimensions, RCP< Matrix > &P) const
void BuildConstantP(RCP< Matrix > &P, RCP< CrsGraph > &prolongatorGraph, RCP< Matrix > &A) const
void ComputeLinearInterpolationStencil(const int numDimensions, const int numInterpolationPoints, const Array< Array< real_type > > coord, Array< real_type > &stencil) const