Thyra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Thyra_DefaultDiagonalLinearOpWithSolve_def.hpp
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Thyra: Interfaces and Support for Abstract Numerical Algorithms
5 // Copyright (2004) 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 Roscoe A. Bartlett (bartlettra@ornl.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
42 #ifndef THYRA_DIAGONAL_LINEAR_OP_WITH_SOLVE_HPP
43 #define THYRA_DIAGONAL_LINEAR_OP_WITH_SOLVE_HPP
44 
45 #include "Thyra_DefaultDiagonalLinearOpWithSolve_decl.hpp"
46 #include "Thyra_DefaultDiagonalLinearOp.hpp"
47 #include "Thyra_MultiVectorStdOps.hpp"
48 #include "Thyra_VectorBase.hpp"
49 #include "Thyra_TestingTools.hpp" // ToDo: I need to have a better way to get eps()!
50 #include "Teuchos_Assert.hpp"
51 
52 
53 namespace Thyra {
54 
55 
56 // Constructors/initializers/accessors
57 
58 
59 template<class Scalar>
61 {}
62 
63 
64 template<class Scalar>
66  const RCP<const VectorBase<Scalar> > &diag
67  )
68 {
69  this->initialize(diag);
70 }
71 
72 
73 // protected
74 
75 
76 // Overridden from LinearOpWithSolveBase
77 
78 
79 template<class Scalar>
80 bool
82  EOpTransp M_trans) const
83 {
85  return (ST::isComplex ? M_trans==NOTRANS || M_trans==TRANS : true);
86 }
87 
88 
89 template<class Scalar>
90 bool
92  EOpTransp M_trans, const SolveMeasureType& /* solveMeasureType */) const
93 {
94  return this->solveSupportsImpl(M_trans); // I am a direct solver!
95 }
96 
97 
98 template<class Scalar>
101  const EOpTransp transp,
102  const MultiVectorBase<Scalar> &B,
103  const Ptr<MultiVectorBase<Scalar> > &X,
104  const Ptr<const SolveCriteria<Scalar> > solveCriteria
105  ) const
106 {
107 
108 #ifdef THYRA_DEBUG
109  TEUCHOS_ASSERT(this->solveSupportsImpl(transp));
110 #else
111  (void)transp;
112 #endif
113 
115 
116  assign(X, ST::zero());
117 
118  const Ordinal numCols = B.domain()->dim();
119  SolveStatus<Scalar> overallSolveStatus;
120 
121  for (Ordinal col_j = 0; col_j < numCols; ++col_j) {
122 
123  const RCP<const VectorBase<Scalar> > b = B.col(col_j);
124  const RCP<VectorBase<Scalar> > x = X->col(col_j);
125 
126  ele_wise_divide( ST::one(), *b, *this->getDiag(), x.ptr() );
127 
128  }
129 
130  SolveStatus<Scalar> solveStatus;
131  solveStatus.solveStatus =
132  (nonnull(solveCriteria) && !solveCriteria->solveMeasureType.useDefault()
135  return solveStatus;
136 
137 }
138 
139 
140 } // end namespace Thyra
141 
142 
143 #endif // THYRA_DIAGONAL_LINEAR_OP_WITH_SOLVE_HPP
void initialize(int *argc, char ***argv)
EOpTransp
Enumeration for determining how a linear operator is applied. `*.
Use the non-transposed operator.
The final solution status is unknown but he solve did not totally fail.
SolveStatus< Scalar > solveImpl(const EOpTransp transp, const MultiVectorBase< Scalar > &B, const Ptr< MultiVectorBase< Scalar > > &X, const Ptr< const SolveCriteria< Scalar > > solveCriteria) const
Use the transposed operator.
Teuchos::Ordinal Ordinal
Type for the dimension of a vector space. `*.
Interface for a collection of column vectors called a multi-vector.
Ptr< T > ptr() const
RCP< const VectorBase< Scalar > > col(Ordinal j) const
Calls colImpl().
static ScalarMag unknownTolerance()
Abstract interface for finite-dimensional dense vectors.
Simple struct for the return status from a solve.
bool solveSupportsSolveMeasureTypeImpl(EOpTransp M_trans, const SolveMeasureType &solveMeasureType) const
ESolveStatus solveStatus
The return status of the solve.
virtual RCP< const VectorSpaceBase< Scalar > > domain() const =0
Return a smart pointer for the domain space for this operator.
bool nonnull(const boost::shared_ptr< T > &p)
ScalarMag achievedTol
The maximum final tolerance actually achieved by the (block) linear solve. A value of unknownToleranc...
The requested solution criteria has likely been achieved.
#define TEUCHOS_ASSERT(assertion_test)
Simple struct that defines the requested solution criteria for a solve.