49 #include "Epetra_SerialComm.h"
50 #include "Epetra_Comm.h"
51 #include "Epetra_Map.h"
52 #include "Epetra_RowMatrix.h"
53 #include "Epetra_CrsMatrix.h"
54 #include "Epetra_Vector.h"
55 #include "Epetra_MultiVector.h"
56 #include "Epetra_Util.h"
57 #include "Teuchos_ParameterList.hpp"
58 #include "Teuchos_RefCountPtr.hpp"
62 using namespace Teuchos;
74 DropTolerance_(1e-12),
75 IsInitialized_(
false),
84 ApplyInverseTime_(0.0),
86 ApplyInverseFlops_(0.0),
142 cerr <<
"Caught an exception while parsing the parameter list" << endl;
143 cerr <<
"This typically means that a parameter was set with the" << endl;
144 cerr <<
"wrong type (for example, int instead of double). " << endl;
145 cerr <<
"please check the documentation for the type required by each parameer." << endl;
165 cout <<
" There are too many processors !!! " << endl;
166 cerr <<
"Ifpack_IKLU can be used with Comm().NumProc() == 1" << endl;
167 cerr <<
"only. This class is a subdomain solver for Ifpack_AdditiveSchwarz," << endl;
168 cerr <<
"and it is currently not meant to be used otherwise." << endl;
180 std::vector<int> RowIndices(Length);
181 std::vector<double> RowValues(Length);
193 &RowValues[0],&RowIndices[0]));
194 for (
int j = 0 ; j < RowNnz ; ++j) {
195 csrA_->
j[count++] = RowIndices[j];
205 cout <<
"AMD Perm (from inside KLU) [" << i <<
"] = " <<
cssS_->
q[i] << endl;
242 bool distributed = (
Comm().
NumProc() > 1)?
true:
false;
243 #if !defined(EPETRA_NO_32BIT_GLOBAL_INDICES) || !defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
256 std::vector<int> RowIndices(Length);
257 std::vector<double> RowValues(Length);
264 &RowValues[0],&RowIndices[0]));
267 cout <<
"The number of nonzeros for this row does not math the expected number of nonzeros!!!" << endl;
269 for (
int j = 0 ; j < RowNnz ; ++j) {
271 csrA_->
x[count++] = RowValues[j];
283 std::vector<int> numEntriesL( NumMyRows_ ), numEntriesU( NumMyRows_ );
285 numEntriesL[i] = ( L_tmp->
p[i+1] - L_tmp->
p[i] );
286 numEntriesU[i] = ( U_tmp->
p[i+1] - U_tmp->
p[i] );
292 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
295 L_->InsertGlobalValues( i, numEntriesL[i], &(L_tmp->
x[L_tmp->
p[i]]), &(L_tmp->
j[L_tmp->
p[i]]) );
296 U_->InsertGlobalValues( i, numEntriesU[i], &(U_tmp->
x[U_tmp->
p[i]]), &(U_tmp->
j[U_tmp->
p[i]]) );
301 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
304 const int MaxNumEntries_L_U =
std::max(
L_->MaxNumEntries(),
U_->MaxNumEntries());
305 std::vector<long long> entries(MaxNumEntries_L_U);
308 std::copy(&(L_tmp->
j[L_tmp->
p[i]]), &(L_tmp->
j[L_tmp->
p[i]]) + numEntriesL[i], entries.begin());
309 L_->InsertGlobalValues( i, numEntriesL[i], &(L_tmp->
x[L_tmp->
p[i]]), &(entries[0]) );
311 std::copy(&(U_tmp->
j[U_tmp->
p[i]]), &(U_tmp->
j[U_tmp->
p[i]]) + numEntriesU[i], entries.begin());
312 U_->InsertGlobalValues( i, numEntriesU[i], &(U_tmp->
x[U_tmp->
p[i]]), &(entries[0]) );
317 throw "Ifpack_IKLU::Compute: GlobalIndices type unknown for SerialMap_";
322 long long MyNonzeros =
L_->NumGlobalNonzeros64() +
U_->NumGlobalNonzeros64();
341 if (X.NumVectors() != Y.NumVectors())
363 for (
int j=0; j<X.NumVectors(); ++j) {
364 Xcopy->ReplaceMyValue( invq[i], j, (*X(j))[i] );
383 for (
int j=0; j<Y.NumVectors(); ++j) {
384 Y.ReplaceMyValue(
csrnN_->
perm[i], j, (*(*Ytemp)(j))[i] );
389 #ifdef IFPACK_FLOPCOUNTERS
408 const int MaxIters,
const double Tol,
427 if (!
Comm().MyPID()) {
429 os <<
"================================================================================" << endl;
430 os <<
"Ifpack_IKLU: " <<
Label() << endl << endl;
434 os <<
"Relax value = " <<
RelaxValue() << endl;
435 os <<
"Condition number estimate = " <<
Condest() << endl;
441 <<
" % of A)" << endl;
442 os <<
"nonzeros / rows = "
446 os <<
"Phase # calls Total Time (s) Total MFlops MFlops/s" << endl;
447 os <<
"----- ------- -------------- ------------ --------" << endl;
450 <<
" 0.0 0.0" << endl;
451 os <<
"Compute() " << std::setw(5) <<
NumCompute()
457 os <<
" " << std::setw(15) << 0.0 << endl;
464 os <<
" " << std::setw(15) << 0.0 << endl;
465 os <<
"================================================================================" << endl;
Ifpack_IKLU(const Epetra_RowMatrix *A)
Ifpack_IKLU constuctor with variable number of indices per row.
const Epetra_RowMatrix & Matrix() const
Returns a reference to the matrix to be preconditioned.
bool IsComputed_
true if this object has been computed
const Epetra_RowMatrix & A_
reference to the matrix to be preconditioned.
double ApplyInverseTime_
Contains the time for all successful calls to ApplyInverse().
virtual const Epetra_Map & RowMatrixRowMap() const =0
csrn * csr_lu(const csr *A, const css *S, double tol)
T & get(ParameterList &l, const std::string &name)
double ElapsedTime(void) const
virtual std::ostream & Print(std::ostream &os) const
Prints basic information on iostream. This function is used by operator<<.
bool IsInitialized_
true if this object has been initialized
csr * csr_spalloc(int m, int n, int nzmax, int values, int triplet)
Teuchos::RefCountPtr< Epetra_SerialComm > SerialComm_
int SetParameters(Teuchos::ParameterList ¶meterlis)
Set parameters using a Teuchos::ParameterList object.
long long NumGlobalNonzeros64() const
int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
virtual int NumApplyInverse() const
Returns the number of calls to ApplyInverse().
int NumMyNonzeros_
Number of local nonzeros.
int NumApplyInverse_
Contains the number of successful call to ApplyInverse().
Epetra_Time Time_
Used for timing purposed.
virtual double InitializeTime() const
Returns the time spent in Initialize().
double DropTolerance_
Discards all elements below this tolerance.
virtual double ComputeTime() const
Returns the time spent in Compute().
Ifpack_CondestType
Ifpack_CondestType: enum to define the type of condition number estimate.
virtual int SumAll(double *PartialSums, double *GlobalSums, int Count) const =0
csr * csrA_
Containers for the matrix storage and permutation.
virtual int MaxNumEntries() const =0
long long GlobalNonzeros_
Global number of nonzeros in L and U factors.
virtual const Epetra_Comm & Comm() const =0
double ApplyInverseFlops_
Contain sthe number of flops for ApplyInverse().
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
double LevelOfFill_
Level-of-fill.
Teuchos::RefCountPtr< Epetra_CrsMatrix > U_
U factor.
virtual double ApplyInverseFlops() const
Returns the number of flops in the application of the preconditioner.
const char * Label() const
Returns the label of this object.
virtual const Epetra_BlockMap & Map() const =0
bool IsInitialized() const
Returns true if the preconditioner has been successfully initialized.
csrn * csr_nfree(csrn *N)
int NumCompute_
Contains the number of successful call to Compute().
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
bool IsComputed() const
If factor is completed, this query returns true, otherwise it returns false.
virtual int NumMyRows() const =0
double Athresh_
Absolute threshold.
int NumInitialize_
Contains the number of successful calls to Initialize().
std::string Ifpack_toString(const int &x)
Converts an integer to std::string.
std::string Label_
Label for this object.
virtual int NumCompute() const
Returns the number of calls to Compute().
double Condest_
Condition number estimate.
double LevelOfFill() const
virtual long long NumGlobalNonzeros64() const =0
virtual ~Ifpack_IKLU()
Ifpack_IKLU Destructor.
double Ifpack_Condest(const Ifpack_Preconditioner &IFP, const Ifpack_CondestType CT, const int MaxIters, const double Tol, Epetra_RowMatrix *Matrix)
double DropTolerance() const
Gets the dropping tolerance.
double Rthresh_
Relative threshold.
virtual int NumProc() const =0
double InitializeTime_
Contains the time for all successful calls to Initialize().
int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Ifpack_IKLU forward/back solve on a Epetra_MultiVector X in Y...
double Relax_
relaxation value
double RelativeThreshold() const
Get relative threshold value.
virtual int ExtractMyRowCopy(int MyRow, int Length, int &NumEntries, double *Values, int *Indices) const =0
bool operator()(const double &x, const double &y)
double RelaxValue() const
Set relative threshold value.
Teuchos::RefCountPtr< Epetra_CrsMatrix > L_
L factor.
bool UseTranspose_
true if transpose has to be used.
const Epetra_Comm & Comm() const
Returns the Epetra_BlockMap object associated with the range of this matrix operator.
virtual int NumInitialize() const
Returns the number of calls to Initialize().
css * csr_sqr(int order, const csr *A)
csrn * csrnN_
Container for the L and U factor.
Teuchos::RefCountPtr< Epetra_Map > SerialMap_
#define IFPACK_CHK_ERR(ifpack_err)
void ResetStartTime(void)
double ComputeTime_
Contains the time for all successful calls to Compute().
double AbsoluteThreshold() const
Get absolute threshold value.
int NumMyRows_
Number of local rows.
int Initialize()
Initialize L and U with values from user matrix A.
virtual int NumMyNonzeros() const =0
virtual long long NumGlobalRows64() const =0
double Condest() const
Returns the computed estimated condition number, or -1.0 if no computed.
void Destroy()
Releases all allocated memory.
int Compute()
Compute IC factor U using the specified graph, diagonal perturbation thresholds and relaxation parame...
virtual double ComputeFlops() const
Returns the number of flops in the computation phase.