Amesos2 - Direct Sparse Solver Interfaces  Version of the Day
Amesos2_MultiVecAdapter_decl.hpp
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // Amesos2: Templated Direct Sparse Solver Package
6 // Copyright 2011 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39 //
40 // ***********************************************************************
41 //
42 // @HEADER
43 
55 #ifndef AMESOS2_MULTIVEC_ADAPTER_DECL_HPP
56 #define AMESOS2_MULTIVEC_ADAPTER_DECL_HPP
57 
58 #include <Teuchos_RCP.hpp>
59 #include <Teuchos_Ptr.hpp>
60 #include <Teuchos_ArrayView.hpp>
61 #include <Tpetra_Map.hpp>
62 
63 #include "Amesos2_TypeDecl.hpp"
64 #include "Amesos2_VectorTraits.hpp"
65 
66 namespace Amesos2 {
67 
68 
175  template <class MV>
176  struct MultiVecAdapter {};
177 
185  template <class MV>
186  Teuchos::RCP<MultiVecAdapter<MV> >
187  createMultiVecAdapter(Teuchos::RCP<MV> mv){
188  using Teuchos::rcp;
189 
190  if(mv.is_null()) return Teuchos::null;
191  return( rcp(new MultiVecAdapter<MV>(mv)) );
192  }
193 
194  template <class MV>
195  Teuchos::RCP<const MultiVecAdapter<MV> >
196  createConstMultiVecAdapter(Teuchos::RCP<const MV> mv){
197  using Teuchos::rcp;
198  using Teuchos::rcp_const_cast;
199 
200  if(mv.is_null()) return Teuchos::null;
201  return( rcp(new MultiVecAdapter<MV>(Teuchos::rcp_const_cast<MV,const MV>(mv))).getConst() );
202  }
203 
204 
206  // Utilities for getting and putting data from MultiVecs //
208 
209  namespace Util {
210 
217  template <typename MV, typename V>
219 
220  typedef typename VectorTraits<V>::ptr_scalar_type ptr_return_type ;
221 
222  static ptr_return_type * get_pointer_to_vector ( const Teuchos::Ptr< MV> &mv ) ;
223 
224  static ptr_return_type * get_pointer_to_vector ( Teuchos::Ptr< MV> &mv ) ;
225 
226  static ptr_return_type * get_pointer_to_vector ( const Teuchos::Ptr< const MV > &mv ) ;
227 
228  static ptr_return_type * get_pointer_to_vector ( Teuchos::Ptr< const MV > &mv ) ;
229  };
230 
231  /*
232  * If the multivector scalar type and the desired scalar tpye are
233  * the same, then we can do a simple straight copy.
234  */
235  template <typename MV>
236  struct same_type_get_copy {
237  static void apply(const Teuchos::Ptr<const MV> mv,
238  const Teuchos::ArrayView<typename MV::scalar_t>& v,
239  const size_t ldx,
240  Teuchos::Ptr<const Tpetra::Map<typename MV::local_ordinal_t, typename MV::global_ordinal_t, typename MV::node_t> > distribution_map,
241  EDistribution distribution );
242  };
243 
244  /*
245  * In the case where the scalar type of the multi-vector and the
246  * corresponding S type are different, then we need to first get a
247  * copy of the scalar values, then convert each one into the S
248  * type before inserting into the vals array.
249  */
250  template <typename MV, typename S>
251  struct diff_type_get_copy {
252  static void apply(const Teuchos::Ptr<const MV> mv,
253  const Teuchos::ArrayView<S>& v,
254  const size_t ldx,
255  Teuchos::Ptr<const Tpetra::Map<typename MV::local_ordinal_t, typename MV::global_ordinal_t, typename MV::node_t> > distribution_map,
256  EDistribution distribution );
257  };
258 
265  template <class MV, typename S>
267  static void
268  do_get (const Teuchos::Ptr<const MV>& mv,
269  const Teuchos::ArrayView<S>& vals,
270  const size_t ldx,
271  Teuchos::Ptr<const Tpetra::Map<typename MV::local_ordinal_t, typename MV::global_ordinal_t, typename MV::node_t> > distribution_map,
272  EDistribution distribution = ROOTED);
273 
274  static void
275  do_get (const Teuchos::Ptr<const MV>& mv,
276  const Teuchos::ArrayView<S>& vals,
277  const size_t ldx,
278  EDistribution distribution,
279  typename MV::global_ordinal_t indexBase = 0);
280 
281  static void
282  do_get (const Teuchos::Ptr<const MV>& mv,
283  const Teuchos::ArrayView<S>& vals,
284  const size_t ldx);
285  };
286 
287  template <class MV, typename KV>
288  struct get_1d_copy_helper_kokkos_view {
289  static void
290  do_get (const Teuchos::Ptr<const MV>& mv,
291  KV& kokkos_vals,
292  const size_t ldx,
293  Teuchos::Ptr<const Tpetra::Map<typename MV::local_ordinal_t, typename MV::global_ordinal_t, typename MV::node_t> > distribution_map,
294  EDistribution distribution = ROOTED);
295 
296  static void
297  do_get (const Teuchos::Ptr<const MV>& mv,
298  KV& kokkos_vals,
299  const size_t ldx,
300  EDistribution distribution,
301  typename MV::global_ordinal_t indexBase = 0);
302 
303  static void
304  do_get (const Teuchos::Ptr<const MV>& mv,
305  KV& kokkos_vals,
306  const size_t ldx);
307  };
308 
309  /*
310  * If the multivector scalar type and the desired scalar tpye are
311  * the same, then we can do a simple straight copy.
312  */
313  template <typename MV>
314  struct same_type_data_put {
315  static void apply(const Teuchos::Ptr<MV>& mv,
316  const Teuchos::ArrayView<typename MV::scalar_t>& data,
317  const size_t ldx,
318  Teuchos::Ptr<const Tpetra::Map<typename MV::local_ordinal_t, typename MV::global_ordinal_t, typename MV::node_t> > distribution_map,
319  EDistribution distribution );
320  };
321 
322  /*
323  * In the case where the scalar type of the multi-vector and the
324  * corresponding S type are different, then we need to first get a
325  * copy of the scalar values, then convert each one into the S
326  * type before inserting into the vals array.
327  */
328  template <typename MV, typename S>
329  struct diff_type_data_put {
330  static void apply(const Teuchos::Ptr<MV>& mv,
331  const Teuchos::ArrayView<S>& data,
332  const size_t ldx,
333  Teuchos::Ptr<const Tpetra::Map<typename MV::local_ordinal_t, typename MV::global_ordinal_t, typename MV::node_t> > distribution_map,
334  EDistribution distribution );
335  };
336 
343  template <class MV, typename S>
345  static void do_put(const Teuchos::Ptr<MV>& mv,
346  const Teuchos::ArrayView<S>& data,
347  const size_t ldx,
348  Teuchos::Ptr<const Tpetra::Map<typename MV::local_ordinal_t, typename MV::global_ordinal_t, typename MV::node_t> > distribution_map,
349  EDistribution distribution = ROOTED);
350 
351  static void do_put(const Teuchos::Ptr<MV>& mv,
352  const Teuchos::ArrayView<S>& data,
353  const size_t ldx,
354  EDistribution distribution, typename MV::global_ordinal_t indexBase = 0);
355 
356  static void do_put(const Teuchos::Ptr<MV>& mv,
357  const Teuchos::ArrayView<S>& data,
358  const size_t ldx);
359  };
360 
361  template <class MV, typename KV>
362  struct put_1d_data_helper_kokkos_view {
363  static void do_put(const Teuchos::Ptr<MV>& mv,
364  KV& kokkos_data,
365  const size_t ldx,
366  Teuchos::Ptr<const Tpetra::Map<typename MV::local_ordinal_t, typename MV::global_ordinal_t, typename MV::node_t> > distribution_map,
367  EDistribution distribution = ROOTED);
368 
369  static void do_put(const Teuchos::Ptr<MV>& mv,
370  KV& kokkos_data,
371  const size_t ldx,
372  EDistribution distribution, typename MV::global_ordinal_t indexBase = 0);
373 
374  static void do_put(const Teuchos::Ptr<MV>& mv,
375  KV& kokkos_data,
376  const size_t ldx);
377  };
378  }
379 } // end namespace Amesos2
380 
383 #ifdef HAVE_AMESOS2_EPETRA
385 #endif
386 
387 #endif // AMESOS2_MULTIVEC_ADAPTER_DECL_HPP
Teuchos::RCP< MultiVecAdapter< MV > > createMultiVecAdapter(Teuchos::RCP< MV > mv)
Factory creation method for MultiVecAdapters.
Definition: Amesos2_MultiVecAdapter_decl.hpp:187
Amesos2::MultiVecAdapter specialization for the Kokkos::View class.
Helper class for getting 1-D copies of multivectors.
Definition: Amesos2_MultiVecAdapter_decl.hpp:266
Helper struct for getting pointers to the MV data - only used when number of vectors = 1 and single M...
Definition: Amesos2_MultiVecAdapter_decl.hpp:218
static void do_get(const Teuchos::Ptr< const MV > &mv, const Teuchos::ArrayView< S > &vals, const size_t ldx, Teuchos::Ptr< const Tpetra::Map< typename MV::local_ordinal_t, typename MV::global_ordinal_t, typename MV::node_t > > distribution_map, EDistribution distribution=ROOTED)
Helper class for getting 1-D copies of multivectors.
Definition: Amesos2_MultiVecAdapter_def.hpp:146
Amesos2::MultiVecAdapter specialization for the Tpetra::MultiVector class.
Amesos2::MultiVecAdapter specialization for the Epetra_MultiVector class.
Enum and other types declarations for Amesos2.
static void do_put(const Teuchos::Ptr< MV > &mv, const Teuchos::ArrayView< S > &data, const size_t ldx, Teuchos::Ptr< const Tpetra::Map< typename MV::local_ordinal_t, typename MV::global_ordinal_t, typename MV::node_t > > distribution_map, EDistribution distribution=ROOTED)
Helper class for putting 1-D data arrays into multivectors.
Definition: Amesos2_MultiVecAdapter_def.hpp:317
EDistribution
Definition: Amesos2_TypeDecl.hpp:123
Helper class for putting 1-D data arrays into multivectors.
Definition: Amesos2_MultiVecAdapter_decl.hpp:344
A templated MultiVector class adapter for Amesos2.
Definition: Amesos2_MultiVecAdapter_decl.hpp:176