53 #ifndef AMESOS2_CHOLMOD_DEF_HPP 54 #define AMESOS2_CHOLMOD_DEF_HPP 56 #include <Teuchos_Tuple.hpp> 57 #include <Teuchos_ParameterList.hpp> 58 #include <Teuchos_StandardParameterEntryValidators.hpp> 67 template <
class Matrix,
class Vector>
69 Teuchos::RCP<const Matrix> A,
70 Teuchos::RCP<Vector> X,
71 Teuchos::RCP<const Vector> B )
79 CHOL::cholmod_start(&data_.c);
80 (&data_.c)->nmethods=9;
81 (&data_.c)->quick_return_if_not_posdef=1;
90 template <
class Matrix,
class Vector>
93 CHOL::cholmod_free_factor (&(data_.L), &(data_.c));
95 CHOL::cholmod_free_dense(&(data_.Y), &data_.c);
96 CHOL::cholmod_free_dense(&(data_.E), &data_.c);
98 CHOL::cholmod_finish(&(data_.c));
101 template<
class Matrix,
class Vector>
105 #ifdef HAVE_AMESOS2_TIMERS 106 Teuchos::TimeMonitor preOrderTimer(this->timers_.preOrderTime_);
108 data_.L = CHOL::cholmod_analyze (&data_.A, &(data_.c));
118 template <
class Matrix,
class Vector>
122 if ( !skip_symfact ){
123 #ifdef HAVE_AMESOS2_TIMERS 124 Teuchos::TimeMonitor symFactTimer(this->timers_.symFactTime_);
126 CHOL::cholmod_resymbol (&data_.A, NULL, 0,
true, data_.L, &(data_.c));
133 skip_symfact =
false;
140 template <
class Matrix,
class Vector>
148 #ifdef HAVE_AMESOS2_DEBUG 149 TEUCHOS_TEST_FOR_EXCEPTION(data_.A.ncol != as<size_t>(this->globalNumCols_),
151 "Error in converting to cholmod_sparse: wrong number of global columns." );
152 TEUCHOS_TEST_FOR_EXCEPTION(data_.A.nrow != as<size_t>(this->globalNumRows_),
154 "Error in converting to cholmod_sparse: wrong number of global rows." );
158 #ifdef HAVE_AMESOS2_TIMERS 159 Teuchos::TimeMonitor numFactTimer(this->timers_.numFactTime_);
162 #ifdef HAVE_AMESOS2_VERBOSE_DEBUG 163 std::cout <<
"Cholmod:: Before numeric factorization" << std::endl;
164 std::cout <<
"nzvals_ : " << nzvals_.toString() << std::endl;
165 std::cout <<
"rowind_ : " << rowind_.toString() << std::endl;
166 std::cout <<
"colptr_ : " << colptr_.toString() << std::endl;
169 CHOL::cholmod_factorize(&data_.A, data_.L, &(data_.c));
171 info = (&data_.c)->status;
178 Teuchos::broadcast(*(this->matrixA_->getComm()), 0, &info);
180 TEUCHOS_TEST_FOR_EXCEPTION(info == 2,
182 "Memory allocation failure in Cholmod factorization");
188 template <
class Matrix,
class Vector>
195 const global_size_type ld_rhs = X->getGlobalLength();
196 const size_t nrhs = X->getGlobalNumVectors();
198 Teuchos::RCP<
const Tpetra::Map<local_ordinal_type,
199 global_ordinal_type, node_type> > map2;
201 const size_t val_store_size = as<size_t>(ld_rhs * nrhs);
202 Teuchos::Array<chol_type> xValues(val_store_size);
203 Teuchos::Array<chol_type> bValues(val_store_size);
206 #ifdef HAVE_AMESOS2_TIMERS 207 Teuchos::TimeMonitor mvConvTimer(this->timers_.vecConvTime_);
208 Teuchos::TimeMonitor redistTimer( this->timers_.vecRedistTime_ );
211 chol_type>::do_get(B, bValues(),
213 ROOTED, this->rowIndexBase_);
221 #ifdef HAVE_AMESOS2_TIMERS 222 Teuchos::TimeMonitor mvConvTimer(this->timers_.vecConvTime_);
225 function_map::cholmod_init_dense(as<int>(this->globalNumRows_),
226 as<int>(nrhs), as<int>(ld_rhs),
227 bValues.getRawPtr(), &data_.b);
228 function_map::cholmod_init_dense(as<int>(this->globalNumRows_),
229 as<int>(nrhs), as<int>(ld_rhs),
230 xValues.getRawPtr(), &data_.x);
235 #ifdef HAVE_AMESOS2_TIMERS 236 Teuchos::TimeMonitor solveTimer(this->timers_.solveTime_);
239 CHOL::cholmod_dense *xtemp = &(data_.x);
241 CHOL::cholmod_solve2(CHOLMOD_A, data_.L, &data_.b, NULL,
242 &xtemp, NULL, &data_.Y, &data_.E, &data_.c);
244 ierr = (&data_.c)->status;
250 Teuchos::broadcast(*(this->getComm()), 0, &ierr);
252 TEUCHOS_TEST_FOR_EXCEPTION(ierr == -2, std::runtime_error,
"Ran out of memory" );
256 #ifdef HAVE_AMESOS2_TIMERS 257 Teuchos::TimeMonitor redistTimer(this->timers_.vecRedistTime_);
263 ROOTED, this->rowIndexBase_);
272 template <
class Matrix,
class Vector>
276 return( this->matrixA_->getGlobalNumRows() == this->matrixA_->getGlobalNumCols() );
280 template <
class Matrix,
class Vector>
285 using Teuchos::getIntegralValue;
286 using Teuchos::ParameterEntryValidator;
288 RCP<const Teuchos::ParameterList> valid_params = getValidParameters_impl();
291 if( parameterList->isParameter(
"Trans") ){
292 RCP<const ParameterEntryValidator> trans_validator = valid_params->getEntry(
"Trans").validator();
293 parameterList->getEntry(
"Trans").setValidator(trans_validator);
297 if( parameterList->isParameter(
"IterRefine") ){
298 RCP<const ParameterEntryValidator> refine_validator = valid_params->getEntry(
"IterRefine").validator();
299 parameterList->getEntry(
"IterRefine").setValidator(refine_validator);
304 if( parameterList->isParameter(
"ColPerm") ){
305 RCP<const ParameterEntryValidator> colperm_validator = valid_params->getEntry(
"ColPerm").validator();
306 parameterList->getEntry(
"ColPerm").setValidator(colperm_validator);
311 (&data_.c)->dbound = parameterList->get<
double>(
"dbound", 0.0);
313 bool prefer_upper = parameterList->get<
bool>(
"PreferUpper",
true);
315 (&data_.c)->prefer_upper = prefer_upper ? 1 : 0;
317 (&data_.c)->print = parameterList->get<
int>(
"print",3);
319 (&data_.c)->nmethods = parameterList->get<
int>(
"nmethods",0);
324 template <
class Matrix,
class Vector>
325 Teuchos::RCP<const Teuchos::ParameterList>
329 using Teuchos::tuple;
330 using Teuchos::ParameterList;
331 using Teuchos::EnhancedNumberValidator;
332 using Teuchos::setStringToIntegralParameter;
333 using Teuchos::stringToIntegralParameterEntryValidator;
335 static Teuchos::RCP<const Teuchos::ParameterList> valid_params;
337 if( is_null(valid_params) ){
338 Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
341 Teuchos::RCP<EnhancedNumberValidator<int> > print_validator
342 = Teuchos::rcp(
new EnhancedNumberValidator<int>(0,5));
344 Teuchos::RCP<EnhancedNumberValidator<int> > nmethods_validator
345 = Teuchos::rcp(
new EnhancedNumberValidator<int>(0,9));
347 pl->set(
"nmethods", 0,
"Specifies the number of different ordering methods to try", nmethods_validator);
349 pl->set(
"print", 3,
"Specifies the verbosity of the print statements", print_validator);
351 pl->set(
"dbound", 0.0,
352 "Specifies the smallest absolute value on the diagonal D for the LDL' factorization");
355 pl->set(
"Equil",
true,
"Whether to equilibrate the system before solve");
357 pl->set(
"PreferUpper",
true,
358 "Specifies whether the matrix will be " 359 "stored in upper triangular form.");
368 template <
class Matrix,
class Vector>
374 #ifdef HAVE_AMESOS2_TIMERS 375 Teuchos::TimeMonitor convTimer(this->timers_.mtxConvTime_);
380 nzvals_.resize(this->globalNumNonZeros_);
381 rowind_.resize(this->globalNumNonZeros_);
382 colptr_.resize(this->globalNumCols_ + 1);
387 #ifdef HAVE_AMESOS2_TIMERS 388 Teuchos::TimeMonitor mtxRedistTimer( this->timers_.mtxRedistTime_ );
391 TEUCHOS_TEST_FOR_EXCEPTION(this->rowIndexBase_ != this->columnIndexBase_,
393 "Row and column maps have different indexbase ");
396 nzvals_(), rowind_(),
397 colptr_(), nnz_ret,
ROOTED,
399 this->rowIndexBase_);
403 TEUCHOS_TEST_FOR_EXCEPTION(nnz_ret != as<int>(this->globalNumNonZeros_),
405 "Did not get the expected number of non-zero vals");
407 function_map::cholmod_init_sparse(as<size_t>(this->globalNumRows_),
408 as<size_t>(this->globalNumCols_),
409 as<size_t>(this->globalNumNonZeros_),
421 template<
class Matrix,
class Vector>
427 #endif // AMESOS2_CHOLMOD_DEF_HPP ~Cholmod()
Destructor.
Definition: Amesos2_Cholmod_def.hpp:91
int preOrdering_impl()
Performs pre-ordering on the matrix to increase efficiency.
Definition: Amesos2_Cholmod_def.hpp:103
Amesos2::SolverCore: A templated interface for interaction with third-party direct sparse solvers...
Definition: Amesos2_SolverCore_decl.hpp:105
EPhase
Used to indicate a phase in the direct solution.
Definition: Amesos2_TypeDecl.hpp:65
int solve_impl(const Teuchos::Ptr< MultiVecAdapter< Vector > > X, const Teuchos::Ptr< const MultiVecAdapter< Vector > > B) const
CHOLMOD specific solve.
Definition: Amesos2_Cholmod_def.hpp:190
Helper class for getting 1-D copies of multivectors.
Definition: Amesos2_MultiVecAdapter_decl.hpp:243
Definition: Amesos2_TypeDecl.hpp:142
Cholmod(Teuchos::RCP< const Matrix > A, Teuchos::RCP< Vector > X, Teuchos::RCP< const Vector > B)
Initialize from Teuchos::RCP.
Definition: Amesos2_Cholmod_def.hpp:68
Definition: Amesos2_AbstractConcreteMatrixAdapter.hpp:48
A generic helper class for getting a CCS representation of a Matrix.
Definition: Amesos2_Util.hpp:580
void setParameters_impl(const Teuchos::RCP< Teuchos::ParameterList > ¶meterList)
Definition: Amesos2_Cholmod_def.hpp:282
Amesos2 CHOLMOD declarations.
A Matrix adapter interface for Amesos2.
Definition: Amesos2_MatrixAdapter_decl.hpp:76
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters_impl() const
Definition: Amesos2_Cholmod_def.hpp:326
Definition: Amesos2_TypeDecl.hpp:127
bool loadA_impl(EPhase current_phase)
Reads matrix data into internal structures.
Definition: Amesos2_Cholmod_def.hpp:370
int symbolicFactorization_impl()
Perform symbolic factorization of the matrix using CHOLMOD.
Definition: Amesos2_Cholmod_def.hpp:120
Amesos2 interface to the SuperLU package.
Definition: Amesos2_Cholmod_decl.hpp:73
Helper class for putting 1-D data arrays into multivectors.
Definition: Amesos2_MultiVecAdapter_decl.hpp:296
A templated MultiVector class adapter for Amesos2.
Definition: Amesos2_MultiVecAdapter_decl.hpp:175
int numericFactorization_impl()
CHOLMOD specific numeric factorization.
Definition: Amesos2_Cholmod_def.hpp:142
bool matrixShapeOK_impl() const
Determines whether the shape of the matrix is OK for this solver.
Definition: Amesos2_Cholmod_def.hpp:274