46 #ifndef MUELU_GEOINTERPFACTORY_DEF_HPP
47 #define MUELU_GEOINTERPFACTORY_DEF_HPP
54 #include <Xpetra_Map.hpp>
55 #include <Xpetra_MapFactory.hpp>
58 #include <Xpetra_MultiVectorFactory.hpp>
60 #include <Xpetra_CrsMatrixWrap.hpp>
67 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
69 GetOStream(
Runtime1) <<
"I constructed a GeoInterpFactory object... Nothing else to do here." << std::endl;
72 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
77 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
79 Input(fineLevel,
"A");
80 Input(fineLevel,
"A00");
81 Input(fineLevel,
"A10");
82 Input(fineLevel,
"A20");
84 Input(fineLevel,
"VElementList");
85 Input(fineLevel,
"PElementList");
86 Input(fineLevel,
"MElementList");
88 Input(coarseLevel,
"VElementList");
89 Input(coarseLevel,
"PElementList");
90 Input(coarseLevel,
"MElementList");
105 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
107 GetOStream(
Runtime1) <<
"Starting 'build' routine...\n";
115 return BuildP(fineLevel, coarseLevel);
118 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
122 GetOStream(
Runtime1) <<
"Starting 'BuildP' routine...\n";
134 GetOStream(
Runtime1) <<
"done getting fine level elements...\n";
135 GetOStream(
Runtime1) <<
"getting coarse level elements...\n";
143 coarseLevel.
Get(
"VElementList", coarseElementVDOFs, coarseLevel.
GetFactoryManager()->GetFactory(
"VElementList").get());
144 coarseLevel.
Get(
"PElementList", coarseElementPDOFs, coarseLevel.
GetFactoryManager()->GetFactory(
"PElementList").get());
145 coarseLevel.
Get(
"MElementList", coarseElementMDOFs, coarseLevel.
GetFactoryManager()->GetFactory(
"MElementList").get());
147 GetOStream(
Runtime1) <<
"computing various numbers...\n";
149 GO totalFineElements = fineElementMDOFs->numRows();
150 LO nFineElements = (int)sqrt(totalFineElements);
151 GO totalCoarseElements = coarseElementMDOFs->numRows();
152 LO nCoarseElements = (int)sqrt(totalCoarseElements);
155 GO nM = (2 * nCoarseElements + 1) * (2 * nCoarseElements + 1);
157 GO nP = (nCoarseElements + 1) * (nCoarseElements + 1);
160 RCP<Matrix> fineA00 = Get<RCP<Matrix> >(fineLevel,
"A00");
161 RCP<Matrix> fineA10 = Get<RCP<Matrix> >(fineLevel,
"A10");
162 RCP<Matrix> fineA20 = Get<RCP<Matrix> >(fineLevel,
"A20");
164 GetOStream(
Runtime1) <<
"creating coarse grid maps...\n";
170 GO fNV = rowMapforPV->getGlobalNumElements();
171 GO fNP = rowMapforPP->getGlobalNumElements();
172 GO fNM = rowMapforPM->getGlobalNumElements();
178 RCP<Matrix> FineA = Factory::Get<RCP<Matrix> >(fineLevel,
"A");
186 GetOStream(
Runtime1) <<
"creating coarse grid matrices...\n";
188 size_t maxEntriesPerRowV = 9,
189 maxEntriesPerRowP = 4,
190 maxEntriesPerRowM = 9;
192 RCP<Matrix> P =
rcp(
new CrsMatrixWrap(rowMapforP, maxEntriesPerRowV));
193 RCP<Matrix> PV =
rcp(
new CrsMatrixWrap(rowMapforPV, maxEntriesPerRowV));
194 RCP<Matrix> PP =
rcp(
new CrsMatrixWrap(rowMapforPP, maxEntriesPerRowP));
195 RCP<Matrix> PM =
rcp(
new CrsMatrixWrap(rowMapforPM, maxEntriesPerRowM));
254 GO fineElement[4] = {0, 1, nFineElements, nFineElements + 1};
256 GetOStream(
Runtime1) <<
"start building matrices...\n";
257 GetOStream(
Runtime1) <<
"nCoarseElements = " << nCoarseElements << std::endl;
259 for (
GO coarseElement = 0; coarseElement < totalCoarseElements; coarseElement++) {
268 if (coarseElement < nCoarseElements) {
271 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 0);
273 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 1);
275 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 4);
279 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[0], 4), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
281 colPtrV[0] = 2 * colPtrM[0];
282 colPtrV[1] = 2 * colPtrM[1];
283 colPtrV[2] = 2 * colPtrM[2];
285 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 8), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
287 colPtrV[0] = 2 * colPtrM[0] + 1;
288 colPtrV[1] = 2 * colPtrM[1] + 1;
289 colPtrV[2] = 2 * colPtrM[2] + 1;
291 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 9), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
294 colPtrP[0] = (*coarseElementPDOFs)(coarseElement, 0);
296 colPtrP[1] = (*coarseElementPDOFs)(coarseElement, 1);
300 PP->insertGlobalValues((*fineElementPDOFs)(fineElement[0], 1), colPtrP.
view(0, nnz), valPtrP.view(0, nnz));
303 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 4);
307 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[0], 1), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
309 colPtrV[0] = 2 * colPtrM[0];
311 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 2), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
313 colPtrV[0] = 2 * colPtrM[0] + 1;
315 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 3), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
318 colPtrP[0] = (*coarseElementPDOFs)(coarseElement, 1);
322 PP->insertGlobalValues((*fineElementPDOFs)(fineElement[1], 1), colPtrP.
view(0, nnz), valPtrP.view(0, nnz));
325 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 0);
327 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 1);
329 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 4);
333 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[1], 4), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
335 colPtrV[0] = 2 * colPtrM[0];
336 colPtrV[1] = 2 * colPtrM[1];
337 colPtrV[2] = 2 * colPtrM[2];
339 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[1], 8), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
341 colPtrV[0] = 2 * colPtrM[0] + 1;
342 colPtrV[1] = 2 * colPtrM[1] + 1;
343 colPtrV[2] = 2 * colPtrM[2] + 1;
345 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[1], 9), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
348 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 1);
352 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[1], 1), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
354 colPtrV[0] = 2 * colPtrM[0];
356 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[1], 2), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
358 colPtrV[0] = 2 * colPtrM[0] + 1;
360 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[1], 3), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
363 if (coarseElement == 0) {
366 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 0);
370 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[0], 0), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
372 colPtrV[0] = 2 * colPtrM[0];
374 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 0), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
376 colPtrV[0] = 2 * colPtrM[0] + 1;
378 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 1), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
381 colPtrP[0] = (*coarseElementPDOFs)(coarseElement, 0);
385 PP->insertGlobalValues((*fineElementPDOFs)(fineElement[0], 0), colPtrP.
view(0, nnz), valPtrP.view(0, nnz));
391 if (coarseElement % (nCoarseElements) == 0) {
394 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 0);
396 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 3);
398 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 7);
402 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[0], 7), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
404 colPtrV[0] = 2 * colPtrM[0];
405 colPtrV[1] = 2 * colPtrM[1];
406 colPtrV[2] = 2 * colPtrM[2];
408 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 14), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
410 colPtrV[0] = 2 * colPtrM[0] + 1;
411 colPtrV[1] = 2 * colPtrM[1] + 1;
412 colPtrV[2] = 2 * colPtrM[2] + 1;
414 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 15), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
417 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 7);
421 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[0], 3), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
423 colPtrV[0] = 2 * colPtrM[0];
425 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 6), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
427 colPtrV[0] = 2 * colPtrM[0] + 1;
429 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 7), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
432 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 0);
434 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 3);
436 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 7);
440 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[2], 7), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
442 colPtrV[0] = 2 * colPtrM[0];
443 colPtrV[1] = 2 * colPtrM[1];
444 colPtrV[2] = 2 * colPtrM[2];
446 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[2], 14), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
448 colPtrV[0] = 2 * colPtrM[0] + 1;
449 colPtrV[1] = 2 * colPtrM[1] + 1;
450 colPtrV[2] = 2 * colPtrM[2] + 1;
452 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[2], 15), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
455 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 3);
459 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[2], 3), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
461 colPtrV[0] = 2 * colPtrM[0];
463 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[2], 6), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
465 colPtrV[0] = 2 * colPtrM[0] + 1;
467 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[2], 7), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
470 colPtrP[0] = (*coarseElementPDOFs)(coarseElement, 0);
472 colPtrP[1] = (*coarseElementPDOFs)(coarseElement, 3);
476 PP->insertGlobalValues((*fineElementPDOFs)(fineElement[0], 3), colPtrP.
view(0, nnz), valPtrP.view(0, nnz));
479 colPtrP[0] = (*coarseElementPDOFs)(coarseElement, 3);
483 PP->insertGlobalValues((*fineElementPDOFs)(fineElement[2], 3), colPtrP.
view(0, nnz), valPtrP.view(0, nnz));
489 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 8);
493 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[3], 0), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
495 colPtrV[0] = 2 * colPtrM[0];
497 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 0), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
499 colPtrV[0] = 2 * colPtrM[0] + 1;
501 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 1), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
504 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 5);
508 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[3], 1), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
510 colPtrV[0] = 2 * colPtrM[0];
512 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 2), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
514 colPtrV[0] = 2 * colPtrM[0] + 1;
516 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 3), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
519 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 2);
523 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[3], 2), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
525 colPtrV[0] = 2 * colPtrM[0];
527 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 4), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
529 colPtrV[0] = 2 * colPtrM[0] + 1;
531 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 5), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
534 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 6);
538 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[3], 3), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
540 colPtrV[0] = 2 * colPtrM[0];
542 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 6), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
544 colPtrV[0] = 2 * colPtrM[0] + 1;
546 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 7), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
549 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 4);
551 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 6);
553 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 8);
557 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[0], 5), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
559 colPtrV[0] = 2 * colPtrM[0];
560 colPtrV[1] = 2 * colPtrM[1];
561 colPtrV[2] = 2 * colPtrM[2];
563 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 10), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
565 colPtrV[0] = 2 * colPtrM[0] + 1;
566 colPtrV[1] = 2 * colPtrM[1] + 1;
567 colPtrV[2] = 2 * colPtrM[2] + 1;
569 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 11), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
572 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 5);
574 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 7);
576 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 8);
580 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[0], 6), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
582 colPtrV[0] = 2 * colPtrM[0];
583 colPtrV[1] = 2 * colPtrM[1];
584 colPtrV[2] = 2 * colPtrM[2];
586 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 12), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
588 colPtrV[0] = 2 * colPtrM[0] + 1;
589 colPtrV[1] = 2 * colPtrM[1] + 1;
590 colPtrV[2] = 2 * colPtrM[2] + 1;
592 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 13), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
595 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 1);
597 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 2);
599 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 5);
603 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[1], 5), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
605 colPtrV[0] = 2 * colPtrM[0];
606 colPtrV[1] = 2 * colPtrM[1];
607 colPtrV[2] = 2 * colPtrM[2];
609 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[1], 10), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
611 colPtrV[0] = 2 * colPtrM[0] + 1;
612 colPtrV[1] = 2 * colPtrM[1] + 1;
613 colPtrV[2] = 2 * colPtrM[2] + 1;
615 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[1], 11), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
618 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 5);
620 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 7);
622 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 8);
626 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[1], 6), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
628 colPtrV[0] = 2 * colPtrM[0];
629 colPtrV[1] = 2 * colPtrM[1];
630 colPtrV[2] = 2 * colPtrM[2];
632 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[1], 12), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
634 colPtrV[0] = 2 * colPtrM[0] + 1;
635 colPtrV[1] = 2 * colPtrM[1] + 1;
636 colPtrV[2] = 2 * colPtrM[2] + 1;
638 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[1], 13), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
641 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 4);
643 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 6);
645 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 8);
649 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[2], 5), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
651 colPtrV[0] = 2 * colPtrM[0];
652 colPtrV[1] = 2 * colPtrM[1];
653 colPtrV[2] = 2 * colPtrM[2];
655 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[2], 10), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
657 colPtrV[0] = 2 * colPtrM[0] + 1;
658 colPtrV[1] = 2 * colPtrM[1] + 1;
659 colPtrV[2] = 2 * colPtrM[2] + 1;
661 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[2], 11), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
664 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 2);
666 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 3);
668 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 6);
672 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[2], 6), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
674 colPtrV[0] = 2 * colPtrM[0];
675 colPtrV[1] = 2 * colPtrM[1];
676 colPtrV[2] = 2 * colPtrM[2];
678 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[2], 12), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
680 colPtrV[0] = 2 * colPtrM[0] + 1;
681 colPtrV[1] = 2 * colPtrM[1] + 1;
682 colPtrV[2] = 2 * colPtrM[2] + 1;
684 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[2], 13), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
687 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 1);
689 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 2);
691 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 5);
695 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[3], 5), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
697 colPtrV[0] = 2 * colPtrM[0];
698 colPtrV[1] = 2 * colPtrM[1];
699 colPtrV[2] = 2 * colPtrM[2];
701 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 10), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
703 colPtrV[0] = 2 * colPtrM[0] + 1;
704 colPtrV[1] = 2 * colPtrM[1] + 1;
705 colPtrV[2] = 2 * colPtrM[2] + 1;
707 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 11), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
710 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 2);
712 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 3);
714 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 6);
718 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[3], 6), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
720 colPtrV[0] = 2 * colPtrM[0];
721 colPtrV[1] = 2 * colPtrM[1];
722 colPtrV[2] = 2 * colPtrM[2];
724 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 12), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
726 colPtrV[0] = 2 * colPtrM[0] + 1;
727 colPtrV[1] = 2 * colPtrM[1] + 1;
728 colPtrV[2] = 2 * colPtrM[2] + 1;
730 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 13), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
733 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 0);
734 valPtrM[0] = 0.140625;
735 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 1);
736 valPtrM[1] = -0.046875;
737 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 2);
738 valPtrM[2] = 0.015625;
739 colPtrM[3] = (*coarseElementMDOFs)(coarseElement, 3);
740 valPtrM[3] = -0.046875;
741 colPtrM[4] = (*coarseElementMDOFs)(coarseElement, 4);
742 valPtrM[4] = 0.28125;
743 colPtrM[5] = (*coarseElementMDOFs)(coarseElement, 5);
744 valPtrM[5] = -0.09375;
745 colPtrM[6] = (*coarseElementMDOFs)(coarseElement, 6);
746 valPtrM[6] = -0.09375;
747 colPtrM[7] = (*coarseElementMDOFs)(coarseElement, 7);
748 valPtrM[7] = 0.28125;
749 colPtrM[8] = (*coarseElementMDOFs)(coarseElement, 8);
753 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[0], 8), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
755 colPtrV[0] = 2 * colPtrM[0];
756 colPtrV[1] = 2 * colPtrM[1];
757 colPtrV[2] = 2 * colPtrM[2];
758 colPtrV[3] = 2 * colPtrM[3];
759 colPtrV[4] = 2 * colPtrM[4];
760 colPtrV[5] = 2 * colPtrM[5];
761 colPtrV[6] = 2 * colPtrM[6];
762 colPtrV[7] = 2 * colPtrM[7];
763 colPtrV[8] = 2 * colPtrM[8];
765 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 16), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
767 colPtrV[0] = 2 * colPtrM[0] + 1;
768 colPtrV[1] = 2 * colPtrM[1] + 1;
769 colPtrV[2] = 2 * colPtrM[2] + 1;
770 colPtrV[3] = 2 * colPtrM[3] + 1;
771 colPtrV[4] = 2 * colPtrM[4] + 1;
772 colPtrV[5] = 2 * colPtrM[5] + 1;
773 colPtrV[6] = 2 * colPtrM[6] + 1;
774 colPtrV[7] = 2 * colPtrM[7] + 1;
775 colPtrV[8] = 2 * colPtrM[8] + 1;
777 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[0], 17), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
780 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 0);
781 valPtrM[0] = -0.046875;
782 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 1);
783 valPtrM[1] = 0.140625;
784 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 2);
785 valPtrM[2] = -0.046875;
786 colPtrM[3] = (*coarseElementMDOFs)(coarseElement, 3);
787 valPtrM[3] = 0.015625;
788 colPtrM[4] = (*coarseElementMDOFs)(coarseElement, 4);
789 valPtrM[4] = 0.28125;
790 colPtrM[5] = (*coarseElementMDOFs)(coarseElement, 5);
791 valPtrM[5] = 0.28125;
792 colPtrM[6] = (*coarseElementMDOFs)(coarseElement, 6);
793 valPtrM[6] = -0.09375;
794 colPtrM[7] = (*coarseElementMDOFs)(coarseElement, 7);
795 valPtrM[7] = -0.09375;
796 colPtrM[8] = (*coarseElementMDOFs)(coarseElement, 8);
800 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[1], 8), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
802 colPtrV[0] = 2 * colPtrM[0];
803 colPtrV[1] = 2 * colPtrM[1];
804 colPtrV[2] = 2 * colPtrM[2];
805 colPtrV[3] = 2 * colPtrM[3];
806 colPtrV[4] = 2 * colPtrM[4];
807 colPtrV[5] = 2 * colPtrM[5];
808 colPtrV[6] = 2 * colPtrM[6];
809 colPtrV[7] = 2 * colPtrM[7];
810 colPtrV[8] = 2 * colPtrM[8];
812 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[1], 16), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
814 colPtrV[0] = 2 * colPtrM[0] + 1;
815 colPtrV[1] = 2 * colPtrM[1] + 1;
816 colPtrV[2] = 2 * colPtrM[2] + 1;
817 colPtrV[3] = 2 * colPtrM[3] + 1;
818 colPtrV[4] = 2 * colPtrM[4] + 1;
819 colPtrV[5] = 2 * colPtrM[5] + 1;
820 colPtrV[6] = 2 * colPtrM[6] + 1;
821 colPtrV[7] = 2 * colPtrM[7] + 1;
822 colPtrV[8] = 2 * colPtrM[8] + 1;
824 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[1], 17), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
827 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 0);
828 valPtrM[0] = -0.046875;
829 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 1);
830 valPtrM[1] = 0.015625;
831 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 2);
832 valPtrM[2] = -0.046875;
833 colPtrM[3] = (*coarseElementMDOFs)(coarseElement, 3);
834 valPtrM[3] = 0.140625;
835 colPtrM[4] = (*coarseElementMDOFs)(coarseElement, 4);
836 valPtrM[4] = -0.09375;
837 colPtrM[5] = (*coarseElementMDOFs)(coarseElement, 5);
838 valPtrM[5] = -0.09375;
839 colPtrM[6] = (*coarseElementMDOFs)(coarseElement, 6);
840 valPtrM[6] = 0.28125;
841 colPtrM[7] = (*coarseElementMDOFs)(coarseElement, 7);
842 valPtrM[7] = 0.28125;
843 colPtrM[8] = (*coarseElementMDOFs)(coarseElement, 8);
847 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[2], 8), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
849 colPtrV[0] = 2 * colPtrM[0];
850 colPtrV[1] = 2 * colPtrM[1];
851 colPtrV[2] = 2 * colPtrM[2];
852 colPtrV[3] = 2 * colPtrM[3];
853 colPtrV[4] = 2 * colPtrM[4];
854 colPtrV[5] = 2 * colPtrM[5];
855 colPtrV[6] = 2 * colPtrM[6];
856 colPtrV[7] = 2 * colPtrM[7];
857 colPtrV[8] = 2 * colPtrM[8];
859 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[2], 16), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
861 colPtrV[0] = 2 * colPtrM[0] + 1;
862 colPtrV[1] = 2 * colPtrM[1] + 1;
863 colPtrV[2] = 2 * colPtrM[2] + 1;
864 colPtrV[3] = 2 * colPtrM[3] + 1;
865 colPtrV[4] = 2 * colPtrM[4] + 1;
866 colPtrV[5] = 2 * colPtrM[5] + 1;
867 colPtrV[6] = 2 * colPtrM[6] + 1;
868 colPtrV[7] = 2 * colPtrM[7] + 1;
869 colPtrV[8] = 2 * colPtrM[8] + 1;
871 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[2], 17), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
874 colPtrM[0] = (*coarseElementMDOFs)(coarseElement, 0);
875 valPtrM[0] = 0.015625;
876 colPtrM[1] = (*coarseElementMDOFs)(coarseElement, 1);
877 valPtrM[1] = -0.046875;
878 colPtrM[2] = (*coarseElementMDOFs)(coarseElement, 2);
879 valPtrM[2] = 0.140625;
880 colPtrM[3] = (*coarseElementMDOFs)(coarseElement, 3);
881 valPtrM[3] = -0.046875;
882 colPtrM[4] = (*coarseElementMDOFs)(coarseElement, 4);
883 valPtrM[4] = -0.09375;
884 colPtrM[5] = (*coarseElementMDOFs)(coarseElement, 5);
885 valPtrM[5] = 0.28125;
886 colPtrM[6] = (*coarseElementMDOFs)(coarseElement, 6);
887 valPtrM[6] = 0.28125;
888 colPtrM[7] = (*coarseElementMDOFs)(coarseElement, 7);
889 valPtrM[7] = -0.09375;
890 colPtrM[8] = (*coarseElementMDOFs)(coarseElement, 8);
894 PM->insertGlobalValues((*fineElementMDOFs)(fineElement[3], 8), colPtrM.
view(0, nnz), valPtrM.view(0, nnz));
896 colPtrV[0] = 2 * colPtrM[0];
897 colPtrV[1] = 2 * colPtrM[1];
898 colPtrV[2] = 2 * colPtrM[2];
899 colPtrV[3] = 2 * colPtrM[3];
900 colPtrV[4] = 2 * colPtrM[4];
901 colPtrV[5] = 2 * colPtrM[5];
902 colPtrV[6] = 2 * colPtrM[6];
903 colPtrV[7] = 2 * colPtrM[7];
904 colPtrV[8] = 2 * colPtrM[8];
906 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 16), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
908 colPtrV[0] = 2 * colPtrM[0] + 1;
909 colPtrV[1] = 2 * colPtrM[1] + 1;
910 colPtrV[2] = 2 * colPtrM[2] + 1;
911 colPtrV[3] = 2 * colPtrM[3] + 1;
912 colPtrV[4] = 2 * colPtrM[4] + 1;
913 colPtrV[5] = 2 * colPtrM[5] + 1;
914 colPtrV[6] = 2 * colPtrM[6] + 1;
915 colPtrV[7] = 2 * colPtrM[7] + 1;
916 colPtrV[8] = 2 * colPtrM[8] + 1;
918 PV->insertGlobalValues((*fineElementVDOFs)(fineElement[3], 17), colPtrV.
view(0, nnz), valPtrM.view(0, nnz));
921 colPtrP[0] = (*coarseElementPDOFs)(coarseElement, 0);
923 colPtrP[1] = (*coarseElementPDOFs)(coarseElement, 1);
925 colPtrP[2] = (*coarseElementPDOFs)(coarseElement, 2);
927 colPtrP[3] = (*coarseElementPDOFs)(coarseElement, 3);
931 PP->insertGlobalValues((*fineElementPDOFs)(fineElement[0], 2), colPtrP.
view(0, nnz), valPtrP.view(0, nnz));
934 colPtrP[0] = (*coarseElementPDOFs)(coarseElement, 1);
936 colPtrP[1] = (*coarseElementPDOFs)(coarseElement, 2);
940 PP->insertGlobalValues((*fineElementPDOFs)(fineElement[1], 2), colPtrP.
view(0, nnz), valPtrP.view(0, nnz));
943 colPtrP[0] = (*coarseElementPDOFs)(coarseElement, 2);
947 PP->insertGlobalValues((*fineElementPDOFs)(fineElement[3], 2), colPtrP.
view(0, nnz), valPtrP.view(0, nnz));
950 colPtrP[0] = (*coarseElementPDOFs)(coarseElement, 2);
952 colPtrP[1] = (*coarseElementPDOFs)(coarseElement, 3);
956 PP->insertGlobalValues((*fineElementPDOFs)(fineElement[3], 3), colPtrP.
view(0, nnz), valPtrP.view(0, nnz));
959 if ((coarseElement + 1) % (nCoarseElements) == 0)
961 fineElement[0] = fineElement[3] + 1;
962 fineElement[1] = fineElement[0] + 1;
963 fineElement[2] = fineElement[0] + nFineElements;
964 fineElement[3] = fineElement[2] + 1;
966 fineElement[0] = fineElement[1] + 1;
967 fineElement[1] = fineElement[0] + 1;
968 fineElement[2] = fineElement[3] + 1;
969 fineElement[3] = fineElement[2] + 1;
974 for (
GO VRow = 0; VRow < fNV; VRow++) {
978 PV->getGlobalRowView(VRow, colPtr, valPtr);
981 P->insertGlobalValues(VRow, colPtr, valPtr);
985 for (
GO PRow = 0; PRow < fNP; PRow++) {
990 PP->getGlobalRowView(PRow, colPtr, valPtr);
993 for (
LO jj = 0; jj < colPtr.
size(); jj++) {
994 newColPtr[jj] += colPtr[jj];
998 P->insertGlobalValues(PRow + fNV, newColPtr.view(0, colPtr.
size()), valPtr);
1002 for (
GO MRow = 0; MRow < fNM; MRow++) {
1007 PM->getGlobalRowView(MRow, colPtr, valPtr);
1010 for (
LO jj = 0; jj < colPtr.
size(); jj++) {
1011 newColPtr[jj] += colPtr[jj];
1015 P->insertGlobalValues(MRow + fNV + fNP, newColPtr.view(0, colPtr.
size()), valPtr);
1019 PV->fillComplete(colMapforPV, rowMapforPV);
1020 PP->fillComplete(colMapforPP, rowMapforPP);
1021 PM->fillComplete(colMapforPM, rowMapforPM);
1025 Set(coarseLevel,
"PV", PV);
1026 Set(coarseLevel,
"PP", PP);
1027 Set(coarseLevel,
"PM", PM);
1028 Set(coarseLevel,
"P", P);
1030 Set(coarseLevel,
"NV", nV);
1031 Set(coarseLevel,
"NP", nP);
1032 Set(coarseLevel,
"NM", nM);
1038 #define MUELU_GEOINTERPFACTORY_SHORT
1039 #endif // MUELU_GEOINTERPFACTORY_DEF_HPP
T & Get(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Get data without decrementing associated storage counter (i.e., read-only access). Usage: Level->Get< RCP<Matrix> >("A", factory) if factory == NULL => use default factory.
void BuildP(Level &fineLevel, Level &coarseLevel) const
Abstract Build method.
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Specifies the data that this class needs, and the factories that generate that data.
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Class that holds all level-specific information.
virtual ~GeoInterpFactory()
Destructor.
const RCP< const FactoryManagerBase > GetFactoryManager()
returns the current factory manager
Description of what is happening (more verbose)
GeoInterpFactory()
Constructor.
static Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > createUniformContigMap(UnderlyingLib lib, global_size_t numElements, const Teuchos::RCP< const Teuchos::Comm< int >> &comm)
ArrayView< T > view(size_type lowerOffset, size_type size) const