56 const Array<Scalar> &t_values
57 ,
typename DataStore<Scalar>::DataStoreVector_t *data_out )
const
62 typedef Teuchos::ScalarTraits<Scalar> ST;
63#ifdef HAVE_RYTHMOS_DEBUG
64 assertInterpolatePreconditions((*nodes_),t_values,data_out);
66 RCP<Teuchos::FancyOStream> out = this->getOStream();
67 Teuchos::OSTab ostab(out,1,
"HI::interpolator");
68 if ( Teuchos::as<int>(this->getVerbLevel()) >= Teuchos::as<int>(Teuchos::VERB_HIGH) ) {
69 *out <<
"(*nodes_):" << std::endl;
70 for (Teuchos::Ordinal i=0 ; i<(*nodes_).size() ; ++i) {
71 *out <<
"(*nodes_)[" << i <<
"] = " << std::endl;
72 (*nodes_)[i].describe(*out,Teuchos::VERB_EXTREME);
74 *out <<
"t_values = " << std::endl;
75 for (Teuchos::Ordinal i=0 ; i<t_values.size() ; ++i) {
76 *out <<
"t_values[" << i <<
"] = " << t_values[i] << std::endl;
78 for (Teuchos::Ordinal i=0; i<data_out->size() ; ++i) {
79 *out <<
"data_out[" << i <<
"] = " << std::endl;
80 (*data_out)[i].describe(*out,Teuchos::VERB_EXTREME);
84 if (t_values.size() == 0) {
88 if ((*nodes_).size() == 1) {
91 DataStore<Scalar> DS((*nodes_)[0]);
92 data_out->push_back(DS);
96 for (
int i=0 ; i<Teuchos::as<int>((*nodes_).size())-1 ; ++i) {
97 const Scalar& t0 = (*nodes_)[i].time;
98 const Scalar& t1 = (*nodes_)[i+1].time;
99 while ((t0 <= t_values[n]) && (t_values[n] <= t1)) {
100 const Scalar& t = t_values[n];
103 DataStore<Scalar> DS((*nodes_)[i]);
104 data_out->push_back(DS);
105 }
else if (t == t1) {
106 DataStore<Scalar> DS((*nodes_)[i+1]);
107 data_out->push_back(DS);
109 RCP<const Thyra::VectorBase<Scalar> > x0 = (*nodes_)[i ].x;
110 RCP<const Thyra::VectorBase<Scalar> > x1 = (*nodes_)[i+1].x;
111 RCP<const Thyra::VectorBase<Scalar> > xdot0 = (*nodes_)[i ].xdot;
112 RCP<const Thyra::VectorBase<Scalar> > xdot1 = (*nodes_)[i+1].xdot;
115 RCP<Thyra::VectorBase<Scalar> > tmp_vec = x0->clone_v();
116 RCP<Thyra::VectorBase<Scalar> > xdot_temp = x1->clone_v();
119 Scalar t_t0 = t - t0;
120 Scalar t_t1 = t - t1;
124 Thyra::Vt_S(xdot_temp.ptr(),Scalar(ST::one()/dt));
125 Thyra::Vp_StV(xdot_temp.ptr(),Scalar(-ST::one()/dt),*x0);
128 DataStore<Scalar> DS;
133 RCP<Thyra::VectorBase<Scalar> > x_vec = x0->clone_v();
134 Thyra::Vp_StV(x_vec.ptr(),t_t0,*xdot0);
135 tmp_t = t_t0*t_t0/dt;
136 Thyra::V_StVpStV(tmp_vec.ptr(),tmp_t,*xdot_temp,Scalar(-ST::one()*tmp_t),*xdot0);
137 Thyra::Vp_V(x_vec.ptr(),*tmp_vec);
138 tmp_t = t_t0*t_t0*t_t1/dt2;
139 Thyra::V_StVpStV(tmp_vec.ptr(),tmp_t,*xdot1,Scalar(-2*tmp_t),*xdot_temp);
140 Thyra::Vp_StV(tmp_vec.ptr(),tmp_t,*xdot0);
141 Thyra::Vp_V(x_vec.ptr(),*tmp_vec);
146 RCP<Thyra::VectorBase<Scalar> > xdot_vec = xdot0->clone_v();
148 Thyra::Vp_StV(xdot_vec.ptr(),Scalar(2*tmp_t),*xdot_temp);
149 Thyra::Vp_StV(xdot_vec.ptr(),Scalar(-2*tmp_t),*xdot0);
150 tmp_t = Scalar((2*t_t0*t_t1+t_t0*t_t0)/dt2);
151 Thyra::V_StVpStV(tmp_vec.ptr(),tmp_t,*xdot1,Scalar(-2*tmp_t),*xdot_temp);
152 Thyra::Vp_StV(tmp_vec.ptr(),tmp_t,*xdot0);
153 Thyra::Vp_V(xdot_vec.ptr(),*tmp_vec);
158 DS.accuracy = (t_t0)*(t_t0)*(t_t1)*(t_t1);
161 data_out->push_back(DS);
164 if (n == Teuchos::as<int>(t_values.size())) {