/usr/include/trilinos/EpetraExt_CrsSingletonFilter_LinearProblem.h is in libtrilinos-epetraext-dev 12.12.1-5.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 | //@HEADER
// ***********************************************************************
//
// EpetraExt: Epetra Extended - Linear Algebra Services Package
// Copyright (2011) Sandia Corporation
//
// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
// the U.S. Government retains certain rights in this software.
//
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are
// met:
//
// 1. Redistributions of source code must retain the above copyright
// notice, this list of conditions and the following disclaimer.
//
// 2. Redistributions in binary form must reproduce the above copyright
// notice, this list of conditions and the following disclaimer in the
// documentation and/or other materials provided with the distribution.
//
// 3. Neither the name of the Corporation nor the names of the
// contributors may be used to endorse or promote products derived from
// this software without specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
//
// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
//
// ***********************************************************************
//@HEADER
#ifndef _EpetraExt_LINEARPROBLEM_CRSSINGLETONFILTER_H_
#define _EpetraExt_LINEARPROBLEM_CRSSINGLETONFILTER_H_
#include "Epetra_ConfigDefs.h"
#include "Epetra_Object.h"
#include "Epetra_CrsMatrix.h"
#include "Epetra_MapColoring.h"
#include "Epetra_SerialDenseVector.h"
#include "EpetraExt_Transform.h"
#include "Teuchos_RCP.hpp"
class Epetra_LinearProblem;
class Epetra_Map;
class Epetra_MultiVector;
class Epetra_Import;
class Epetra_Export;
class Epetra_IntVector;
namespace EpetraExt {
//! Epetra_CrsSingletonFilter: A class for explicitly eliminating matrix rows and columns.
/*! The Epetra_CrsSingletonFilter class takes an existing Epetra_LinearProblem object, analyzes
it structure and explicitly eliminates singleton rows and columns from the matrix and appropriately
modifies the RHS and LHS of the linear problem. The result of this process is a reduced system of equations
that is itself an Epetra_LinearProblem object. The reduced system can then be solved using any solver that
is understands an Epetra_LinearProblem. The solution for the full system is obtained by calling ComputeFullSolution().
Singleton rows are defined to be rows that have a single nonzero entry in the matrix. The equation associated with
this row can be explicitly eliminated because it involved only one variable. For example if row i has a single nonzero
value in column j, call it A(i,j), we can explicitly solve for x(j) = b(i)/A(i,j), where b(i) is the ith entry of the RHS
and x(j) is the jth entry of the LHS.
Singleton columns are defined to be columns that have a single nonzero entry in the matrix. The variable associated
with this column is fully dependent, meaning that the solution for all other variables does not depend on it. If this
entry is A(i,j) then the ith row and jth column can be removed from the system and x(j) can be solved after the solution
for all other variables is determined.
By removing singleton rows and columns, we can often produce a reduced system that is smaller and far less dense, and in
general having better numerical properties.
The basic procedure for using this class is as follows:
<ol>
<li> Construct full problem: Construct and Epetra_LinearProblem containing the "full" matrix, RHS and LHS. This is
done outside of Epetra_CrsSingletonFilter class.
Presumably, you have some reason to believe that this system may contain singletons.
<li> Construct an Epetra_CrsSingletonFilter instance: Constructor needs no arguments.
<li> Analyze matrix: Invoke the Analyze() method, passing in the Epetra_RowMatrix object from your full linear
problem mentioned in the first step above.
<li> Go/No Go decision to construct reduced problem:
Query the results of the Analyze method using the SingletonsDetected() method. This method
returns "true" if there were singletons found in the matrix. You can also query any of the other methods
in the Filter Statistics section to determine if you want to proceed with the construction of the reduced system.
<li> Construct reduced problem:
If, in the previous step, you determine that you want to proceed with the construction of the reduced problem,
you should next call the ConstructReducedProblem() method, passing in the full linear problem object from the first
step. This method will use the information from the Analyze() method to construct a reduce problem that has
explicitly eliminated the singleton rows, solved for the corresponding LHS values and updated the RHS. This
step will also remove singleton columns from the reduced system. Once the solution of the reduced problem is
is computed (via any solver that understands an Epetra_LinearProblem), you should call the ComputeFullSolution()
method to compute the LHS values assocaited with the singleton columns.
<li> Solve reduced problem: Obtain a pointer to the reduced problem using the ReducedProblem() method.
Using the solver of your choice, solve the reduced system.
<li> Compute solution to full problem: Once the solution the reduced problem is determined, the ComputeFullSolution()
method will place the reduced solution values into the appropriate locations of the full solution LHS and then
compute the values associated with column singletons. At this point, you have a complete solution to the original
full problem.
<li> Solve a subsequent full problem that differs from the original problem only in values: It is often the case that the
structure of a problem will be the same for a sequence of linear problems. In this case, the UpdateReducedProblem()
method can be useful. After going through the above process one time, if you have a linear problem that is structural
\e identical to the previous problem, you can minimize memory and time costs by using the UpdateReducedProblem()
method, passing in the subsequent problem. Once you have called the UpdateReducedProblem() method, you can then
solve the reduce problem problem as you wish, and then compute the full solution as before. The pointer generated
by ReducedProblem() will not change when UpdateReducedProblem() is called.
</ol>
*/
class LinearProblem_CrsSingletonFilter : public SameTypeTransform<Epetra_LinearProblem> {
public:
//@{ \name Constructors/Destructor.
//! Epetra_CrsSingletonFilter default constructor.
LinearProblem_CrsSingletonFilter( bool verbose = false );
//! Epetra_CrsSingletonFilter Destructor
virtual ~LinearProblem_CrsSingletonFilter();
//@}
///
NewTypeRef operator()( OriginalTypeRef orig );
///
bool analyze( OriginalTypeRef orig );
///
NewTypeRef construct();
///
bool fwd();
///
bool rvs();
//@{ \name Analyze methods.
//! Analyze the input matrix, removing row/column pairs that have singletons.
/*! Analyzes the user's input matrix to determine rows and columns that should be explicitly
eliminated to create the reduced system. Look for rows and columns that have single entries.
These rows/columns
can easily be removed from the problem.
The results of calling this method are two MapColoring objects accessible via RowMapColors() and
ColMapColors() accessor methods. All rows/columns that would be eliminated in the reduced system
have a color of 1 in the corresponding RowMapColors/ColMapColors object. All kept rows/cols have a
color of 0.
*/
int Analyze(Epetra_RowMatrix * FullMatrix);
//! Returns true if singletons were detected in this matrix (must be called after Analyze() to be effective).
bool SingletonsDetected() const {if (!AnalysisDone_) return(false); else return(NumSingletons()>0);};
//@}
//@{ \name Reduce methods.
//! Return a reduced linear problem based on results of Analyze().
/*! Creates a new Epetra_LinearProblem object based on the results of the Analyze phase. A pointer
to the reduced problem is obtained via a call to ReducedProblem().
\return Error code, set to 0 if no error.
*/
int ConstructReducedProblem(Epetra_LinearProblem * Problem);
//! Update a reduced linear problem using new values.
/*! Updates an existing Epetra_LinearProblem object using new matrix, LHS and RHS values. The matrix
structure must be \e identical to the matrix that was used to construct the original reduced problem.
\return Error code, set to 0 if no error.
*/
int UpdateReducedProblem(Epetra_LinearProblem * Problem);
//@}
//@{ \name Methods to construct Full System Solution.
//! Compute a solution for the full problem using the solution of the reduced problem, put in LHS of FullProblem().
/*! After solving the reduced linear system, this method can be called to compute the
solution to the original problem, assuming the solution for the reduced system is valid. The solution of the
unreduced, original problem will be in the LHS of the original Epetra_LinearProblem.
*/
int ComputeFullSolution();
//@}
//@{ \name Filter Statistics.
//! Return number of rows that contain a single entry, returns -1 if Analysis not performed yet.
int NumRowSingletons() const {return(NumGlobalRowSingletons_);};
//! Return number of columns that contain a single entry that are \e not associated with singleton row, returns -1 if Analysis not performed yet.
int NumColSingletons() const {return(NumGlobalColSingletons_);};
//! Return total number of singletons detected, returns -1 if Analysis not performed yet.
/*! Return total number of singletons detected across all processors. This method will not return a
valid result until after the Analyze() method is called. The dimension of the reduced system can
be computed by subtracting this number from dimension of full system.
\warning This method returns -1 if Analyze() method has not been called.
*/
int NumSingletons() const {return(NumColSingletons()+NumRowSingletons());};
//! Returns ratio of reduced system to full system dimensions, returns -1.0 if reduced problem not constructed.
double RatioOfDimensions() const {return(RatioOfDimensions_);};
//! Returns ratio of reduced system to full system nonzero count, returns -1.0 if reduced problem not constructed.
double RatioOfNonzeros() const {return(RatioOfNonzeros_);};
//@}
//@{ \name Attribute Access Methods.
//! Returns pointer to the original unreduced Epetra_LinearProblem.
Epetra_LinearProblem * FullProblem() const {return(FullProblem_);};
//! Returns pointer to the derived reduced Epetra_LinearProblem.
Epetra_LinearProblem * ReducedProblem() const {return(ReducedProblem_.get());};
//! Returns pointer to Epetra_CrsMatrix from full problem.
Epetra_RowMatrix * FullMatrix() const {return(FullMatrix_);};
//! Returns pointer to Epetra_CrsMatrix from full problem.
Epetra_CrsMatrix * ReducedMatrix() const {return(ReducedMatrix_.get());};
//! Returns pointer to Epetra_MapColoring object: color 0 rows are part of reduced system.
Epetra_MapColoring * RowMapColors() const {return(RowMapColors_);};
//! Returns pointer to Epetra_MapColoring object: color 0 columns are part of reduced system.
Epetra_MapColoring * ColMapColors() const {return(ColMapColors_);};
//! Returns pointer to Epetra_Map describing the reduced system row distribution.
Epetra_Map * ReducedMatrixRowMap() const {return(ReducedMatrixRowMap_);};
//! Returns pointer to Epetra_Map describing the reduced system column distribution.
Epetra_Map * ReducedMatrixColMap() const {return(ReducedMatrixColMap_);};
//! Returns pointer to Epetra_Map describing the domain map for the reduced system.
Epetra_Map * ReducedMatrixDomainMap() const {return(ReducedMatrixDomainMap_);};
//! Returns pointer to Epetra_Map describing the range map for the reduced system.
Epetra_Map * ReducedMatrixRangeMap() const {return(ReducedMatrixRangeMap_);};
//@}
protected:
// This pointer will be zero if full matrix is not a CrsMatrix.
Epetra_CrsMatrix * FullCrsMatrix() const {return(FullCrsMatrix_);};
const Epetra_Map & FullMatrixRowMap() const {return(FullMatrix()->RowMatrixRowMap());};
const Epetra_Map & FullMatrixColMap() const {return(FullMatrix()->RowMatrixColMap());};
const Epetra_Map & FullMatrixDomainMap() const {return((FullMatrix()->OperatorDomainMap()));};
const Epetra_Map & FullMatrixRangeMap() const {return((FullMatrix()->OperatorRangeMap()));};
void InitializeDefaults();
int ComputeEliminateMaps();
int Setup(Epetra_LinearProblem * Problem);
int InitFullMatrixAccess();
int GetRow(int Row, int & NumIndices, int * & Indices);
int GetRow(int Row, int & NumIndices, double * & Values, int * & Indices);
#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
int GetRowGCIDs(int Row, int & NumIndices, double * & Values, int * & GlobalIndices);
#endif
#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
int GetRowGCIDs(int Row, int & NumIndices, double * & Values, long long * & GlobalIndices);
#endif
int CreatePostSolveArrays(const Epetra_IntVector & RowIDs,
const Epetra_MapColoring & RowMapColors,
const Epetra_IntVector & ColProfiles,
const Epetra_IntVector & NewColProfiles,
const Epetra_IntVector & ColHasRowWithSingleton);
int ConstructRedistributeExporter(Epetra_Map * SourceMap, Epetra_Map * TargetMap,
Epetra_Export * & RedistributeExporter,
Epetra_Map * & RedistributeMap);
Epetra_LinearProblem * FullProblem_;
Teuchos::RCP<Epetra_LinearProblem> ReducedProblem_;
Epetra_RowMatrix * FullMatrix_;
Epetra_CrsMatrix * FullCrsMatrix_;
Teuchos::RCP<Epetra_CrsMatrix> ReducedMatrix_;
Epetra_MultiVector * ReducedRHS_;
Epetra_MultiVector * ReducedLHS_;
Epetra_Map * ReducedMatrixRowMap_;
Epetra_Map * ReducedMatrixColMap_;
Epetra_Map * ReducedMatrixDomainMap_;
Epetra_Map * ReducedMatrixRangeMap_;
Epetra_Map * OrigReducedMatrixDomainMap_;
Epetra_Import * Full2ReducedRHSImporter_;
Epetra_Import * Full2ReducedLHSImporter_;
Epetra_Export * RedistributeDomainExporter_;
int * ColSingletonRowLIDs_;
int * ColSingletonColLIDs_;
int * ColSingletonPivotLIDs_;
double * ColSingletonPivots_;
int AbsoluteThreshold_;
double RelativeThreshold_;
int NumMyRowSingletons_;
int NumMyColSingletons_;
int NumGlobalRowSingletons_;
int NumGlobalColSingletons_;
double RatioOfDimensions_;
double RatioOfNonzeros_;
bool HaveReducedProblem_;
bool UserDefinedEliminateMaps_;
bool AnalysisDone_;
bool SymmetricElimination_;
Epetra_MultiVector * tempExportX_;
Epetra_MultiVector * tempX_;
Epetra_MultiVector * tempB_;
Epetra_MultiVector * RedistributeReducedLHS_;
int * Indices_int_;
#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
long long * Indices_LL_;
#endif
Epetra_SerialDenseVector Values_;
Epetra_MapColoring * RowMapColors_;
Epetra_MapColoring * ColMapColors_;
bool FullMatrixIsCrsMatrix_;
int MaxNumMyEntries_;
bool verbose_;
private:
//! Copy constructor (defined as private so it is unavailable to user).
LinearProblem_CrsSingletonFilter(const LinearProblem_CrsSingletonFilter & Problem){};
template<typename int_type>
int TConstructReducedProblem(Epetra_LinearProblem * Problem);
template<typename int_type>
int TUpdateReducedProblem(Epetra_LinearProblem * Problem);
template<typename int_type>
int TConstructRedistributeExporter(Epetra_Map * SourceMap, Epetra_Map * TargetMap,
Epetra_Export * & RedistributeExporter,
Epetra_Map * & RedistributeMap);
};
} //namespace EpetraExt
#endif /* _EpetraExt_LINEARPROBLEM_CRSSINGLETONFILTER_H_ */
|