Ifpack_DiagonalFilter.cpp

00001 /*@HEADER
00002 // ***********************************************************************
00003 //
00004 //       Ifpack: Object-Oriented Algebraic Preconditioner Package
00005 //                 Copyright (2002) Sandia Corporation
00006 //
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 //
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00025 //
00026 // ***********************************************************************
00027 //@HEADER
00028 */
00029 
00030 #include "Ifpack_ConfigDefs.h"
00031 #include "Ifpack_DiagonalFilter.h"
00032 #include "Epetra_ConfigDefs.h"
00033 #include "Epetra_RowMatrix.h"
00034 #include "Epetra_Comm.h"
00035 #include "Epetra_Map.h"
00036 #include "Epetra_MultiVector.h"
00037 #include "Epetra_Vector.h"
00038 
00039 //==============================================================================
00040 Ifpack_DiagonalFilter::Ifpack_DiagonalFilter(const Teuchos::RefCountPtr<Epetra_RowMatrix>& Matrix,
00041                          double AbsoluteThreshold,
00042                          double RelativeThreshold) :
00043   A_(Matrix),
00044   AbsoluteThreshold_(AbsoluteThreshold),
00045   RelativeThreshold_(RelativeThreshold)
00046 {
00047   Epetra_Time Time(Comm());
00048   
00049   pos_.resize(NumMyRows());
00050   val_.resize(NumMyRows());
00051   
00052   vector<int> Indices(MaxNumEntries());
00053   vector<double> Values(MaxNumEntries());
00054   int NumEntries;
00055   
00056   for (int MyRow = 0 ; MyRow < NumMyRows() ; ++MyRow) {
00057     
00058     pos_[MyRow] = -1;
00059     val_[MyRow] = 0.0;
00060     int ierr = A_->ExtractMyRowCopy(MyRow, MaxNumEntries(), NumEntries,
00061                     &Values[0], &Indices[0]);
00062     assert (ierr == 0);
00063     
00064     for (int i = 0 ; i < NumEntries ; ++i) {
00065       if (Indices[i] == MyRow) {
00066     pos_[MyRow] = i;
00067     val_[MyRow] = Values[i] * (RelativeThreshold_ - 1) +
00068       AbsoluteThreshold_ * EPETRA_SGN(Values[i]);
00069       }
00070       break;
00071     }
00072   }
00073   cout << "TIME = " << Time.ElapsedTime() << endl;
00074 }
00075 
00076 //==============================================================================
00077 int Ifpack_DiagonalFilter::
00078 ExtractMyRowCopy(int MyRow, int Length, int& NumEntries, 
00079          double* Values, int* Indices) const
00080 {
00081 
00082   IFPACK_CHK_ERR(A_->ExtractMyRowCopy(MyRow, Length, NumEntries,
00083                      Values,Indices));
00084 
00085   if (pos_[MyRow] != -1)
00086     Values[pos_[MyRow]] += val_[MyRow];
00087 
00088   return(0);
00089 }
00090 
00091 //==============================================================================
00092 int Ifpack_DiagonalFilter::
00093 Multiply(bool TransA, const Epetra_MultiVector& X, 
00094      Epetra_MultiVector& Y) const
00095 {
00096 
00097   if (X.NumVectors() != Y.NumVectors())
00098     IFPACK_CHK_ERR(-2);
00099 
00100   IFPACK_CHK_ERR(A_->Multiply(TransA, X, Y));
00101 
00102   for (int v = 0 ; v < X.NumVectors() ; ++v)
00103     for (int i = 0 ; i < NumMyRows() ; ++i)
00104       Y[v][i] += val_[i] * X[v][i];
00105       
00106 
00107   return(0);
00108 }

Generated on Wed May 12 21:30:18 2010 for IFPACK by  doxygen 1.4.7