00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039 #ifndef NOX_EPETRA_BROYDENOPERATOR_H
00040 #define NOX_EPETRA_BROYDENOPERATOR_H
00041
00042 #include "NOX.H"
00043 #include "NOX_Common.H"
00044 #include "NOX_Epetra.H"
00045 #include "Teuchos_ParameterList.hpp"
00046 #include "NOX_Solver_PrePostOperator.H"
00047 #include "NOX_Abstract_PrePostOperator.H"
00048
00049 #include "Epetra_Vector.h"
00050 #include "Epetra_CrsGraph.h"
00051 #include "Epetra_CrsMatrix.h"
00052
00053 namespace NOX {
00054
00055 namespace Epetra {
00056
00095 class BroydenOperator : public Epetra_RowMatrix,
00096 public NOX::Abstract::PrePostOperator,
00097 public NOX::Epetra::Interface::Jacobian,
00098 public NOX::Epetra::Interface::Preconditioner
00099 {
00100
00101 public:
00102
00103 class ReplacementInterface
00104 {
00105 public:
00106
00107 enum FILL_TYPE { JACOBIAN ,
00108 PRECONDITIONER
00109 };
00110
00111 virtual ~ReplacementInterface()
00112 { }
00113
00114 virtual Teuchos::RCP<const Epetra_CrsMatrix> getReplacementValuesMatrix( const Epetra_Vector & x, FILL_TYPE ) = 0;
00115
00116 };
00117
00119 BroydenOperator(Teuchos::ParameterList & nlParams,
00120 const Teuchos::RCP<NOX::Utils>& utils,
00121 Epetra_Vector & solnVec,
00122 const Teuchos::RCP<Epetra_CrsMatrix>& broydMat0,
00123 bool verbose = false);
00124
00126 BroydenOperator(const BroydenOperator &);
00127
00129 virtual ~BroydenOperator();
00130
00132 virtual const char* Label () const;
00133
00135 virtual int SetUseTranspose(bool UseTranspose);
00136
00138 virtual int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00139
00141 virtual int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00142
00144 virtual bool UseTranspose() const;
00145
00147 virtual bool HasNormInf() const;
00148
00150 virtual const Epetra_Map & OperatorDomainMap() const;
00151
00153 virtual const Epetra_Map & OperatorRangeMap() const;
00154
00156 virtual bool Filled() const;
00157
00159 virtual int NumMyRowEntries(int MyRow, int & NumEntries) const;
00160
00162 virtual int MaxNumEntries() const;
00163
00165 virtual int ExtractMyRowCopy(int MyRow, int Length, int & NumEntries, double *Values, int * Indices) const;
00166
00168 virtual int ExtractDiagonalCopy(Epetra_Vector & Diagonal) const;
00169
00171 virtual int Multiply(bool TransA, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00172
00174 virtual int Solve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00175
00177 virtual int InvRowSums(Epetra_Vector& x) const;
00178
00180 virtual int LeftScale(const Epetra_Vector& x);
00181
00183 virtual int InvColSums(Epetra_Vector& x) const;
00184
00186 virtual int RightScale(const Epetra_Vector& x);
00187
00189 virtual double NormInf() const;
00190
00192 virtual double NormOne() const;
00193
00195 virtual int NumGlobalNonzeros() const;
00196
00198 virtual int NumGlobalRows() const;
00199
00201 virtual int NumGlobalCols() const;
00202
00204 virtual int NumGlobalDiagonals() const;
00205
00207 virtual int NumMyNonzeros() const;
00208
00210 virtual int NumMyRows() const;
00211
00213 virtual int NumMyCols() const;
00214
00216 virtual int NumMyDiagonals() const;
00217
00219 virtual bool LowerTriangular() const;
00220
00222 virtual bool UpperTriangular() const;
00223
00225 virtual const Epetra_Comm & Comm() const;
00226
00228 virtual const Epetra_Map & RowMatrixRowMap() const;
00229
00231 virtual const Epetra_Map & RowMatrixColMap() const;
00232
00234 virtual const Epetra_Import * RowMatrixImporter() const;
00235
00237 virtual const Epetra_BlockMap& Map() const;
00238
00240 void setStepVector ( Epetra_Vector & vec );
00241
00243 void setStepVector ( NOX::Epetra::Vector & vec );
00244
00246 void setYieldVector( Epetra_Vector & vec );
00247
00249 void setYieldVector( NOX::Epetra::Vector & vec );
00250
00252 bool computeSparseBroydenUpdate();
00253
00255 void removeEntriesFromBroydenUpdate( const Epetra_CrsGraph & graph);
00256
00258 const Epetra_CrsMatrix & getBroydenMatrix()
00259 { return *crsMatrix; };
00260
00262 void resetBroydenMatrix( const Epetra_CrsMatrix & mat );
00263
00265 void addReplacementInterface(ReplacementInterface * i )
00266 { replacementInterfaces.push_back(i); }
00267
00275
00276 virtual bool isStep() const;
00277 virtual bool isYield() const;
00278 virtual bool isBroyden() const;
00279
00281
00282 #ifdef HAVE_NOX_DEBUG
00283
00284 void outputActiveEntries();
00285 #endif
00286
00287
00288 protected:
00289
00291 virtual bool initialize( Teuchos::ParameterList & nlParams, const Epetra_Vector & x );
00292
00294 virtual bool computeJacobian( const Epetra_Vector & x, Epetra_Operator& Jac);
00295
00297 virtual bool computePreconditioner( const Epetra_Vector &,
00298 Epetra_Operator& Prec,
00299 Teuchos::ParameterList * params = 0 );
00300
00302 virtual void runPreSolve( const NOX::Solver::Generic & solver);
00303
00305 virtual void runPreIterate( const NOX::Solver::Generic & solver);
00306
00308 virtual void runPostIterate( const NOX::Solver::Generic & solver);
00309
00311 virtual void runPostSolve( const NOX::Solver::Generic & solver);
00312
00314 void replaceBroydenMatrixValues( const Epetra_CrsMatrix & mat);
00315
00316 protected:
00317
00318 bool verbose;
00319
00320 Teuchos::RCP<NOX::Epetra::Vector> stepVec ;
00321
00322 Teuchos::RCP<NOX::Epetra::Vector> yieldVec ;
00323
00324 Teuchos::RCP<NOX::Epetra::Vector> workVec ;
00325
00326 Teuchos::RCP<NOX::Epetra::Vector> oldX ;
00327
00328 Teuchos::RCP<NOX::Epetra::Vector> oldF ;
00329
00330 Teuchos::RCP<Epetra_CrsMatrix> crsMatrix;
00331
00332 Teuchos::RCP<NOX::Epetra::Interface::Jacobian> jacIntPtr;
00333 Teuchos::RCP<Epetra_CrsMatrix> jacMatrixPtr;
00334
00335 Teuchos::RCP<NOX::Epetra::Interface::Preconditioner> precIntPtr;
00336 Teuchos::RCP<Epetra_CrsMatrix> precMatrixPtr;
00337
00339 Teuchos::ParameterList & nlParams;
00340
00342 const Teuchos::RCP<NOX::Utils>& utils;
00343
00345 NOX::Solver::PrePostOperator prePostOperator;
00346
00348 string label;
00349
00354 bool isValidStep;
00355 bool isValidYield;
00356 bool isValidBroyden;
00358
00360 vector<bool> entriesRemoved;
00361
00363 map< int, list<int> > retainedEntries;
00364
00366 vector<ReplacementInterface *> replacementInterfaces;
00367
00368 };
00369 }
00370 }
00371
00372 #endif