48 #ifndef _ZOLTAN2_SPHYNXALGORITHM_HPP_
49 #define _ZOLTAN2_SPHYNXALGORITHM_HPP_
75 #include "Zoltan2Sphynx_config.h"
83 #include "AnasaziLOBPCGSolMgr.hpp"
84 #include "AnasaziBlockDavidsonSolMgr.hpp"
85 #include "AnasaziGeneralizedDavidsonSolMgr.hpp"
86 #include "AnasaziBlockKrylovSchurSolMgr.hpp"
87 #include "AnasaziRandomizedSolMgr.hpp"
88 #include "AnasaziBasicEigenproblem.hpp"
89 #include "AnasaziTpetraAdapter.hpp"
91 #include "BelosLinearProblem.hpp"
92 #include "BelosTpetraOperator.hpp"
94 #include "Ifpack2_Factory.hpp"
96 #include "Teuchos_TimeMonitor.hpp"
98 #ifdef HAVE_ZOLTAN2SPHYNX_MUELU
99 #include "MueLu_CreateTpetraPreconditioner.hpp"
104 template <
typename Adapter>
118 using graph_t = Tpetra::CrsGraph<lno_t, gno_t, node_t>;
119 using matrix_t = Tpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t>;
120 using mvector_t = Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t>;
121 using op_t = Tpetra::Operator<scalar_t, lno_t, gno_t, node_t>;
124 typedef Anasazi::MultiVecTraits<scalar_t,mvector_t>
MVT;
131 Sphynx(
const RCP<const Environment> &env,
132 const RCP<Teuchos::ParameterList> ¶ms,
133 const RCP<Teuchos::ParameterList> &sphynxParams,
134 const RCP<
const Comm<int> > &comm,
136 env_(env), params_(params), sphynxParams_(sphynxParams), comm_(comm), adapter_(adapter)
140 const Teuchos::ParameterEntry *pe = params_->getEntryPtr(
"num_global_parts");
141 numGlobalParts_ = pe->getValue<
int>(&numGlobalParts_);
142 if(numGlobalParts_ > 1){
144 Teuchos::TimeMonitor t(*Teuchos::TimeMonitor::getNewTimer(
"Sphynx::Laplacian"));
147 pe = sphynxParams_->getEntryPtr(
"sphynx_verbosity");
149 verbosity_ = pe->getValue<
int>(&verbosity_);
152 pe = sphynxParams_->getEntryPtr(
"sphynx_skip_preprocessing");
154 skipPrep_ = pe->getValue<
bool>(&skipPrep_);
159 graph_ = adapter_->getUserGraph();
161 throw std::runtime_error(
"\nSphynx Error: Preprocessing has not been implemented yet.\n");
187 template<
typename problem_t>
190 template<
typename problem_t>
193 template<
typename problem_t>
196 template<
typename problem_t>
204 std::vector<int> strides);
207 std::vector<const weight_t *>
weights,
208 std::vector<int> strides,
224 auto rowOffsets = graph_->getLocalGraphDevice().row_map;
225 auto rowOffsets_h = Kokkos::create_mirror_view(rowOffsets);
226 Kokkos::deep_copy(rowOffsets_h, rowOffsets);
229 const size_t numGlobalEntries = graph_->getGlobalNumEntries();
230 const size_t numLocalRows = graph_->getLocalNumRows();
231 const size_t numGlobalRows = graph_->getGlobalNumRows();
235 for(
size_t i = 0; i < numLocalRows; i++){
236 if(rowOffsets_h(i+1) - rowOffsets_h(i) - 1 > localMax)
237 localMax = rowOffsets_h(i+1) - rowOffsets_h(i) - 1;
241 size_t globalMax = 0;
242 Teuchos::reduceAll<int, size_t> (*comm_, Teuchos::REDUCE_MAX, 1, &localMax, &globalMax);
244 double avg =
static_cast<double>(numGlobalEntries-numGlobalRows)/numGlobalRows;
245 double maxOverAvg =
static_cast<double>(globalMax)/avg;
249 if(maxOverAvg > 10) {
255 if(comm_->getRank() == 0) {
256 std::cout <<
"Regularity of Graph ----------------" << std::endl;
257 std::cout <<
" Maximum degree: " << globalMax << std::endl;
258 std::cout <<
" Average degree: " << avg << std::endl;
259 std::cout <<
" Max/Avg: " << globalMax/avg << std::endl;
260 std::cout <<
" Regular graph?: " << !irregular_ << std::endl;
280 precType_ =
"jacobi";
281 #ifdef HAVE_ZOLTAN2SPHYNX_MUELU
286 const Teuchos::ParameterEntry *pe = sphynxParams_->getEntryPtr(
"sphynx_preconditioner_type");
288 precType_ = pe->getValue<std::string>(&precType_);
289 if(precType_ !=
"muelu" && precType_ !=
"jacobi" && precType_ !=
"polynomial")
290 throw std::runtime_error(
"\nSphynx Error: " + precType_ +
" is not recognized as a preconditioner.\n"
291 +
" Possible values: muelu (if enabled), jacobi, and polynomial\n");
294 solverType_ = sphynxParams_->get(
"sphynx_eigensolver",
"LOBPCG");
295 TEUCHOS_TEST_FOR_EXCEPTION(!(solverType_ ==
"LOBPCG" || solverType_ ==
"randomized" || solverType_ ==
"BlockDavidson" || solverType_ ==
"GeneralizedDavidson" || solverType_ ==
"BlockKrylovSchur" ),
296 std::invalid_argument,
"Sphynx: sphynx_eigensolver must be set to LOBPCG, randomized, BlockDavidson, GeneralizedDavidson, or BlockKrylovSchur.");
301 if(precType_ ==
"polynomial")
308 pe = sphynxParams_->getEntryPtr(
"sphynx_problem_type");
310 std::string probType =
"";
311 probType = pe->getValue<std::string>(&probType);
313 if(probType ==
"combinatorial")
315 else if(probType ==
"generalized")
317 else if(probType ==
"normalized")
320 throw std::runtime_error(
"\nSphynx Error: " + probType +
" is not recognized as a problem type.\n"
321 +
" Possible values: combinatorial, generalized, and normalized.\n");
327 if(!irregular_ && (precType_ ==
"jacobi" || precType_ ==
"polynomial"))
332 pe = sphynxParams_->getEntryPtr(
"sphynx_tolerance");
334 tolerance_ = pe->getValue<
scalar_t>(&tolerance_);
343 pe = sphynxParams_->getEntryPtr(
"sphynx_initial_guess");
345 std::string initialGuess =
"";
346 initialGuess = pe->getValue<std::string>(&initialGuess);
348 if(initialGuess ==
"random")
350 else if(initialGuess ==
"constants")
353 throw std::runtime_error(
"\nSphynx Error: " + initialGuess +
" is not recognized as an option for initial guess.\n"
354 +
" Possible values: random and constants.\n");
366 if(solverType_ ==
"randomized")
380 auto rowOffsets = graph_->getLocalGraphHost().row_map;
383 Teuchos::RCP<matrix_t> degMat(
new matrix_t (graph_->getRowMap(),
389 lno_t numRows =
static_cast<lno_t>(graph_->getLocalNumRows());
392 for (
lno_t i = 0; i < numRows; ++i) {
393 val[0] = rowOffsets(i+1) - rowOffsets(i) - 1;
395 degMat->insertLocalValues(i, 1, val, ind);
400 degMat->fillComplete(graph_->getDomainMap(), graph_->getRangeMap());
412 using range_policy = Kokkos::RangePolicy<
413 typename node_t::device_type::execution_space, Kokkos::IndexType<lno_t>>;
414 using values_view_t = Kokkos::View<scalar_t*, typename node_t::device_type>;
415 using offset_view_t = Kokkos::View<size_t*, typename node_t::device_type>;
417 const size_t numEnt = graph_->getLocalNumEntries();
418 const size_t numRows = graph_->getLocalNumRows();
421 values_view_t newVal (Kokkos::view_alloc(
"CombLapl::val", Kokkos::WithoutInitializing), numEnt);
422 Kokkos::deep_copy(newVal, -1);
425 offset_view_t diagOffsets(Kokkos::view_alloc(
"Diag Offsets", Kokkos::WithoutInitializing), numRows);
426 graph_->getLocalDiagOffsets(diagOffsets);
429 auto rowOffsets = graph_->getLocalGraphDevice().row_map;
432 Kokkos::parallel_for(
"Combinatorial Laplacian", range_policy(0, numRows),
433 KOKKOS_LAMBDA(
const lno_t i){
434 newVal(rowOffsets(i) + diagOffsets(i)) = rowOffsets(i+1) - rowOffsets(i) - 1;
440 Teuchos::RCP<matrix_t> laplacian (
new matrix_t(graph_, newVal));
441 laplacian->fillComplete (graph_->getDomainMap(), graph_->getRangeMap());
460 using range_policy = Kokkos::RangePolicy<
461 typename node_t::device_type::execution_space, Kokkos::IndexType<lno_t>>;
462 using values_view_t = Kokkos::View<scalar_t*, typename node_t::device_type>;
463 using offset_view_t = Kokkos::View<size_t*, typename node_t::device_type>;
464 using vector_t = Tpetra::Vector<scalar_t, lno_t, gno_t, node_t>;
465 using dual_view_t =
typename vector_t::dual_view_type;
466 using KAT = Kokkos::ArithTraits<scalar_t>;
468 const size_t numEnt = graph_->getLocalNumEntries();
469 const size_t numRows = graph_->getLocalNumRows();
472 values_view_t newVal (Kokkos::view_alloc(
"NormLapl::val", Kokkos::WithoutInitializing), numEnt);
474 Kokkos::deep_copy(newVal, 1);
477 Kokkos::deep_copy(newVal, -1);
481 dual_view_t dv (
"MV::DualView", numRows, 1);
482 auto deginvsqrt = dv.d_view;
485 offset_view_t diagOffsets(Kokkos::view_alloc(
"Diag Offsets", Kokkos::WithoutInitializing), numRows);
486 graph_->getLocalDiagOffsets(diagOffsets);
489 auto rowOffsets = graph_->getLocalGraphDevice().row_map;
493 Kokkos::parallel_for(
"Combinatorial Laplacian", range_policy(0, numRows),
494 KOKKOS_LAMBDA(
const lno_t i){
495 newVal(rowOffsets(i) + diagOffsets(i)) = rowOffsets(i+1) - rowOffsets(i) - 1;
496 deginvsqrt(i,0) = 1.0/KAT::sqrt(rowOffsets(i+1) - rowOffsets(i) - 1);
502 Teuchos::RCP<matrix_t> laplacian (
new matrix_t(graph_, newVal));
503 laplacian->fillComplete (graph_->getDomainMap(), graph_->getRangeMap());
506 vector_t degInvSqrt(graph_->getRowMap(), dv);
509 laplacian->leftScale(degInvSqrt);
510 laplacian->rightScale(degInvSqrt);
522 const Teuchos::RCP<const Environment> env_;
523 Teuchos::RCP<Teuchos::ParameterList> params_;
524 Teuchos::RCP<Teuchos::ParameterList> sphynxParams_;
525 const Teuchos::RCP<const Teuchos::Comm<int> > comm_;
526 const Teuchos::RCP<const Adapter> adapter_;
530 Teuchos::RCP<const graph_t> graph_;
531 Teuchos::RCP<matrix_t> laplacian_;
532 Teuchos::RCP<const matrix_t> degMatrix_;
533 Teuchos::RCP<mvector_t> eigenVectors_;
536 std::string precType_;
537 std::string solverType_;
554 template <
typename Adapter>
557 eigenVectors_ = userEvects;
566 template <
typename Adapter>
570 if(numGlobalParts_ == 1) {
572 size_t numRows =adapter_->getUserGraph()->getLocalNumRows();
573 Teuchos::ArrayRCP<part_t> parts(numRows,0);
574 solution->setParts(parts);
581 int numEigenVectors = (int) log2(numGlobalParts_)+1;
584 if(eigenVectors_ == Teuchos::null){
589 if(computedNumEv <= 1 && (solverType_ ==
"LOBPCG" || solverType_ ==
"GeneralizedDavidson" || solverType_ ==
"BlockDavidson" || solverType_ ==
"BlockKrylovSchur")) {
591 std::runtime_error(
"\nAnasazi Error: LOBPCGSolMgr::solve() returned unconverged.\n"
592 "Sphynx Error: LOBPCG could not compute any eigenvectors.\n"
593 " Increase either max iters or tolerance.\n");
598 computedNumEv = (int) eigenVectors_->getNumVectors();
599 if(computedNumEv <= numEigenVectors) {
601 std::runtime_error(
"\nSphynx Error: Number of eigenvectors given by user\n"
602 " is less than number of Eigenvectors needed for partition." );
612 std::vector<const weight_t *>
weights;
613 std::vector<int> wstrides;
625 template <
typename Adapter>
629 Teuchos::TimeMonitor t(*Teuchos::TimeMonitor::getNewTimer(
"Sphynx::Anasazi"));
633 std::string which = (solverType_ ==
"randomized" ?
"LM" :
"SR");
634 std::string ortho =
"SVQB";
635 bool relConvTol =
false;
636 int maxIterations = sphynxParams_->get(
"sphynx_max_iterations",1000);
637 int blockSize = sphynxParams_->get(
"sphynx_block_size",numEigenVectors);
638 int orthoFreq = sphynxParams_->get(
"sphynx_ortho_freq", 0);
639 int resFreq = sphynxParams_->get(
"sphynx_res_freq", 0);
640 bool isHermitian =
true;
641 bool relLockTol =
false;
643 bool useFullOrtho = sphynxParams_->get(
"sphynx_use_full_ortho",
true);
650 int anasaziVerbosity = Anasazi::Errors + Anasazi::Warnings;
652 anasaziVerbosity += Anasazi::FinalSummary + Anasazi::TimingDetails;
654 anasaziVerbosity += Anasazi::IterationDetails;
656 anasaziVerbosity += Anasazi::StatusTestDetails + Anasazi::OrthoDetails
660 Teuchos::ParameterList anasaziParams;
661 anasaziParams.set(
"Verbosity", anasaziVerbosity);
662 anasaziParams.set(
"Which", which);
663 anasaziParams.set(
"Convergence Tolerance", tolerance_);
664 anasaziParams.set(
"Maximum Iterations", maxIterations);
665 anasaziParams.set(
"Block Size", blockSize);
666 anasaziParams.set(
"Relative Convergence Tolerance", relConvTol);
667 anasaziParams.set(
"Orthogonalization", ortho);
668 anasaziParams.set(
"Use Locking", lock);
669 anasaziParams.set(
"Relative Locking Tolerance", relLockTol);
670 anasaziParams.set(
"Full Ortho", useFullOrtho);
671 anasaziParams.set(
"Orthogonalization Frequency", orthoFreq);
672 anasaziParams.set(
"Residual Frequency", resFreq);
674 if(solverType_ ==
"GeneralizedDavidson" || solverType_ ==
"BlockKrylovSchur" || solverType_ ==
"BlockDavidson"){
675 anasaziParams.set(
"Num Blocks", maxIterations+1 );
676 anasaziParams.set(
"Maximum Restarts", 0 );
677 anasaziParams.set(
"Maximum Subspace Dimension", (maxIterations+1)*blockSize );
681 auto map = laplacian_->getRangeMap();
682 Teuchos::RCP<mvector_t> ivec(
new mvector_t(map, numEigenVectors));
686 Anasazi::MultiVecTraits<scalar_t, mvector_t>::MvRandom(*ivec);
687 ivec->getVectorNonConst(0)->putScalar(1.);
693 ivec->getVectorNonConst(0)->putScalar(1.);
694 for (
int j = 1; j < numEigenVectors; j++)
695 ivec->getVectorNonConst(j)->putScalar(0.);
697 gno_t blkSize = map->getGlobalNumElements() / numEigenVectors;
698 TEUCHOS_TEST_FOR_EXCEPTION(blkSize <= 0, std::runtime_error,
"Blocksize too small for \"constants\" initial guess. Try \"random\".");
700 for (
size_t lid = 0; lid < ivec->getLocalLength(); lid++) {
701 gno_t gid = map->getGlobalElement(lid);
702 for (
int j = 1; j < numEigenVectors; j++){
703 if (((j-1)*blkSize <= gid) && (j*blkSize > gid))
704 ivec->replaceLocalValue(lid,j,1.);
710 using problem_t = Anasazi::BasicEigenproblem<scalar_t, mvector_t, op_t>;
711 Teuchos::RCP<problem_t> problem (
new problem_t(laplacian_, ivec));
712 problem->setHermitian(isHermitian);
713 problem->setNEV(numEigenVectors);
715 if(solverType_ !=
"randomized"){
722 bool boolret = problem->setProblem();
723 if (boolret !=
true) {
724 throw std::runtime_error(
"\nAnasazi::BasicEigenproblem::setProblem() returned with error.\n");
727 Teuchos::RCP<Anasazi::SolverManager<scalar_t, mvector_t, op_t>> solver;
729 if(solverType_ ==
"LOBPCG"){
730 solver = Teuchos::rcp(
new Anasazi::LOBPCGSolMgr<scalar_t, mvector_t, op_t>(problem, anasaziParams));
732 else if(solverType_ ==
"BlockDavidson"){
733 solver = Teuchos::rcp(
new Anasazi::BlockDavidsonSolMgr<scalar_t, mvector_t, op_t>(problem, anasaziParams));
735 else if(solverType_ ==
"GeneralizedDavidson"){
736 solver = Teuchos::rcp(
new Anasazi::GeneralizedDavidsonSolMgr<scalar_t, mvector_t, op_t>(problem, anasaziParams));
738 else if(solverType_ ==
"BlockKrylovSchur"){
739 solver = Teuchos::rcp(
new Anasazi::BlockKrylovSchurSolMgr<scalar_t, mvector_t, op_t>(problem, anasaziParams));
742 solver = Teuchos::rcp(
new Anasazi::Experimental::RandomizedSolMgr<scalar_t, mvector_t, op_t>(problem, anasaziParams));
746 if (verbosity_ > 0 && comm_->getRank() == 0)
747 anasaziParams.print(std::cout);
750 Anasazi::ReturnType returnCode = solver->solve();
753 iter = solver->getNumIters();
757 using solution_t = Anasazi::Eigensolution<scalar_t, mvector_t>;
758 solution_t sol = problem->getSolution();
759 eigenVectors_ = sol.Evecs;
760 int numev = sol.numVecs;
761 std::vector<Anasazi::Value<double>> eigenValues_ = sol.Evals;
764 if (verbosity_ > 0 && comm_->getRank() == 0) {
765 std::cout << std::endl;
766 std::cout <<
"ANASAZI SUMMARY" << std::endl;
767 std::cout <<
"Failed to converge: " << numfailed << std::endl;
768 std::cout <<
"No of iterations : " << iter << std::endl;
769 std::cout <<
"Solve time: " << std::endl;
770 std::cout <<
"No of comp. vecs. : " << numev << std::endl;
773 std::cout <<
"Solver type: " << solverType_ << std::endl;
775 if(solverType_ ==
"randomized") {
776 std::vector<double> normR(numev);
780 Teuchos::SerialDenseMatrix<int,double> T (numev, numev);
782 for (
int i = 0; i < numev; ++i) T(i,i) = eigenValues_[i].realpart;
783 laplacian_->apply(*eigenVectors_, Aevec);
784 MVT::MvTimesMatAddMv(-1.0, *eigenVectors_, T, 1.0, Aevec);
785 MVT::MvNorm(Aevec, normR);
788 if(comm_->getRank() == 0 && verbosity_ > 0) {
789 std::cout << std::endl <<
"Solver manager returned "
790 << (returnCode == Anasazi::Converged ?
"converged." :
"unconverged.")
791 << std::endl << std::endl
792 << std::setw(16) <<
"Eigenvalue"
793 << std::setw(18) <<
"Direct Residual"
795 <<
"------------------------------------------------------" << std::endl;
797 for (
int i = 0; i < numev; ++i) {
798 std::cout << std::setw(16) << 2.0-eigenValues_[i].realpart
799 << std::setw(18) << normR[i] / eigenValues_[i].realpart
802 std::cout <<
"------------------------------------------------------" << std::endl;
811 template <
typename Adapter>
812 template <
typename problem_t>
815 if(comm_->getRank() == 0) std::cout <<
"precType_ is: " << precType_ << std::endl;
817 if(precType_ ==
"muelu") {
820 else if(precType_ ==
"polynomial") {
823 else if(precType_ ==
"jacobi") {
830 template <
typename Adapter>
831 template <
typename problem_t>
834 #ifdef HAVE_ZOLTAN2SPHYNX_MUELU
835 Teuchos::ParameterList paramList;
837 paramList.set(
"verbosity",
"none");
838 else if(verbosity_ == 1)
839 paramList.set(
"verbosity",
"low");
840 else if(verbosity_ == 2)
841 paramList.set(
"verbosity",
"medium");
842 else if(verbosity_ == 3)
843 paramList.set(
"verbosity",
"high");
844 else if(verbosity_ >= 4)
845 paramList.set(
"verbosity",
"extreme");
847 paramList.set(
"smoother: type",
"CHEBYSHEV");
848 Teuchos::ParameterList smootherParamList;
849 smootherParamList.set(
"chebyshev: degree", 3);
850 smootherParamList.set(
"chebyshev: ratio eigenvalue", 7.0);
851 smootherParamList.set(
"chebyshev: eigenvalue max iterations", irregular_ ? 100 : 10);
852 paramList.set(
"smoother: params", smootherParamList);
853 paramList.set(
"use kokkos refactor",
true);
854 paramList.set(
"transpose: use implicit",
true);
858 paramList.set(
"multigrid algorithm",
"unsmoothed");
860 paramList.set(
"coarse: type",
"CHEBYSHEV");
861 Teuchos::ParameterList coarseParamList;
862 coarseParamList.set(
"chebyshev: degree", 3);
863 coarseParamList.set(
"chebyshev: ratio eigenvalue", 7.0);
864 paramList.set(
"coarse: params", coarseParamList);
866 paramList.set(
"max levels", 5);
867 paramList.set(
"aggregation: drop tol", 0.40);
870 using prec_t = MueLu::TpetraOperator<scalar_t, lno_t, gno_t, node_t>;
871 Teuchos::RCP<prec_t> prec = MueLu::CreateTpetraPreconditioner<
874 problem->setPrec(prec);
877 throw std::runtime_error(
"\nSphynx Error: MueLu requested but not compiled into Trilinos.\n");
882 template <
typename Adapter>
883 template <
typename problem_t>
886 int verbosity2 = Belos::Errors;
888 verbosity2 += Belos::Warnings;
889 else if(verbosity_ == 2)
890 verbosity2 += Belos::Warnings + Belos::FinalSummary;
891 else if(verbosity_ == 3)
892 verbosity2 += Belos::Warnings + Belos::FinalSummary + Belos::TimingDetails;
893 else if(verbosity_ >= 4)
894 verbosity2 += Belos::Warnings + Belos::FinalSummary + Belos::TimingDetails
895 + Belos::StatusTestDetails;
897 Teuchos::ParameterList paramList;
898 paramList.set(
"Polynomial Type",
"Roots");
899 paramList.set(
"Orthogonalization",
"ICGS");
900 paramList.set(
"Maximum Degree", laplacian_->getGlobalNumRows() > 100 ? 25 : 5);
901 paramList.set(
"Polynomial Tolerance", 1.0e-6 );
902 paramList.set(
"Verbosity", verbosity2 );
903 paramList.set(
"Random RHS",
false );
904 paramList.set(
"Outer Solver",
"");
905 paramList.set(
"Timer Label",
"Belos Polynomial Solve" );
908 using lproblem_t = Belos::LinearProblem<scalar_t, mvector_t, op_t>;
909 Teuchos::RCP<lproblem_t> innerPolyProblem(
new lproblem_t());
910 innerPolyProblem->setOperator(laplacian_);
912 using btop_t = Belos::TpetraOperator<scalar_t, lno_t, gno_t, node_t>;
913 Teuchos::RCP<btop_t> polySolver(
new btop_t(innerPolyProblem,
914 Teuchos::rcpFromRef(paramList),
916 problem->setPrec(polySolver);
920 template <
typename Adapter>
921 template <
typename problem_t>
925 Teuchos::RCP<Ifpack2::Preconditioner<scalar_t, lno_t, gno_t, node_t>> prec;
926 std::string precType =
"RELAXATION";
928 prec = Ifpack2::Factory::create<matrix_t> (precType, laplacian_);
930 Teuchos::ParameterList precParams;
931 precParams.set(
"relaxation: type",
"Jacobi");
932 precParams.set(
"relaxation: fix tiny diagonal entries",
true);
933 precParams.set(
"relaxation: min diagonal value", 1.0e-8);
935 prec->setParameters(precParams);
939 problem->setPrec(prec);
945 template <
typename Adapter>
951 Teuchos::Array<size_t> columns (computedNumEv-1);
952 for (
int j = 0; j < computedNumEv-1; ++j) {
955 coordinates = eigenVectors->subCopy (columns());
956 coordinates->setCopyOrView (Teuchos::View);
964 template <
typename Adapter>
966 std::vector<int> strides)
969 int numWeights = adapter_->getNumWeightsPerVertex();
970 int numConstraints = numWeights > 0 ? numWeights : 1;
972 size_t myNumVertices = adapter_->getLocalNumVertices();
974 for(
int j = 0; j < numConstraints; j++)
975 weights[j] =
new weight_t[myNumVertices];
979 const gno_t *columnId;
982 if(numWeights == 0) {
985 adapter_->getEdgesView(offset, columnId);
986 for (
size_t i = 0; i < myNumVertices; i++)
987 weights[0][i] = offset[i+1] - offset[i] - 1;
989 vecweights.push_back(weights[0]);
990 strides.push_back(1);
995 for(
int j = 0; j < numConstraints; j++) {
997 if(adapter_->useDegreeAsVertexWeight(j)) {
999 adapter_->getEdgesView(offset, columnId);
1000 for (
size_t i = 0; i < myNumVertices; i++)
1001 weights[j][i] = offset[i+1] - offset[i];
1006 adapter_->getVertexWeightsView(wgt, stride, j);
1008 for (
size_t i = 0; i < myNumVertices; i++)
1009 weights[j][i] = wgt[i];
1012 vecweights.push_back(weights[j]);
1013 strides.push_back(1);
1023 template <
typename Adapter>
1025 std::vector<const weight_t *>
weights,
1026 std::vector<int> strides,
1030 Teuchos::TimeMonitor t(*Teuchos::TimeMonitor::getNewTimer(
"Sphynx::MJ"));
1037 size_t myNumVertices = coordinates->getLocalLength();
1040 Teuchos::RCP<mvector_adapter_t> adapcoordinates(
new mvector_adapter_t(coordinates,
1043 Teuchos::RCP<const base_adapter_t> baseAdapter =
1044 Teuchos::rcp(dynamic_cast<const base_adapter_t *>(adapcoordinates.get()),
false);
1050 Teuchos::RCP<const Comm<int>> comm2 = comm_;
1051 Teuchos::RCP<mj_t> mj(
new mj_t(env_, comm2, baseAdapter));
1054 Teuchos::RCP<solution_t> vectorsolution(
new solution_t(env_, comm2, 1, mj));
1055 mj->partition(vectorsolution);
1058 Teuchos::ArrayRCP<part_t> parts(myNumVertices);
1059 for(
size_t i = 0; i < myNumVertices; i++) parts[i] = vectorsolution->getPartListView()[i];
1060 solution->setParts(parts);
Zoltan2::BaseAdapter< userTypes_t > base_adapter_t
void partition(const Teuchos::RCP< PartitioningSolution< Adapter > > &solution)
Sphynx(const RCP< const Environment > &env, const RCP< Teuchos::ParameterList > ¶ms, const RCP< Teuchos::ParameterList > &sphynxParams, const RCP< const Comm< int > > &comm, const RCP< const XpetraCrsGraphAdapter< graph_t > > &adapter)
Anasazi::MultiVecTraits< scalar_t, mvector_t > MVT
void computeWeights(std::vector< const weight_t * > vecweights, std::vector< int > strides)
void setMueLuPreconditioner(Teuchos::RCP< problem_t > &problem)
std::bitset< NUM_MODEL_FLAGS > modelFlag_t
void eigenvecsToCoords(Teuchos::RCP< mvector_t > &eigenVectors, int computedNumEv, Teuchos::RCP< mvector_t > &coordinates)
map_t::global_ordinal_type gno_t
void computeDegreeMatrix()
typename Adapter::node_t node_t
void setPolynomialPreconditioner(Teuchos::RCP< problem_t > &problem)
Provides access for Zoltan2 to Xpetra::CrsGraph data.
Tpetra::MultiVector< scalar_t, lno_t, gno_t, node_t > mvector_t
void setJacobiPreconditioner(Teuchos::RCP< problem_t > &problem)
Teuchos::RCP< matrix_t > computeCombinatorialLaplacian()
void setUserEigenvectors(const Teuchos::RCP< mvector_t > &userEvects)
Defines the XpetraMultiVectorAdapter.
Defines XpetraCrsGraphAdapter class.
SparseMatrixAdapter_t::part_t part_t
Contains the Multi-jagged algorthm.
Adapter::scalar_t scalar_t
A PartitioningSolution is a solution to a partitioning problem.
Tpetra::CrsGraph< lno_t, gno_t, node_t > graph_t
void determineRegularity()
Tpetra::CrsMatrix< scalar_t, lno_t, gno_t, node_t > matrix_t
void MJwrapper(const Teuchos::RCP< const mvector_t > &coordinates, std::vector< const weight_t * > weights, std::vector< int > strides, const Teuchos::RCP< PartitioningSolution< Adapter >> &solution)
Algorithm defines the base class for all algorithms.
map_t::local_ordinal_type lno_t
An adapter for Xpetra::MultiVector.
int AnasaziWrapper(const int numEigenVectors)
Tpetra::Operator< scalar_t, lno_t, gno_t, node_t > op_t
Defines the CoordinateModel classes.
typename Adapter::offset_t offset_t
void setPreconditioner(Teuchos::RCP< problem_t > &problem)
typename Adapter::scalar_t weight_t
Teuchos::RCP< matrix_t > computeNormalizedLaplacian(bool AHat=false)
Multi Jagged coordinate partitioning algorithm.