50 #ifndef _ZOLTAN2_XPETRAMULTIVECTORADAPTER_HPP_ 51 #define _ZOLTAN2_XPETRAMULTIVECTORADAPTER_HPP_ 58 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA) 59 #include <Xpetra_EpetraMultiVector.hpp> 61 #include <Xpetra_TpetraMultiVector.hpp> 82 template <
typename User>
86 #ifndef DOXYGEN_SHOULD_SKIP_THIS 94 typedef User userCoord_t;
96 typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
97 typedef Xpetra::TpetraMultiVector<
121 std::vector<const scalar_t *> &
weights, std::vector<int> &weightStrides);
139 ids = map_->getNodeElementList().getRawPtr();
146 env_->localInputAssertion(__FILE__, __LINE__,
"invalid weight index",
149 weights_[idx].getStridedList(length, weights, stride);
158 void getEntriesView(
const scalar_t *&elements,
int &stride,
int idx=0)
const;
160 template <
typename Adapter>
164 template <
typename Adapter>
170 RCP<const User> invector_;
171 RCP<const x_mvector_t> vector_;
172 RCP<const Xpetra::Map<lno_t, gno_t, node_t> > map_;
173 RCP<Environment> env_;
177 ArrayRCP<StridedData<lno_t, scalar_t> > weights_;
184 template <
typename User>
186 const RCP<const User> &invector,
187 std::vector<const scalar_t *> &
weights, std::vector<int> &weightStrides):
188 invector_(invector), vector_(), map_(),
190 numWeights_(weights.size()), weights_(weights.size())
194 RCP<x_mvector_t> tmp =
196 vector_ = rcp_const_cast<
const x_mvector_t>(tmp);
197 map_ = vector_->getMap();
198 base_ = map_->getIndexBase();
200 size_t length = vector_->getLocalLength();
202 if (length > 0 && numWeights_ > 0){
204 for (
int w=0; w < numWeights_; w++){
205 if (weightStrides.size())
206 stride = weightStrides[w];
207 ArrayRCP<const scalar_t> wgtV(weights[w], 0, stride*length,
false);
208 weights_[w] = input_t(wgtV, stride);
215 template <
typename User>
217 const RCP<const User> &invector):
218 invector_(invector), vector_(), map_(),
220 numWeights_(0), weights_()
222 RCP<x_mvector_t> tmp =
224 vector_ = rcp_const_cast<
const x_mvector_t>(tmp);
225 map_ = vector_->getMap();
226 base_ = map_->getIndexBase();
230 template <
typename User>
232 const scalar_t *&elements,
int &stride,
int idx)
const 237 if (map_->lib() == Xpetra::UseTpetra){
238 const xt_mvector_t *tvector =
239 dynamic_cast<const xt_mvector_t *
>(vector_.get());
241 vecsize = tvector->getLocalLength();
243 ArrayRCP<const scalar_t> data = tvector->getData(idx);
244 elements = data.get();
247 else if (map_->lib() == Xpetra::UseEpetra){
248 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA) 249 typedef Xpetra::EpetraMultiVectorT<gno_t,node_t> xe_mvector_t;
250 const xe_mvector_t *evector =
251 dynamic_cast<const xe_mvector_t *
>(vector_.get());
253 vecsize = evector->getLocalLength();
255 ArrayRCP<const double> data = evector->getData(idx);
259 elements =
reinterpret_cast<const scalar_t *
>(data.get());
262 throw std::logic_error(
"Epetra requested, but Trilinos is not " 263 "built with Epetra");
267 throw std::logic_error(
"invalid underlying lib");
272 template <
typename User>
273 template <
typename Adapter>
275 const User &in, User *&out,
280 ArrayRCP<gno_t> importList;
284 (solution,
this, importList);
290 importList.getRawPtr());
296 template <
typename User>
297 template <
typename Adapter>
299 const User &in, RCP<User> &out,
304 ArrayRCP<gno_t> importList;
308 (solution,
this, importList);
314 importList.getRawPtr());
InputTraits< User >::scalar_t scalar_t
Helper functions for Partitioning Problems.
fast typical checks for valid arguments
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
InputTraits< User >::gno_t gno_t
XpetraMultiVectorAdapter(const RCP< const User > &invector, std::vector< const scalar_t * > &weights, std::vector< int > &weightStrides)
Constructor.
Defines the VectorAdapter interface.
static RCP< User > doMigration(const User &from, size_t numLocalRows, const gno_t *myNewRows)
Migrate the object Given a user object and a new row distribution, create and return a new user objec...
Traits of Xpetra classes, including migration method.
static RCP< User > convertToXpetra(const RCP< User > &a)
Convert the object to its Xpetra wrapped version.
size_t getImportList(const PartitioningSolution< SolutionAdapter > &solution, const DataAdapter *const data, ArrayRCP< typename DataAdapter::gno_t > &imports)
From a PartitioningSolution, get a list of IDs to be imported. Assumes part numbers in PartitioningSo...
int getNumEntriesPerID() const
Return the number of vectors (typically one).
A PartitioningSolution is a solution to a partitioning problem.
void applyPartitioningSolution(const User &in, User *&out, const PartitioningSolution< Adapter > &solution) const
VectorAdapter defines the interface for vector input.
The StridedData class manages lists of weights or coordinates.
~XpetraMultiVectorAdapter()
Destructor.
void getEntriesView(const scalar_t *&elements, int &stride, int idx=0) const
Provide a pointer to the elements of the specified vector.
InputTraits< User >::part_t part_t
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
int getNumWeightsPerID() const
Returns the number of weights per object. Number of weights per object should be zero or greater...
An adapter for Xpetra::MultiVector.
void getIDsView(const gno_t *&ids) const
size_t getLocalNumIDs() const
Returns the number of objects on this process.
void getWeightsView(const scalar_t *&weights, int &stride, int idx) const
Provide pointer to a weight array with stride.
This file defines the StridedData class.