Teko Version of the Day
Loading...
Searching...
No Matches
Teko_DiagnosticLinearOp.cpp
1/*
2// @HEADER
3//
4// ***********************************************************************
5//
6// Teko: A package for block and physics based preconditioning
7// Copyright 2010 Sandia Corporation
8//
9// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
10// the U.S. Government retains certain rights in this software.
11//
12// Redistribution and use in source and binary forms, with or without
13// modification, are permitted provided that the following conditions are
14// met:
15//
16// 1. Redistributions of source code must retain the above copyright
17// notice, this list of conditions and the following disclaimer.
18//
19// 2. Redistributions in binary form must reproduce the above copyright
20// notice, this list of conditions and the following disclaimer in the
21// documentation and/or other materials provided with the distribution.
22//
23// 3. Neither the name of the Corporation nor the names of the
24// contributors may be used to endorse or promote products derived from
25// this software without specific prior written permission.
26//
27// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
28// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
29// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
30// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
31// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
32// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
33// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
34// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
35// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
36// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
37// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
38//
39// Questions? Contact Eric C. Cyr (eccyr@sandia.gov)
40//
41// ***********************************************************************
42//
43// @HEADER
44
45*/
46
48
49#include "Teuchos_TimeMonitor.hpp"
50
51#include "Thyra_MultiVectorStdOps.hpp"
52
53namespace Teko {
54
59DiagnosticLinearOp::DiagnosticLinearOp(const Teuchos::RCP<std::ostream> & ostrm, const ModifiableLinearOp & A,const std::string & diagnosticString)
60 : outputStream_(ostrm), wrapOpA_(A), wrapOpA_lo_(A), diagString_(diagnosticString), timer_(diagnosticString)
61{
62}
63
64DiagnosticLinearOp::DiagnosticLinearOp(const Teuchos::RCP<std::ostream> & ostrm, const LinearOp & A,const std::string & diagnosticString)
65 : outputStream_(ostrm), wrapOpA_(Teuchos::null), wrapOpA_lo_(A), diagString_(diagnosticString), timer_(diagnosticString)
66{
67}
68
73DiagnosticLinearOp::DiagnosticLinearOp(const Teuchos::RCP<std::ostream> & ostrm,const LinearOp & fwdOp,
74 const ModifiableLinearOp & A,const std::string & diagnosticString)
75 : outputStream_(ostrm), wrapOpA_(A), wrapOpA_lo_(A), fwdOp_(fwdOp), diagString_(diagnosticString), timer_(diagnosticString)
76{
77}
78
80{
81 double elapsedTime = totalTime();
82 int applications = numApplications();
83
84 (*outputStream_) << "DiagnosticLinearOp \"" << diagString_ << "\": "
85 << "elapsed = " << elapsedTime << ", "
86 << "applications = " << applications << ", ";
87 if(applications>0)
88 (*outputStream_) << "timer/app = " << elapsedTime / double(applications) << std::endl;
89 else
90 (*outputStream_) << "timer/app = " << "none" << std::endl;
91}
92
105void DiagnosticLinearOp::implicitApply(const MultiVector & x, MultiVector & y,
106 const double alpha, const double beta) const
107{
108 Teko_DEBUG_SCOPE("DiagnosticLinearOp::implicityApply",10);
109
110 // start timer on construction, end on destruction
111 Teuchos::TimeMonitor monitor(timer_,false);
112
113 MultiVector z; // for temporary storage dealing with nozero beta
114 if(beta!=0.0)
115 z = deepcopy(y);
116
117 wrapOpA_lo_->apply(Thyra::NOTRANS,*x,y.ptr(),alpha,beta);
118
119 // print residual if there is a fwd Op
120 bool printResidual = (fwdOp_!=Teuchos::null);
121 if(printResidual) {
122 // compute residual
123 MultiVector residual = Teko::deepcopy(x);
124 // fwdOp_->apply(Thyra::NOTRANS,*y,residual.ptr(),-1.0,1.0);
125
126 fwdOp_->apply(Thyra::NOTRANS,*y,residual.ptr(),-1.0,alpha);
127 if(beta!=0.0)
128 fwdOp_->apply(Thyra::NOTRANS,*z,residual.ptr(),beta,1.0);
129
130 // calculate norms
131 std::vector<double> norms(y->domain()->dim()); // size of column count
132 std::vector<double> rhsNorms(x->domain()->dim()); // size of column count
133 Thyra::norms_2<double>(*residual,Teuchos::arrayViewFromVector(norms));
134 Thyra::norms_2<double>(*x,Teuchos::arrayViewFromVector(rhsNorms));
135
136 // print out residual norms
137 (*outputStream_) << "DiagnosticLinearOp \"" << diagString_ << "\": residual = [";
138 for(std::size_t i=0;i<norms.size();++i)
139 (*outputStream_) << " " << std::scientific << std::setprecision(4) << norms[i]/rhsNorms[i];// << " (" <<rhsNorms[i]<<") ";
140 (*outputStream_) << " ]" << std::endl;
141
142 residualNorm_ = norms[0];
143 }
144}
145
146} // end namespace Teko
MultiVector deepcopy(const MultiVector &v)
Perform a deep copy of the vector.
virtual ~DiagnosticLinearOp()
Destructor prints out timing information about this operator.
virtual void implicitApply(const MultiVector &x, MultiVector &y, const double alpha=1.0, const double beta=0.0) const
Perform a matrix vector multiply with this operator.