IFPACK  Development
 All Classes Files Functions Variables Enumerations Friends Pages
Ifpack_ICT.cpp
1 /*@HEADER
2 // ***********************************************************************
3 //
4 // Ifpack: Object-Oriented Algebraic Preconditioner Package
5 // Copyright (2002) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 //@HEADER
41 */
42 
43 #include "Ifpack_ConfigDefs.h"
44 #include "Ifpack_Preconditioner.h"
45 #include "Ifpack_ICT.h"
46 #include "Ifpack_Condest.h"
47 #include "Ifpack_Utils.h"
48 #include "Ifpack_HashTable.h"
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"
59 #include <functional>
60 
61 //==============================================================================
62 // FIXME: allocate Comm_ and Time_ the first Initialize() call
63 Ifpack_ICT::Ifpack_ICT(const Epetra_RowMatrix* A) :
64  A_(*A),
65  Comm_(A_.Comm()),
66  Condest_(-1.0),
67  Athresh_(0.0),
68  Rthresh_(1.0),
69  LevelOfFill_(1.0),
70  DropTolerance_(0.0),
71  Relax_(0.0),
72  IsInitialized_(false),
73  IsComputed_(false),
74  UseTranspose_(false),
75  NumMyRows_(0),
76  NumInitialize_(0),
77  NumCompute_(0),
78  NumApplyInverse_(0),
79  InitializeTime_(0.0),
80  ComputeTime_(0.0),
81  ApplyInverseTime_(0.0),
82  ComputeFlops_(0.0),
83  ApplyInverseFlops_(0.0),
84  Time_(Comm()),
85  GlobalNonzeros_(0)
86 {
87  // do nothing here
88 }
89 
90 //==============================================================================
92 {
93  Destroy();
94 }
95 
96 //==============================================================================
97 void Ifpack_ICT::Destroy()
98 {
99  IsInitialized_ = false;
100  IsComputed_ = false;
101 }
102 
103 //==========================================================================
104 int Ifpack_ICT::SetParameters(Teuchos::ParameterList& List)
105 {
106 
107  try
108  {
109  LevelOfFill_ = List.get("fact: ict level-of-fill",LevelOfFill_);
110  Athresh_ = List.get("fact: absolute threshold", Athresh_);
111  Rthresh_ = List.get("fact: relative threshold", Rthresh_);
112  Relax_ = List.get("fact: relax value", Relax_);
113  DropTolerance_ = List.get("fact: drop tolerance", DropTolerance_);
114 
115  // set label
116  Label_ = "ICT (fill=" + Ifpack_toString(LevelOfFill())
117  + ", athr=" + Ifpack_toString(AbsoluteThreshold())
118  + ", rthr=" + Ifpack_toString(RelativeThreshold())
119  + ", relax=" + Ifpack_toString(RelaxValue())
120  + ", droptol=" + Ifpack_toString(DropTolerance())
121  + ")";
122 
123  return(0);
124  }
125  catch (...)
126  {
127  cerr << "Caught an exception while parsing the parameter list" << endl;
128  cerr << "This typically means that a parameter was set with the" << endl;
129  cerr << "wrong type (for example, int instead of double). " << endl;
130  cerr << "please check the documentation for the type required by each parameer." << endl;
131  IFPACK_CHK_ERR(-1);
132  }
133 }
134 
135 //==========================================================================
137 {
138  // clean data if present
139  Destroy();
140 
141  Time_.ResetStartTime();
142 
143  // matrix must be square. Check only on one processor
144  if (Comm().NumProc() == 1 && Matrix().NumMyRows() != Matrix().NumMyCols())
145  IFPACK_CHK_ERR(-2);
146 
147  NumMyRows_ = Matrix().NumMyRows();
148 
149  // nothing else to do here
150  IsInitialized_ = true;
151  ++NumInitialize_;
152  InitializeTime_ += Time_.ElapsedTime();
153 
154  return(0);
155 }
156 
157 //==========================================================================
158 template<typename int_type>
159 int Ifpack_ICT::TCompute()
160 {
161  if (!IsInitialized())
162  IFPACK_CHK_ERR(Initialize());
163 
164  Time_.ResetStartTime();
165  IsComputed_ = false;
166 
167  NumMyRows_ = A_.NumMyRows();
168  int Length = A_.MaxNumEntries();
169  std::vector<int> RowIndices(Length);
170  std::vector<double> RowValues(Length);
171 
172  bool distributed = (Comm().NumProc() > 1)?true:false;
173 
174 #if !defined(EPETRA_NO_32BIT_GLOBAL_INDICES) || !defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
175  if (distributed)
176  {
177  SerialComm_ = Teuchos::rcp(new Epetra_SerialComm);
178 #if !defined(EPETRA_NO_32BIT_GLOBAL_INDICES)
179  if(A_.RowMatrixRowMap().GlobalIndicesInt())
180  SerialMap_ = Teuchos::rcp(new Epetra_Map(NumMyRows_, 0, *SerialComm_));
181  else
182 #endif
183 #if !defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
184  if(A_.RowMatrixRowMap().GlobalIndicesLongLong())
185  SerialMap_ = Teuchos::rcp(new Epetra_Map((long long) NumMyRows_, (long long) 0, *SerialComm_));
186  else
187 #endif
188  throw "Ifpack_ICT::TCompute: Global indices unknown.";
189  assert (SerialComm_.get() != 0);
190  assert (SerialMap_.get() != 0);
191  }
192  else
193  SerialMap_ = Teuchos::rcp(const_cast<Epetra_Map*>(&A_.RowMatrixRowMap()), false);
194 #endif
195 
196  int RowNnz;
197 #ifdef IFPACK_FLOPCOUNTERS
198  double flops = 0.0;
199 #endif
200 
201  H_ = Teuchos::rcp(new Epetra_CrsMatrix(Copy,*SerialMap_,0));
202  if (H_.get() == 0)
203  IFPACK_CHK_ERR(-5); // memory allocation error
204 
205  // get A(0,0) element and insert it (after sqrt)
206  IFPACK_CHK_ERR(A_.ExtractMyRowCopy(0,Length,RowNnz,
207  &RowValues[0],&RowIndices[0]));
208 
209  // skip off-processor elements
210  if (distributed)
211  {
212  int count = 0;
213  for (int i = 0 ;i < RowNnz ; ++i)
214  {
215  if (RowIndices[i] < NumMyRows_){
216  RowIndices[count] = RowIndices[i];
217  RowValues[count] = RowValues[i];
218  ++count;
219  }
220  else
221  continue;
222  }
223  RowNnz = count;
224  }
225 
226  // modify diagonal
227  double diag_val = 0.0;
228  for (int i = 0 ;i < RowNnz ; ++i) {
229  if (RowIndices[i] == 0) {
230  double& v = RowValues[i];
231  diag_val = AbsoluteThreshold() * EPETRA_SGN(v) +
232  RelativeThreshold() * v;
233  break;
234  }
235  }
236 
237  diag_val = sqrt(diag_val);
238  int_type diag_idx = 0;
239  EPETRA_CHK_ERR(H_->InsertGlobalValues(0,1,&diag_val, &diag_idx));
240 
241  // The 10 is just a small constant to limit collisons as the actual keys
242  // we store are the indices and not integers
243  // [0..A_.MaxNumEntries()*LevelofFill()].
244  TIfpack_HashTable<int_type> Hash( 10 * A_.MaxNumEntries() * LevelOfFill(), 1);
245 
246  // start factorization for line 1
247  for (int row_i = 1 ; row_i < NumMyRows_ ; ++row_i) {
248 
249  // get row `row_i' of the matrix
250  IFPACK_CHK_ERR(A_.ExtractMyRowCopy(row_i,Length,RowNnz,
251  &RowValues[0],&RowIndices[0]));
252 
253  // skip off-processor elements
254  if (distributed)
255  {
256  int count = 0;
257  for (int i = 0 ;i < RowNnz ; ++i)
258  {
259  if (RowIndices[i] < NumMyRows_){
260  RowIndices[count] = RowIndices[i];
261  RowValues[count] = RowValues[i];
262  ++count;
263  }
264  else
265  continue;
266  }
267  RowNnz = count;
268  }
269 
270  // number of nonzeros in this row are defined as the nonzeros
271  // of the matrix, plus the level of fill
272  int LOF = (int)(LevelOfFill() * RowNnz);
273  if (LOF == 0) LOF = 1;
274 
275  // convert line `row_i' into hash for fast access
276  Hash.reset();
277 
278  double h_ii = 0.0;
279  for (int i = 0 ; i < RowNnz ; ++i) {
280  if (RowIndices[i] == row_i) {
281  double& v = RowValues[i];
282  h_ii = AbsoluteThreshold() * EPETRA_SGN(v) + RelativeThreshold() * v;
283  }
284  else if (RowIndices[i] < row_i)
285  {
286  Hash.set(RowIndices[i], RowValues[i], true);
287  }
288  }
289 
290  // form element (row_i, col_j)
291  // I start from the first row that has a nonzero column
292  // index in row_i.
293  for (int col_j = RowIndices[0] ; col_j < row_i ; ++col_j) {
294 
295  double h_ij = 0.0, h_jj = 0.0;
296  // note: get() returns 0.0 if col_j is not found
297  h_ij = Hash.get(col_j);
298 
299  // get pointers to row `col_j'
300  int_type* ColIndices;
301  double* ColValues;
302  int ColNnz;
303  int_type col_j_GID = (int_type) H_->RowMap().GID64(col_j);
304  H_->ExtractGlobalRowView(col_j_GID, ColNnz, ColValues, ColIndices);
305 
306  for (int k = 0 ; k < ColNnz ; ++k) {
307  int_type col_k = ColIndices[k];
308 
309  if (col_k == col_j)
310  h_jj = ColValues[k];
311  else {
312  double xxx = Hash.get(col_k);
313  if (xxx != 0.0)
314  {
315  h_ij -= ColValues[k] * xxx;
316 #ifdef IFPACK_FLOPCOUNTERS
317  flops += 2.0;
318 #endif
319  }
320  }
321  }
322 
323  h_ij /= h_jj;
324 
325  if (IFPACK_ABS(h_ij) > DropTolerance_)
326  {
327  Hash.set(col_j, h_ij);
328  }
329 
330 #ifdef IFPACK_FLOPCOUNTERS
331  // only approx
332  ComputeFlops_ += 2.0 * flops + 1.0;
333 #endif
334  }
335 
336  int size = Hash.getNumEntries();
337 
338  std::vector<double> AbsRow(size);
339  int count = 0;
340 
341  // +1 because I use the extra position for diagonal in insert
342  std::vector<int_type> keys(size + 1);
343  std::vector<double> values(size + 1);
344 
345  Hash.arrayify(&keys[0], &values[0]);
346 
347  for (int i = 0 ; i < size ; ++i)
348  {
349  AbsRow[i] = IFPACK_ABS(values[i]);
350  }
351  count = size;
352 
353  double cutoff = 0.0;
354  if (count > LOF) {
355  nth_element(AbsRow.begin(), AbsRow.begin() + LOF, AbsRow.begin() + count,
356 
357  std::greater<double>());
358  cutoff = AbsRow[LOF];
359  }
360 
361  for (int i = 0 ; i < size ; ++i)
362  {
363  h_ii -= values[i] * values[i];
364  }
365 
366  if (h_ii < 0.0) h_ii = 1e-12;;
367 
368  h_ii = sqrt(h_ii);
369 
370 #ifdef IFPACK_FLOPCOUNTERS
371  // only approx, + 1 == sqrt
372  ComputeFlops_ += 2 * size + 1;
373 #endif
374 
375  double DiscardedElements = 0.0;
376 
377  count = 0;
378  for (int i = 0 ; i < size ; ++i)
379  {
380  if (IFPACK_ABS(values[i]) > cutoff)
381  {
382  values[count] = values[i];
383  keys[count] = keys[i];
384  ++count;
385  }
386  else
387  DiscardedElements += values[i];
388  }
389 
390  if (RelaxValue() != 0.0) {
391  DiscardedElements *= RelaxValue();
392  h_ii += DiscardedElements;
393  }
394 
395  values[count] = h_ii;
396  keys[count] = (int_type) H_->RowMap().GID64(row_i);
397  ++count;
398 
399  H_->InsertGlobalValues((int_type) H_->RowMap().GID64(row_i), count, &(values[0]), (int_type*)&(keys[0]));
400  }
401 
402  IFPACK_CHK_ERR(H_->FillComplete());
403 
404 #if 0
405  // to check the complete factorization
406  Epetra_Vector LHS(Matrix().RowMatrixRowMap());
407  Epetra_Vector RHS1(Matrix().RowMatrixRowMap());
408  Epetra_Vector RHS2(Matrix().RowMatrixRowMap());
409  Epetra_Vector RHS3(Matrix().RowMatrixRowMap());
410  LHS.Random();
411 
412  Matrix().Multiply(false,LHS,RHS1);
413  H_->Multiply(true,LHS,RHS2);
414  H_->Multiply(false,RHS2,RHS3);
415 
416  RHS1.Update(-1.0, RHS3, 1.0);
417  cout << endl;
418  cout << RHS1;
419 #endif
420  long long MyNonzeros = H_->NumGlobalNonzeros64();
421  Comm().SumAll(&MyNonzeros, &GlobalNonzeros_, 1);
422 
423  IsComputed_ = true;
424 #ifdef IFPACK_FLOPCOUNTERS
425  double TotalFlops; // sum across all the processors
426  A_.Comm().SumAll(&flops, &TotalFlops, 1);
427  ComputeFlops_ += TotalFlops;
428 #endif
429  ++NumCompute_;
430  ComputeTime_ += Time_.ElapsedTime();
431 
432  return(0);
433 
434 }
435 
437 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
438  if(A_.RowMatrixRowMap().GlobalIndicesInt()) {
439  return TCompute<int>();
440  }
441  else
442 #endif
443 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
444  if(A_.RowMatrixRowMap().GlobalIndicesLongLong()) {
445  return TCompute<long long>();
446  }
447  else
448 #endif
449  throw "Ifpack_ICT::Compute: GlobalIndices type unknown for A_";
450 }
451 
452 //=============================================================================
453 int Ifpack_ICT::ApplyInverse(const Epetra_MultiVector& X,
454  Epetra_MultiVector& Y) const
455 {
456 
457  if (!IsComputed())
458  IFPACK_CHK_ERR(-3); // compute preconditioner first
459 
460  if (X.NumVectors() != Y.NumVectors())
461  IFPACK_CHK_ERR(-2); // Return error: X and Y not the same size
462 
463  Time_.ResetStartTime();
464 
465  // AztecOO gives X and Y pointing to the same memory location,
466  // need to create an auxiliary vector, Xcopy
467  Teuchos::RefCountPtr<const Epetra_MultiVector> Xcopy;
468  if (X.Pointers()[0] == Y.Pointers()[0])
469  Xcopy = Teuchos::rcp( new Epetra_MultiVector(X) );
470  else
471  Xcopy = Teuchos::rcp( &X, false );
472 
473  // NOTE: H_ is based on SerialMap_, while Xcopy is based
474  // on A.Map()... which are in general different. However, Solve()
475  // does not seem to care... which is fine with me.
476  //
477  EPETRA_CHK_ERR(H_->Solve(false,false,false,*Xcopy,Y));
478  EPETRA_CHK_ERR(H_->Solve(false,true,false,Y,Y));
479 
480 #ifdef IFPACK_FLOPCOUNTERS
481  // these are global flop count
482  ApplyInverseFlops_ += 4.0 * GlobalNonzeros_;
483 #endif
484 
485  ++NumApplyInverse_;
486  ApplyInverseTime_ += Time_.ElapsedTime();
487 
488  return(0);
489 }
490 //=============================================================================
491 // This function finds X such that LDU Y = X or U(trans) D L(trans) Y = X for multiple RHS
492 int Ifpack_ICT::Apply(const Epetra_MultiVector& X,
493  Epetra_MultiVector& Y) const
494 {
495 
496  IFPACK_CHK_ERR(-98);
497 }
498 
499 //=============================================================================
500 double Ifpack_ICT::Condest(const Ifpack_CondestType CT,
501  const int MaxIters, const double Tol,
502  Epetra_RowMatrix* Matrix_in)
503 {
504  if (!IsComputed()) // cannot compute right now
505  return(-1.0);
506 
507  // NOTE: this is computing the *local* condest
508  if (Condest_ == -1.0)
509  Condest_ = Ifpack_Condest(*this, CT, MaxIters, Tol, Matrix_in);
510 
511  return(Condest_);
512 }
513 
514 //=============================================================================
515 std::ostream&
516 Ifpack_ICT::Print(std::ostream& os) const
517 {
518  if (!Comm().MyPID()) {
519  os << endl;
520  os << "================================================================================" << endl;
521  os << "Ifpack_ICT: " << Label() << endl << endl;
522  os << "Level-of-fill = " << LevelOfFill() << endl;
523  os << "Absolute threshold = " << AbsoluteThreshold() << endl;
524  os << "Relative threshold = " << RelativeThreshold() << endl;
525  os << "Relax value = " << RelaxValue() << endl;
526  os << "Condition number estimate = " << Condest() << endl;
527  os << "Global number of rows = " << Matrix().NumGlobalRows64() << endl;
528  if (IsComputed_) {
529  os << "Number of nonzeros of H = " << H_->NumGlobalNonzeros64() << endl;
530  os << "nonzeros / rows = "
531  << 1.0 * H_->NumGlobalNonzeros64() / H_->NumGlobalRows64() << endl;
532  }
533  os << endl;
534  os << "Phase # calls Total Time (s) Total MFlops MFlops/s" << endl;
535  os << "----- ------- -------------- ------------ --------" << endl;
536  os << "Initialize() " << std::setw(5) << NumInitialize()
537  << " " << std::setw(15) << InitializeTime()
538  << " 0.0 0.0" << endl;
539  os << "Compute() " << std::setw(5) << NumCompute()
540  << " " << std::setw(15) << ComputeTime()
541  << " " << std::setw(15) << 1.0e-6 * ComputeFlops();
542  if (ComputeTime() != 0.0)
543  os << " " << std::setw(15) << 1.0e-6 * ComputeFlops() / ComputeTime() << endl;
544  else
545  os << " " << std::setw(15) << 0.0 << endl;
546  os << "ApplyInverse() " << std::setw(5) << NumApplyInverse()
547  << " " << std::setw(15) << ApplyInverseTime()
548  << " " << std::setw(15) << 1.0e-6 * ApplyInverseFlops();
549  if (ApplyInverseTime() != 0.0)
550  os << " " << std::setw(15) << 1.0e-6 * ApplyInverseFlops() / ApplyInverseTime() << endl;
551  else
552  os << " " << std::setw(15) << 0.0 << endl;
553  os << "================================================================================" << endl;
554  os << endl;
555  }
556 
557 
558  return(os);
559 }
double RelaxValue() const
Returns the relaxation value.
Definition: Ifpack_ICT.h:317
virtual ostream & Print(std::ostream &os) const
Prints basic information on iostream. This function is used by operator<<.
Definition: Ifpack_ICT.cpp:516
double LevelOfFill() const
Returns the level-of-fill.
Definition: Ifpack_ICT.h:299
const Epetra_Comm & Comm() const
Returns the Epetra_BlockMap object associated with the range of this matrix operator.
Definition: Ifpack_ICT.h:224
int SetParameters(Teuchos::ParameterList &parameterlis)
Set parameters using a Teuchos::ParameterList object.
Definition: Ifpack_ICT.cpp:104
bool IsInitialized() const
Returns true is the preconditioner has been successfully initialized.
Definition: Ifpack_ICT.h:116
int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Ifpack_ICT forward/back solve on a Epetra_MultiVector X in Y...
Definition: Ifpack_ICT.cpp:453
string Ifpack_toString(const int &x)
Convertes an integer to string.
int Initialize()
Initialize L and U with values from user matrix A.
Definition: Ifpack_ICT.cpp:136
double DropTolerance() const
Returns the drop threshold.
Definition: Ifpack_ICT.h:323
Ifpack_ICT(const Epetra_RowMatrix *A)
Ifpack_ICT constuctor with variable number of indices per row.
Definition: Ifpack_ICT.cpp:63
virtual int NumInitialize() const
Returns the number of calls to Initialize().
Definition: Ifpack_ICT.h:242
virtual int NumCompute() const
Returns the number of calls to Compute().
Definition: Ifpack_ICT.h:248
double Condest() const
Returns the computed condition number estimate, or -1.0 if not computed.
Definition: Ifpack_ICT.h:175
virtual double ComputeTime() const
Returns the time spent in Compute().
Definition: Ifpack_ICT.h:266
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
Definition: Ifpack_ICT.h:272
virtual double ComputeFlops() const
Returns the number of flops in all applications of Compute().
Definition: Ifpack_ICT.h:284
bool IsComputed() const
If factor is completed, this query returns true, otherwise it returns false.
Definition: Ifpack_ICT.h:142
double RelativeThreshold() const
Returns the relative threshold.
Definition: Ifpack_ICT.h:311
const Epetra_RowMatrix & Matrix() const
Returns a reference to the matrix to be preconditioned.
Definition: Ifpack_ICT.h:110
virtual double ApplyInverseFlops() const
Returns the number of flops in all applications of ApplyInverse().
Definition: Ifpack_ICT.h:290
virtual double InitializeTime() const
Returns the time spent in Initialize().
Definition: Ifpack_ICT.h:260
int Compute()
Compute IC factor U using the specified graph, diagonal perturbation thresholds and relaxation parame...
Definition: Ifpack_ICT.cpp:436
virtual int NumApplyInverse() const
Returns the number of calls to ApplyInverse().
Definition: Ifpack_ICT.h:254
double AbsoluteThreshold() const
Returns the absolute threshold.
Definition: Ifpack_ICT.h:305
virtual ~Ifpack_ICT()
Ifpack_ICT Destructor.
Definition: Ifpack_ICT.cpp:91