52 #ifndef AMESOS2_KLU2_DEF_HPP
53 #define AMESOS2_KLU2_DEF_HPP
55 #include <Teuchos_Tuple.hpp>
56 #include <Teuchos_ParameterList.hpp>
57 #include <Teuchos_StandardParameterEntryValidators.hpp>
65 template <
class Matrix,
class Vector>
67 Teuchos::RCP<const Matrix> A,
68 Teuchos::RCP<Vector> X,
69 Teuchos::RCP<const Vector> B )
72 , is_contiguous_(true)
74 ::KLU2::klu_defaults<klu2_dtype, local_ordinal_type> (&(data_.common_)) ;
75 data_.symbolic_ = NULL;
76 data_.numeric_ = NULL;
83 template <
class Matrix,
class Vector>
91 if (data_.symbolic_ != NULL)
92 ::KLU2::klu_free_symbolic<klu2_dtype, local_ordinal_type>
93 (&(data_.symbolic_), &(data_.common_)) ;
94 if (data_.numeric_ != NULL)
95 ::KLU2::klu_free_numeric<klu2_dtype, local_ordinal_type>
96 (&(data_.numeric_), &(data_.common_)) ;
109 template <
class Matrix,
class Vector>
112 return (this->root_ && (this->matrixA_->getComm()->getSize() == 1) && is_contiguous_);
115 template<
class Matrix,
class Vector>
121 #ifdef HAVE_AMESOS2_TIMERS
122 Teuchos::TimeMonitor preOrderTimer(this->timers_.preOrderTime_);
129 template <
class Matrix,
class Vector>
133 if (data_.symbolic_ != NULL) {
134 ::KLU2::klu_free_symbolic<klu2_dtype, local_ordinal_type>
135 (&(data_.symbolic_), &(data_.common_)) ;
138 if ( single_proc_optimization() ) {
139 host_ordinal_type_array host_row_ptr_view;
140 host_ordinal_type_array host_cols_view;
141 this->matrixA_->returnRowPtr_kokkos_view(host_row_ptr_view);
142 this->matrixA_->returnColInd_kokkos_view(host_cols_view);
143 data_.symbolic_ = ::KLU2::klu_analyze<klu2_dtype, local_ordinal_type>
144 ((local_ordinal_type)this->globalNumCols_, host_row_ptr_view.data(),
145 host_cols_view.data(), &(data_.common_)) ;
149 data_.symbolic_ = ::KLU2::klu_analyze<klu2_dtype, local_ordinal_type>
150 ((local_ordinal_type)this->globalNumCols_, host_col_ptr_view_.data(),
151 host_rows_view_.data(), &(data_.common_)) ;
159 template <
class Matrix,
class Vector>
173 #ifdef HAVE_AMESOS2_TIMERS
174 Teuchos::TimeMonitor numFactTimer(this->timers_.numFactTime_);
177 if (data_.numeric_ != NULL) {
178 ::KLU2::klu_free_numeric<klu2_dtype, local_ordinal_type>
179 (&(data_.numeric_), &(data_.common_));
182 if ( single_proc_optimization() ) {
183 host_ordinal_type_array host_row_ptr_view;
184 host_ordinal_type_array host_cols_view;
185 this->matrixA_->returnRowPtr_kokkos_view(host_row_ptr_view);
186 this->matrixA_->returnColInd_kokkos_view(host_cols_view);
187 this->matrixA_->returnValues_kokkos_view(host_nzvals_view_);
188 klu2_dtype * pValues = function_map::convert_scalar(host_nzvals_view_.data());
189 data_.numeric_ = ::KLU2::klu_factor<klu2_dtype, local_ordinal_type>
190 (host_row_ptr_view.data(), host_cols_view.data(), pValues,
191 data_.symbolic_, &(data_.common_));
194 klu2_dtype * pValues = function_map::convert_scalar(host_nzvals_view_.data());
195 data_.numeric_ = ::KLU2::klu_factor<klu2_dtype, local_ordinal_type>
196 (host_col_ptr_view_.data(), host_rows_view_.data(), pValues,
197 data_.symbolic_, &(data_.common_));
205 if(data_.numeric_ ==
nullptr) {
212 this->setNnzLU( as<size_t>((data_.numeric_)->lnz) + as<size_t>((data_.numeric_)->unz) );
219 Teuchos::broadcast(*(this->matrixA_->getComm()), 0, &info);
221 TEUCHOS_TEST_FOR_EXCEPTION(info > 0, std::runtime_error,
222 "KLU2 numeric factorization failed");
227 template <
class Matrix,
class Vector>
236 const global_size_type ld_rhs = this->root_ ? X->getGlobalLength() : 0;
237 const size_t nrhs = X->getGlobalNumVectors();
242 #ifdef HAVE_AMESOS2_TIMERS
243 Teuchos::TimeMonitor mvConvTimer(this->timers_.vecConvTime_);
244 Teuchos::TimeMonitor redistTimer( this->timers_.vecRedistTime_ );
246 const bool initialize_data =
true;
247 const bool do_not_initialize_data =
false;
248 if ( single_proc_optimization() && nrhs == 1 ) {
250 bDidAssignB = Util::get_1d_copy_helper_kokkos_view<MultiVecAdapter<Vector>,
251 host_solve_array_t>::do_get(initialize_data, B, bValues_, as<size_t>(ld_rhs));
253 bDidAssignX = Util::get_1d_copy_helper_kokkos_view<MultiVecAdapter<Vector>,
254 host_solve_array_t>::do_get(do_not_initialize_data, X, xValues_, as<size_t>(ld_rhs));
258 bDidAssignB = Util::get_1d_copy_helper_kokkos_view<MultiVecAdapter<Vector>,
259 host_solve_array_t>::do_get(initialize_data, B, bValues_,
261 (is_contiguous_ ==
true) ? ROOTED : CONTIGUOUS_AND_ROOTED,
262 this->rowIndexBase_);
265 bDidAssignX = Util::get_1d_copy_helper_kokkos_view<MultiVecAdapter<Vector>,
266 host_solve_array_t>::do_get(do_not_initialize_data, X, xValues_,
268 (is_contiguous_ ==
true) ? ROOTED : CONTIGUOUS_AND_ROOTED,
269 this->rowIndexBase_);
278 Kokkos::deep_copy(xValues_, bValues_);
286 klu2_dtype * pxValues = function_map::convert_scalar(xValues_.data());
287 klu2_dtype * pbValues = function_map::convert_scalar(bValues_.data());
291 TEUCHOS_TEST_FOR_EXCEPTION(pbValues ==
nullptr,
292 std::runtime_error,
"Amesos2 Runtime Error: b_vector returned null ");
294 TEUCHOS_TEST_FOR_EXCEPTION(pxValues ==
nullptr,
295 std::runtime_error,
"Amesos2 Runtime Error: x_vector returned null ");
298 if ( single_proc_optimization() && nrhs == 1 ) {
299 #ifdef HAVE_AMESOS2_TIMERS
300 Teuchos::TimeMonitor solveTimer(this->timers_.solveTime_);
308 ::KLU2::klu_tsolve2<klu2_dtype, local_ordinal_type>
309 (data_.symbolic_, data_.numeric_,
310 (local_ordinal_type)this->globalNumCols_,
311 (local_ordinal_type)nrhs,
312 pbValues, pxValues, &(data_.common_)) ;
315 ::KLU2::klu_solve2<klu2_dtype, local_ordinal_type>
316 (data_.symbolic_, data_.numeric_,
317 (local_ordinal_type)this->globalNumCols_,
318 (local_ordinal_type)nrhs,
319 pbValues, pxValues, &(data_.common_)) ;
330 #ifdef HAVE_AMESOS2_TIMERS
331 Teuchos::TimeMonitor solveTimer(this->timers_.solveTime_);
338 if ( single_proc_optimization() ) {
339 ::KLU2::klu_tsolve<klu2_dtype, local_ordinal_type>
340 (data_.symbolic_, data_.numeric_,
341 (local_ordinal_type)this->globalNumCols_,
342 (local_ordinal_type)nrhs,
343 pxValues, &(data_.common_)) ;
346 ::KLU2::klu_solve<klu2_dtype, local_ordinal_type>
347 (data_.symbolic_, data_.numeric_,
348 (local_ordinal_type)this->globalNumCols_,
349 (local_ordinal_type)nrhs,
350 pxValues, &(data_.common_)) ;
358 if ( single_proc_optimization() ) {
359 ::KLU2::klu_solve<klu2_dtype, local_ordinal_type>
360 (data_.symbolic_, data_.numeric_,
361 (local_ordinal_type)this->globalNumCols_,
362 (local_ordinal_type)nrhs,
363 pxValues, &(data_.common_)) ;
366 ::KLU2::klu_tsolve<klu2_dtype, local_ordinal_type>
367 (data_.symbolic_, data_.numeric_,
368 (local_ordinal_type)this->globalNumCols_,
369 (local_ordinal_type)nrhs,
370 pxValues, &(data_.common_)) ;
379 #ifdef HAVE_AMESOS2_TIMERS
380 Teuchos::TimeMonitor redistTimer( this->timers_.vecRedistTime_ );
383 Util::put_1d_data_helper_kokkos_view<
386 (is_contiguous_ ==
true) ? ROOTED : CONTIGUOUS_AND_ROOTED,
387 this->rowIndexBase_);
395 template <
class Matrix,
class Vector>
402 return( this->matrixA_->getGlobalNumRows() == this->matrixA_->getGlobalNumCols() );
406 template <
class Matrix,
class Vector>
411 using Teuchos::getIntegralValue;
412 using Teuchos::ParameterEntryValidator;
414 RCP<const Teuchos::ParameterList> valid_params = getValidParameters_impl();
416 transFlag_ = this->control_.useTranspose_ ? 1: 0;
418 if( parameterList->isParameter(
"Trans") ){
419 RCP<const ParameterEntryValidator> trans_validator = valid_params->getEntry(
"Trans").validator();
420 parameterList->getEntry(
"Trans").setValidator(trans_validator);
422 transFlag_ = getIntegralValue<int>(*parameterList,
"Trans");
425 if( parameterList->isParameter(
"IsContiguous") ){
426 is_contiguous_ = parameterList->get<
bool>(
"IsContiguous");
431 template <
class Matrix,
class Vector>
432 Teuchos::RCP<const Teuchos::ParameterList>
436 using Teuchos::tuple;
437 using Teuchos::ParameterList;
438 using Teuchos::setStringToIntegralParameter;
440 static Teuchos::RCP<const Teuchos::ParameterList> valid_params;
442 if( is_null(valid_params) )
444 Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
446 pl->set(
"Equil",
true,
"Whether to equilibrate the system before solve, does nothing now");
447 pl->set(
"IsContiguous",
true,
"Whether GIDs contiguous");
449 setStringToIntegralParameter<int>(
"Trans",
"NOTRANS",
450 "Solve for the transpose system or not",
451 tuple<string>(
"NOTRANS",
"TRANS",
"CONJ"),
452 tuple<string>(
"Solve with transpose",
453 "Do not solve with transpose",
454 "Solve with the conjugate transpose"),
464 template <
class Matrix,
class Vector>
470 if(current_phase == SOLVE)
return(
false);
472 if ( single_proc_optimization() ) {
478 #ifdef HAVE_AMESOS2_TIMERS
479 Teuchos::TimeMonitor convTimer(this->timers_.mtxConvTime_);
484 host_nzvals_view_ = host_value_type_array(
485 Kokkos::ViewAllocateWithoutInitializing(
"host_nzvals_view_"), this->globalNumNonZeros_);
486 host_rows_view_ = host_ordinal_type_array(
487 Kokkos::ViewAllocateWithoutInitializing(
"host_rows_view_"), this->globalNumNonZeros_);
488 host_col_ptr_view_ = host_ordinal_type_array(
489 Kokkos::ViewAllocateWithoutInitializing(
"host_col_ptr_view_"), this->globalNumRows_ + 1);
492 local_ordinal_type nnz_ret = 0;
494 #ifdef HAVE_AMESOS2_TIMERS
495 Teuchos::TimeMonitor mtxRedistTimer( this->timers_.mtxRedistTime_ );
500 ::do_get(this->matrixA_.ptr(), host_nzvals_view_, host_rows_view_, host_col_ptr_view_, nnz_ret,
501 (is_contiguous_ ==
true) ? ROOTED : CONTIGUOUS_AND_ROOTED,
503 this->rowIndexBase_);
507 TEUCHOS_TEST_FOR_EXCEPTION( nnz_ret != as<local_ordinal_type>(this->globalNumNonZeros_),
509 "Did not get the expected number of non-zero vals");
518 template<
class Matrix,
class Vector>
524 #endif // AMESOS2_KLU2_DEF_HPP
Amesos2::SolverCore: A templated interface for interaction with third-party direct sparse solvers...
Definition: Amesos2_SolverCore_decl.hpp:105
KLU2(Teuchos::RCP< const Matrix > A, Teuchos::RCP< Vector > X, Teuchos::RCP< const Vector > B)
Initialize from Teuchos::RCP.
Definition: Amesos2_KLU2_def.hpp:66
A generic helper class for getting a CCS representation of a Matrix.
Definition: Amesos2_Util.hpp:648
int solve_impl(const Teuchos::Ptr< MultiVecAdapter< Vector > > X, const Teuchos::Ptr< const MultiVecAdapter< Vector > > B) const
KLU2 specific solve.
Definition: Amesos2_KLU2_def.hpp:229
EPhase
Used to indicate a phase in the direct solution.
Definition: Amesos2_TypeDecl.hpp:65
Amesos2 KLU2 declarations.
bool loadA_impl(EPhase current_phase)
Reads matrix data into internal structures.
Definition: Amesos2_KLU2_def.hpp:466
~KLU2()
Destructor.
Definition: Amesos2_KLU2_def.hpp:84
void setParameters_impl(const Teuchos::RCP< Teuchos::ParameterList > ¶meterList)
Definition: Amesos2_KLU2_def.hpp:408
int symbolicFactorization_impl()
Perform symbolic factorization of the matrix using KLU2.
Definition: Amesos2_KLU2_def.hpp:131
int preOrdering_impl()
Performs pre-ordering on the matrix to increase efficiency.
Definition: Amesos2_KLU2_def.hpp:117
bool matrixShapeOK_impl() const
Determines whether the shape of the matrix is OK for this solver.
Definition: Amesos2_KLU2_def.hpp:397
A Matrix adapter interface for Amesos2.
Definition: Amesos2_MatrixAdapter_decl.hpp:76
int numericFactorization_impl()
KLU2 specific numeric factorization.
Definition: Amesos2_KLU2_def.hpp:161
Amesos2 interface to the KLU2 package.
Definition: Amesos2_KLU2_decl.hpp:72
bool single_proc_optimization() const
can we optimize size_type and ordinal_type for straight pass through, also check that is_contiguous_ ...
Definition: Amesos2_KLU2_def.hpp:111
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters_impl() const
Definition: Amesos2_KLU2_def.hpp:433
A templated MultiVector class adapter for Amesos2.
Definition: Amesos2_MultiVecAdapter_decl.hpp:176