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
93 typedef User userCoord_t;
95 typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
96 typedef Xpetra::TpetraMultiVector<
116 std::vector<const scalar_t *> &
weights, std::vector<int> &weightStrides);
134 ids = map_->getLocalElementList().getRawPtr();
138 Kokkos::View<const gno_t *, typename node_t::device_type> &ids)
const {
139 if (map_->lib() == Xpetra::UseTpetra) {
140 using device_type =
typename node_t::device_type;
141 const xt_mvector_t *tvector =
142 dynamic_cast<const xt_mvector_t *
>(vector_.get());
150 ids = Kokkos::create_mirror_view_and_copy(device_type(),
151 tvector->getTpetra_MultiVector()->getMap()->getMyGlobalIndices());
153 else if (map_->lib() == Xpetra::UseEpetra) {
154 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
158 throw std::logic_error(
"Epetra requested, but Trilinos is not "
159 "built with Epetra");
163 throw std::logic_error(
"getIDsKokkosView called but not on Tpetra or Epetra!");
171 if(idx<0 || idx >= numWeights_)
173 std::ostringstream emsg;
174 emsg << __FILE__ <<
":" << __LINE__
175 <<
" Invalid weight index " << idx << std::endl;
176 throw std::runtime_error(emsg.str());
180 weights_[idx].getStridedList(length, weights, stride);
184 typename node_t::device_type> &wgt)
const {
185 typedef Kokkos::View<scalar_t**, typename node_t::device_type> view_t;
186 wgt = view_t(
"wgts", vector_->getLocalLength(), numWeights_);
187 typename view_t::HostMirror host_wgt = Kokkos::create_mirror_view(wgt);
188 for(
int idx = 0; idx < numWeights_; ++idx) {
192 weights_[idx].getStridedList(length, weights, stride);
193 size_t fill_index = 0;
194 for(
size_t n = 0; n < length; n += stride) {
195 host_wgt(fill_index++,idx) = weights[n];
198 Kokkos::deep_copy(wgt, host_wgt);
211 Kokkos::View<
scalar_t **, Kokkos::LayoutLeft,
212 typename node_t::device_type> & elements)
const;
214 template <
typename Adapter>
218 template <
typename Adapter>
224 RCP<const User> invector_;
225 RCP<const x_mvector_t> vector_;
226 RCP<const Xpetra::Map<lno_t, gno_t, node_t> > map_;
229 ArrayRCP<StridedData<lno_t, scalar_t> > weights_;
236 template <
typename User>
238 const RCP<const User> &invector,
239 std::vector<const scalar_t *> &
weights, std::vector<int> &weightStrides):
240 invector_(invector), vector_(), map_(),
241 numWeights_(weights.size()), weights_(weights.size())
246 RCP<x_mvector_t> tmp =
248 vector_ = rcp_const_cast<
const x_mvector_t>(tmp);
252 map_ = vector_->getMap();
254 size_t length = vector_->getLocalLength();
256 if (length > 0 && numWeights_ > 0){
258 for (
int w=0; w < numWeights_; w++){
259 if (weightStrides.size())
260 stride = weightStrides[w];
261 ArrayRCP<const scalar_t> wgtV(weights[w], 0, stride*length,
false);
262 weights_[w] = input_t(wgtV, stride);
269 template <
typename User>
271 const RCP<const User> &invector):
272 invector_(invector), vector_(), map_(),
273 numWeights_(0), weights_()
276 RCP<x_mvector_t> tmp =
278 vector_ = rcp_const_cast<
const x_mvector_t>(tmp);
282 map_ = vector_->getMap();
286 template <
typename User>
288 const scalar_t *&elements,
int &stride,
int idx)
const
293 if (map_->lib() == Xpetra::UseTpetra){
294 const xt_mvector_t *tvector =
295 dynamic_cast<const xt_mvector_t *
>(vector_.get());
297 vecsize = tvector->getLocalLength();
299 ArrayRCP<const scalar_t> data = tvector->getData(idx);
300 elements = data.get();
303 else if (map_->lib() == Xpetra::UseEpetra){
304 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
305 typedef Xpetra::EpetraMultiVectorT<gno_t,node_t> xe_mvector_t;
306 const xe_mvector_t *evector =
307 dynamic_cast<const xe_mvector_t *
>(vector_.get());
309 vecsize = evector->getLocalLength();
311 ArrayRCP<const double> data = evector->getData(idx);
315 elements =
reinterpret_cast<const scalar_t *
>(data.get());
318 throw std::logic_error(
"Epetra requested, but Trilinos is not "
319 "built with Epetra");
323 throw std::logic_error(
"invalid underlying lib");
328 template <
typename User>
331 Kokkos::View<scalar_t **, Kokkos::LayoutLeft, typename node_t::device_type> & elements)
const
333 if (map_->lib() == Xpetra::UseTpetra){
334 const xt_mvector_t *tvector =
335 dynamic_cast<const xt_mvector_t *
>(vector_.get());
337 Kokkos::View<scalar_t **, Kokkos::LayoutLeft, typename node_t::device_type> view2d =
338 tvector->getTpetra_MultiVector()->template getLocalView<typename node_t::device_type>(Tpetra::Access::ReadWrite);
344 else if (map_->lib() == Xpetra::UseEpetra){
345 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
346 typedef Xpetra::EpetraMultiVectorT<gno_t,node_t> xe_mvector_t;
347 const xe_mvector_t *evector =
348 dynamic_cast<const xe_mvector_t *
>(vector_.get());
350 Kokkos::View<scalar_t **, Kokkos::LayoutLeft, typename node_t::device_type>
351 (
"elements", evector->getLocalLength(), evector->getNumVectors());
352 if(evector->getLocalLength() > 0) {
353 for(
size_t idx = 0; idx < evector->getNumVectors(); ++idx) {
354 const scalar_t * ptr;
356 getEntriesView(ptr, stride, idx);
357 for(
size_t n = 0; n < evector->getLocalLength(); ++n) {
358 elements(n, idx) = ptr[n];
363 throw std::logic_error(
"Epetra requested, but Trilinos is not "
364 "built with Epetra");
368 throw std::logic_error(
"getEntriesKokkosView called but not using Tpetra or Epetra!");
373 template <
typename User>
374 template <
typename Adapter>
376 const User &in, User *&out,
381 ArrayRCP<gno_t> importList;
385 (solution,
this, importList);
391 importList.getRawPtr());
397 template <
typename User>
398 template <
typename Adapter>
400 const User &in, RCP<User> &out,
405 ArrayRCP<gno_t> importList;
409 (solution,
this, importList);
415 importList.getRawPtr());
Helper functions for Partitioning Problems.
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
map_t::global_ordinal_type 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.
typename InputTraits< User >::part_t part_t
void getIDsKokkosView(Kokkos::View< const gno_t *, typename node_t::device_type > &ids) const
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...
virtual void getIDsKokkosView(ConstIdsDeviceView &ids) const
Provide a Kokkos view to this process' identifiers.
int getNumEntriesPerID() const
Return the number of vectors.
typename InputTraits< User >::node_t node_t
A PartitioningSolution is a solution to a partitioning problem.
typename InputTraits< User >::gno_t gno_t
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.
void getEntriesView(const scalar_t *&elements, int &stride, int idx=0) const
Provide a pointer to the elements of the specified vector.
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
void getWeightsKokkos2dView(Kokkos::View< scalar_t **, typename node_t::device_type > &wgt) const
size_t getLocalNumIDs() const
Returns the number of objects on this process.
typename BaseAdapter< User >::scalar_t scalar_t
typename InputTraits< User >::lno_t lno_t
void getWeightsView(const scalar_t *&weights, int &stride, int idx) const
Provide pointer to a weight array with stride.
void getEntriesKokkosView(Kokkos::View< scalar_t **, Kokkos::LayoutLeft, typename node_t::device_type > &elements) const
This file defines the StridedData class.
Zoltan2::BasicUserTypes< zscalar_t, zlno_t, zgno_t > user_t