Rythmos - Transient Integration for Differential Equations  Version of the Day
 All Classes Functions Variables Typedefs Pages
Rythmos_PointwiseInterpolationBufferAppender.hpp
1 //@HEADER
2 // ***********************************************************************
3 //
4 // Rythmos Package
5 // Copyright (2006) 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 // This library is free software; you can redistribute it and/or modify
11 // it under the terms of the GNU Lesser General Public License as
12 // published by the Free Software Foundation; either version 2.1 of the
13 // License, or (at your option) any later version.
14 //
15 // This library is distributed in the hope that it will be useful, but
16 // WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 // Lesser General Public License for more details.
19 //
20 // You should have received a copy of the GNU Lesser General Public
21 // License along with this library; if not, write to the Free Software
22 // Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23 // USA
24 // Questions? Contact Todd S. Coffey (tscoffe@sandia.gov)
25 //
26 // ***********************************************************************
27 //@HEADER
28 
29 #ifndef RYTHMOS_POINTWISE_INTERPOLATION_BUFFER_APPENDER_HPP
30 #define RYTHMOS_POINTWISE_INTERPOLATION_BUFFER_APPENDER_HPP
31 
32 #include "Rythmos_InterpolationBufferAppenderBase.hpp"
33 #include "Teuchos_ParameterListAcceptorDefaultBase.hpp"
34 #include "Teuchos_VerboseObjectParameterListHelpers.hpp"
35 #include "Rythmos_InterpolationBufferHelpers.hpp"
36 
37 
38 namespace Rythmos {
39 
40 
44 template<class Scalar>
46  : virtual public InterpolationBufferAppenderBase<Scalar>,
47  virtual public Teuchos::ParameterListAcceptorDefaultBase
48 {
49 public:
50 
52  typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType ScalarMag;
53 
57  void append(
58  const InterpolationBufferBase<Scalar>& interpBuffSource,
59  const TimeRange<Scalar>& range,
60  const Ptr<InterpolationBufferBase<Scalar> > &interpBuffSink
61  );
62 
65 
67  void describe(
68  Teuchos::FancyOStream &out,
69  const Teuchos::EVerbosityLevel verbLevel
70  ) const;
71 
73 
76 
78  void setParameterList(const RCP<ParameterList> &paramList);
79 
81  RCP<const ParameterList> getValidParameters() const;
82 
84 
85 };
86 
87 
88 
93 template<class Scalar>
94 RCP<PointwiseInterpolationBufferAppender<Scalar> >
96 {
97  return Teuchos::rcp(new PointwiseInterpolationBufferAppender<Scalar>);
98 }
99 
100 
101 //
102 // Implementations
103 //
104 
105 
106 template<class Scalar>
108  const InterpolationBufferBase<Scalar>& interpBuffSource,
109  const TimeRange<Scalar>& appendRange,
110  const Ptr<InterpolationBufferBase<Scalar> > &interpBuffSink
111  )
112 {
113  TEUCHOS_ASSERT( !is_null(interpBuffSink) );
114 #ifdef HAVE_RYTHMOS_DEBUG
115  this->assertAppendPreconditions(interpBuffSource,appendRange,*interpBuffSink);
116 #endif // HAVE_RYTHMOS_DEBUG
117 
118  RCP<Teuchos::FancyOStream> out = this->getOStream();
119  Teuchos::OSTab ostab(out,1,"PointwiseInterpolationBufferAppender::append");
120  if ( Teuchos::as<int>(this->getVerbLevel()) >= Teuchos::as<int>(Teuchos::VERB_HIGH) ) {
121  *out << "Interpolation Buffer source range = [" << interpBuffSource.getTimeRange().lower() << "," <<
122  interpBuffSource.getTimeRange().upper() << "]" << std::endl;
123  *out << "Append range = [" << appendRange.lower() << "," << appendRange.upper() << "]" << std::endl;
124  *out << "Interpolation Buffer sink range = [" << interpBuffSink->getTimeRange().lower() << "," <<
125  interpBuffSink->getTimeRange().upper() << "]" << std::endl;
126  }
127  // Set up appendRange correctly to be either (] or [):
128  RCP<const TimeRange<Scalar> > correctedAppendRange = Teuchos::rcp(&appendRange,false);
129  if (compareTimeValues<Scalar>(interpBuffSink->getTimeRange().upper(),appendRange.lower()) == 0) {
130  // adding to end of buffer
131  correctedAppendRange = Teuchos::rcp(new TimeRange_oc<Scalar>(appendRange));
132  if ( Teuchos::as<int>(this->getVerbLevel()) >= Teuchos::as<int>(Teuchos::VERB_HIGH) ) {
133  *out << "Corrected append range = (" << correctedAppendRange->lower() << "," <<
134  correctedAppendRange->upper() << "]" << std::endl;
135  }
136  }
137  else if (compareTimeValues<Scalar>(interpBuffSink->getTimeRange().lower(),appendRange.upper()) == 0) {
138  // adding to beginning of buffer
139  correctedAppendRange = Teuchos::rcp(new TimeRange_co<Scalar>(appendRange));
140  if ( Teuchos::as<int>(this->getVerbLevel()) >= Teuchos::as<int>(Teuchos::VERB_HIGH) ) {
141  *out << "Corrected append range = [" << correctedAppendRange->lower() << "," <<
142  correctedAppendRange->upper() << ")" << std::endl;
143  }
144  }
145 
146  Array<Scalar> time_vec_in;
147  interpBuffSource.getNodes(&time_vec_in);
148 
149  Array<Scalar> time_vec;
150  selectPointsInTimeRange(time_vec_in,*correctedAppendRange,Teuchos::outArg(time_vec));
151  if ( Teuchos::as<int>(this->getVerbLevel()) >= Teuchos::as<int>(Teuchos::VERB_HIGH) ) {
152  *out << "Selected points for appending to sink buffer: " << time_vec << std::endl;
153  }
154 
155  Array<RCP<const Thyra::VectorBase<Scalar> > > x_vec;
156  Array<RCP<const Thyra::VectorBase<Scalar> > > xdot_vec;
157  Array<ScalarMag> accuracy_vec;
158  interpBuffSource.getPoints(time_vec, &x_vec, &xdot_vec, &accuracy_vec);
159 
160  if ( Teuchos::as<int>(this->getVerbLevel()) >= Teuchos::as<int>(Teuchos::VERB_HIGH) ) {
161  *out << "Sink buffer range before addPoints = [" << interpBuffSink->getTimeRange().lower() << "," <<
162  interpBuffSink->getTimeRange().upper() << "]" << std::endl;
163  }
164 
165  interpBuffSink->addPoints(time_vec, x_vec, xdot_vec);
166 
167  if ( Teuchos::as<int>(this->getVerbLevel()) >= Teuchos::as<int>(Teuchos::VERB_HIGH) ) {
168  *out << "Sink buffer range after addPoints = [" << interpBuffSink->getTimeRange().lower() << "," <<
169  interpBuffSink->getTimeRange().upper() << "]" << std::endl;
170  }
171 
172 }
173 
174 
175 template<class Scalar>
177  Teuchos::FancyOStream &out,
178  const Teuchos::EVerbosityLevel verbLevel
179  ) const
180 {
181  using Teuchos::as;
182  if (
183  (as<int>(verbLevel) == as<int>(Teuchos::VERB_DEFAULT))
184  || (as<int>(verbLevel) >= as<int>(Teuchos::VERB_LOW))
185  )
186  {
187  out << this->description() << std::endl;
188  }
189 }
190 
191 
192 template<class Scalar>
194  const RCP<ParameterList> &paramList
195  )
196 {
197  TEUCHOS_TEST_FOR_EXCEPT( is_null(paramList) );
198  paramList->validateParameters(*this->getValidParameters());
199  Teuchos::readVerboseObjectSublist(&*paramList,this);
200  setMyParamList(paramList);
201 }
202 
203 
204 template<class Scalar>
205 RCP<const ParameterList>
207 {
208  static RCP<Teuchos::ParameterList> validPL;
209  if (is_null(validPL)) {
210  RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
211  Teuchos::setupVerboseObjectSublist(&*pl);
212  validPL = pl;
213  }
214  return (validPL);
215 }
216 
217 
218 } // namespace Rythmos
219 
220 
221 #endif //RYTHMOS_POINTWISE_INTERPOLATION_BUFFER_APPENDER_HPP
Concrete InterplationBufferAppender subclass that just transfers notes without any regard for accurac...
virtual void getPoints(const Array< Scalar > &time_vec, Array< RCP< const Thyra::VectorBase< Scalar > > > *x_vec, Array< RCP< const Thyra::VectorBase< Scalar > > > *xdot_vec, Array< ScalarMag > *accuracy_vec) const =0
Get values from the buffer at different time points.
virtual void getNodes(Array< Scalar > *time_vec) const =0
Get interpolation nodes.
void append(const InterpolationBufferBase< Scalar > &interpBuffSource, const TimeRange< Scalar > &range, const Ptr< InterpolationBufferBase< Scalar > > &interpBuffSink)
Concrete implementation that simply copies the nodal points between the interpolation buffers...
RCP< PointwiseInterpolationBufferAppender< Scalar > > pointwiseInterpolationBufferAppender()
Nonmember constructor function.
Base class for an interpolation buffer.
Base class for strategy objects that append data from one InterplationBufferBase object to another...
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
virtual TimeRange< Scalar > getTimeRange() const =0
Return the range of time values where interpolation calls can be performed.