Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Stokhos_KokkosCrsMatrixMPVectorUnitTest.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Stokhos Package
5 // Copyright (2009) 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 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Eric T. Phipps (etphipp@sandia.gov).
38 //
39 // ***********************************************************************
40 // @HEADER
41 
44 
46 #include "Stokhos_Ensemble_Sizes.hpp"
49 
50 // For computing DeviceConfig
51 #include "Kokkos_Core.hpp"
52 
53 // Helper functions
54 template< typename IntType >
55 inline
56 IntType map_fem_graph_coord( const IntType& N,
57  const IntType& i,
58  const IntType& j,
59  const IntType& k )
60 {
61  return k + N * ( j + N * i );
62 }
63 
64 template < typename ordinal >
65 inline
66 ordinal generate_fem_graph( ordinal N,
67  std::vector< std::vector<ordinal> >& graph )
68 {
69  graph.resize( N * N * N, std::vector<ordinal>() );
70 
71  ordinal total = 0;
72 
73  for ( int i = 0; i < (int) N; ++i ) {
74  for ( int j = 0; j < (int) N; ++j ) {
75  for ( int k = 0; k < (int) N; ++k ) {
76 
77  const ordinal row = map_fem_graph_coord((int)N,i,j,k);
78 
79  graph[row].reserve(27);
80 
81  for ( int ii = -1; ii < 2; ++ii ) {
82  for ( int jj = -1; jj < 2; ++jj ) {
83  for ( int kk = -1; kk < 2; ++kk ) {
84  if ( 0 <= i + ii && i + ii < (int) N &&
85  0 <= j + jj && j + jj < (int) N &&
86  0 <= k + kk && k + kk < (int) N ) {
87  ordinal col = map_fem_graph_coord((int)N,i+ii,j+jj,k+kk);
88 
89  graph[row].push_back(col);
90  }
91  }}}
92  total += graph[row].size();
93  }}}
94 
95  return total;
96 }
97 
98 template <typename scalar, typename ordinal>
99 inline
100 scalar generate_matrix_coefficient( const ordinal nFEM,
101  const ordinal nStoch,
102  const ordinal iRowFEM,
103  const ordinal iColFEM,
104  const ordinal iStoch )
105 {
106  const scalar A_fem = ( 10.0 + scalar(iRowFEM) / scalar(nFEM) ) +
107  ( 5.0 + scalar(iColFEM) / scalar(nFEM) );
108 
109  const scalar A_stoch = ( 1.0 + scalar(iStoch) / scalar(nStoch) );
110 
111  return A_fem + A_stoch;
112  //return 1.0;
113 }
114 
115 template <typename scalar, typename ordinal>
116 inline
117 scalar generate_vector_coefficient( const ordinal nFEM,
118  const ordinal nStoch,
119  const ordinal iColFEM,
120  const ordinal iStoch )
121 {
122  const scalar X_fem = 100.0 + scalar(iColFEM) / scalar(nFEM);
123  const scalar X_stoch = 1.0 + scalar(iStoch) / scalar(nStoch);
124  return X_fem + X_stoch;
125  //return 1.0;
126 }
127 
128 // Reasonable tolerances for common precisions
129 template <typename Scalar> struct ScalarTol {};
130 template <> struct ScalarTol<float> { static float tol() { return 1e-4; } };
131 template <> struct ScalarTol<double> { static double tol() { return 1e-10; } };
132 
133 // Compare two rank-2 views for equality, to given precision
134 template <typename array_type, typename scalar_type>
135 bool compare_rank_2_views(const array_type& y,
136  const array_type& y_exp,
137  const scalar_type rel_tol,
138  const scalar_type abs_tol,
140 {
141  typedef typename array_type::size_type size_type;
142  typename array_type::HostMirror hy = Kokkos::create_mirror_view(y);
143  typename array_type::HostMirror hy_exp = Kokkos::create_mirror_view(y_exp);
144  Kokkos::deep_copy(hy, y);
145  Kokkos::deep_copy(hy_exp, y_exp);
146 
147  size_type num_rows = y.extent(0);
148  size_type num_cols = y.extent(1);
149  bool success = true;
150  for (size_type i=0; i<num_rows; ++i) {
151  for (size_type j=0; j<num_cols; ++j) {
152  scalar_type diff = std::abs( hy(i,j) - hy_exp(i,j) );
153  scalar_type tol = rel_tol*std::abs(hy_exp(i,j)) + abs_tol;
154  bool s = diff < tol;
155  out << "y_expected(" << i << "," << j << ") - "
156  << "y(" << i << "," << j << ") = " << hy_exp(i,j)
157  << " - " << hy(i,j) << " == "
158  << diff << " < " << tol << " : ";
159  if (s)
160  out << "passed";
161  else
162  out << "failed";
163  out << std::endl;
164  success = success && s;
165  }
166  }
167 
168  return success;
169 }
170 
171 // Helper function to build a diagonal matrix
172 template <typename MatrixType>
173 MatrixType
175  typename MatrixType::ordinal_type mp_vector_size) {
176  typedef typename MatrixType::ordinal_type ordinal_type;
177  typedef typename MatrixType::StaticCrsGraphType matrix_graph_type;
178  typedef typename MatrixType::values_type matrix_values_type;
179 
180  std::vector< std::vector<ordinal_type> > graph(nrow);
181  for (ordinal_type i=0; i<nrow; ++i)
182  graph[i] = std::vector<ordinal_type>(1, i);
183  ordinal_type graph_length = nrow;
184 
185  matrix_graph_type matrix_graph =
186  Kokkos::create_staticcrsgraph<matrix_graph_type>("graph", graph);
187  matrix_values_type matrix_values =
188  matrix_values_type("values", graph_length, mp_vector_size);
189 
190  MatrixType matrix("matrix", nrow, matrix_values, matrix_graph);
191  return matrix;
192 }
193 
194 //
195 // Tests
196 //
197 
198 // Kernel to set diagonal of a matrix to prescribed values
199 template <typename MatrixType>
202  typedef typename MatrixType::size_type size_type;
205 
206  const MatrixType m_matrix;
207  ReplaceDiagonalValuesKernel(const MatrixType matrix) : m_matrix(matrix) {};
208 
209  // Replace diagonal entry for row 'i' with a value
210  KOKKOS_INLINE_FUNCTION
211  void operator() (const size_type i) const {
212  const ordinal_type row = i;
213  const ordinal_type col = i;
214  value_type val = value_type(row);
215  m_matrix.replaceValues(row, &col, 1, &val, false, true);
216  }
217 
218  // Kernel launch
219  static void apply(const MatrixType matrix) {
220  const size_type nrow = matrix.numRows();
221  Kokkos::parallel_for( nrow, ReplaceDiagonalValuesKernel(matrix) );
222  }
223 
224  // Check the result is as expected
225  static bool check(const MatrixType matrix,
226  Teuchos::FancyOStream& out) {
227  typedef typename MatrixType::values_type matrix_values_type;
228  typename matrix_values_type::HostMirror host_matrix_values =
229  Kokkos::create_mirror_view(matrix.values);
230  Kokkos::deep_copy(host_matrix_values, matrix.values);
231  const ordinal_type nrow = matrix.numRows();
232  const ordinal_type vec_size = Kokkos::dimension_scalar(host_matrix_values);
233  bool success = true;
234  for (ordinal_type row=0; row<nrow; ++row) {
235  bool s = compareVecs(host_matrix_values(row),
236  "matrix_values(row)",
237  value_type(vec_size, row),
238  "value_type(row)",
239  0.0, 0.0, out);
240  success = success && s;
241  }
242  return success;
243  }
244 };
245 
246 // Kernel to add values to the diagonal of a matrix
247 template <typename MatrixType>
250  typedef typename MatrixType::size_type size_type;
253 
254  const MatrixType m_matrix;
255  AddDiagonalValuesKernel(const MatrixType matrix) : m_matrix(matrix) {};
256 
257  // Replace diagonal entry for row 'i' with a value
258  KOKKOS_INLINE_FUNCTION
259  void operator() (const size_type i) const {
260  const ordinal_type row = i;
261  const ordinal_type col = i;
262  value_type val = value_type(row);
263  m_matrix.sumIntoValues(row, &col, 1, &val, false, true);
264  }
265 
266  // Kernel launch
267  static void apply(const MatrixType matrix) {
268  const size_type nrow = matrix.numRows();
269  Kokkos::parallel_for( nrow, AddDiagonalValuesKernel(matrix) );
270  }
271 
272  // Check the result is as expected
273  static bool check(const MatrixType matrix,
274  Teuchos::FancyOStream& out) {
275  typedef typename MatrixType::values_type matrix_values_type;
276  typename matrix_values_type::HostMirror host_matrix_values =
277  Kokkos::create_mirror_view(matrix.values);
278  Kokkos::deep_copy(host_matrix_values, matrix.values);
279  const ordinal_type nrow = matrix.numRows();
280  const ordinal_type vec_size = Kokkos::dimension_scalar(host_matrix_values);
281  bool success = true;
282  for (ordinal_type row=0; row<nrow; ++row) {
283  bool s = compareVecs(host_matrix_values(row),
284  "matrix_values(row)",
285  value_type(vec_size, row),
286  "value_type(row)",
287  0.0, 0.0, out);
288  success = success && s;
289  }
290  return success;
291  }
292 };
293 
294 // Kernel to add values to the diagonal of a matrix where each thread
295 // adds to the same row (checks atomic really works)
296 template <typename MatrixType>
299  typedef typename MatrixType::size_type size_type;
302 
303  const MatrixType m_matrix;
304  AddDiagonalValuesAtomicKernel(const MatrixType matrix) : m_matrix(matrix) {};
305 
306  // Replace diagonal entry for row 'i' with a value
307  KOKKOS_INLINE_FUNCTION
308  void operator() (const size_type i) const {
309  const ordinal_type row = 0;
310  const ordinal_type col = 0;
312  m_matrix.sumIntoValues(row, &col, 1, &val, false, true);
313  }
314 
315  // Kernel launch
316  static void apply(const MatrixType matrix) {
317  const size_type nrow = matrix.numRows();
318  Kokkos::parallel_for( nrow, AddDiagonalValuesAtomicKernel(matrix) );
319  }
320 
321  // Check the result is as expected
322  static bool check(const MatrixType matrix,
323  Teuchos::FancyOStream& out) {
324  typedef typename MatrixType::values_type matrix_values_type;
325  typename matrix_values_type::HostMirror host_matrix_values =
326  Kokkos::create_mirror_view(matrix.values);
327  Kokkos::deep_copy(host_matrix_values, matrix.values);
328  const ordinal_type nrow = matrix.numRows();
329  const ordinal_type vec_size = Kokkos::dimension_scalar(host_matrix_values);
330  bool success = true;
331  for (ordinal_type row=0; row<nrow; ++row) {
332  value_type val;
333  if (row == 0)
334  val = value_type( vec_size, nrow*(nrow-1)/2 );
335  else
336  val = value_type( vec_size, 0.0 );
337  bool s = compareVecs(host_matrix_values(row),
338  "matrix_values(row)",
339  val,
340  "val",
341  0.0, 0.0, out);
342  success = success && s;
343  }
344  return success;
345  }
346 };
347 
348 const unsigned VectorSize = STOKHOS_DEFAULT_ENSEMBLE_SIZE;
349 
351  Kokkos_CrsMatrix_MP, ReplaceValues, MatrixScalar )
352 {
353  typedef typename MatrixScalar::ordinal_type Ordinal;
354  typedef typename MatrixScalar::execution_space Device;
355  typedef KokkosSparse::CrsMatrix<MatrixScalar,Ordinal,Device> Matrix;
356 
357  // Build diagonal matrix
358  Ordinal nrow = 10;
359  Matrix matrix = buildDiagonalMatrix<Matrix>(nrow, VectorSize);
360 
361  // Launch our kernel
363  kernel::apply(matrix);
364 
365  // Check the result
366  success = kernel::check(matrix, out);
367 }
368 
370  Kokkos_CrsMatrix_MP, SumIntoValues, MatrixScalar )
371 {
372  typedef typename MatrixScalar::ordinal_type Ordinal;
373  typedef typename MatrixScalar::execution_space Device;
374  typedef KokkosSparse::CrsMatrix<MatrixScalar,Ordinal,Device> Matrix;
375 
376  // Build diagonal matrix
377  Ordinal nrow = 10;
378  Matrix matrix = buildDiagonalMatrix<Matrix>(nrow, VectorSize);
379 
380  // Launch our kernel
381  typedef AddDiagonalValuesKernel<Matrix> kernel;
382  kernel::apply(matrix);
383 
384  // Check the result
385  success = kernel::check(matrix, out);
386 }
387 
389  Kokkos_CrsMatrix_MP, SumIntoValuesAtomic, MatrixScalar )
390 {
391  typedef typename MatrixScalar::ordinal_type Ordinal;
392  typedef typename MatrixScalar::execution_space Device;
393  typedef KokkosSparse::CrsMatrix<MatrixScalar,Ordinal,Device> Matrix;
394 
395  // Build diagonal matrix
396  Ordinal nrow = 10;
397  Matrix matrix = buildDiagonalMatrix<Matrix>(nrow, VectorSize);
398 
399  // Launch our kernel
401  kernel::apply(matrix);
402 
403  // Check the result
404  success = kernel::check(matrix, out);
405 }
406 
407 
408 // Some Stokhos unit tests use View types with a static rank
409 // In this case, dimensions should only be passed for the dynamic rank(s)
410 // These routines are specialized to select the appropriate View ctor
411 template < class ViewType, class OrdinalType, size_t I >
413  static ViewType create_view( const std::string & name, const OrdinalType & fem_length, const OrdinalType & stoch_length ) {
414  return ViewType(name, fem_length, stoch_length);
415  }
416 };
417 
418 template <class ViewType, class OrdinalType>
419 struct RankTypeSelector <ViewType,OrdinalType,1> {
420  static ViewType create_view( const std::string & name, const OrdinalType & fem_length, const OrdinalType & stoch_length ) {
421  (void) stoch_length; // unused if dyn_rank == 1; cast to void to silence compiler warnings
422  return ViewType(name, fem_length);
423  }
424 };
425 
426 template <class ViewType, class OrdinalType>
427 struct RankTypeSelector <ViewType,OrdinalType,0> {
428  static ViewType create_view( const std::string & name, const OrdinalType & fem_length, const OrdinalType & stoch_length ) {
429  (void) stoch_length; // unused if dyn_rank == 1; cast to void to silence compiler warnings
430  (void) fem_length; // unused if dyn_rank == 1; cast to void to silence compiler warnings
431  return ViewType(name);
432  }
433 };
434 
435 
436 template <typename VectorType, typename Multiply>
438  const typename VectorType::ordinal_type stoch_length,
439  KokkosSparse::DeviceConfig dev_config,
440  Multiply multiply_op,
442 {
443  typedef typename VectorType::ordinal_type ordinal_type;
444  typedef typename VectorType::value_type scalar_type;
445  typedef typename VectorType::storage_type storage_type;
447  typedef Kokkos::LayoutRight Layout;
448  typedef Kokkos::View< VectorType*, Layout, execution_space > block_vector_type;
449  typedef KokkosSparse::CrsMatrix< VectorType, ordinal_type, execution_space > block_matrix_type;
450  typedef typename block_matrix_type::StaticCrsGraphType matrix_graph_type;
451  typedef typename block_matrix_type::values_type matrix_values_type;
452 
453  // Check ensemble_length == storage_type::static_size for static storage
455  storage_type::is_static && storage_type::static_size != stoch_length,
456  std::logic_error,
457  "Static storage size must equal ensemble size");
458 
459  // Generate FEM graph:
460  ordinal_type fem_length = nGrid * nGrid * nGrid;
461  std::vector< std::vector<ordinal_type> > fem_graph;
462  ordinal_type fem_graph_length = generate_fem_graph( nGrid, fem_graph );
463 
464  //------------------------------
465  // Generate input multivector:
466 
467  // FIXME: Experimental view needs to be fixed so that construct is called
468  // when not initializing
469  // block_vector_type x =
470  // block_vector_type(Kokkos::ViewAllocateWithoutInitializing("x"), fem_length, stoch_length);
471  // block_vector_type y =
472  // block_vector_type(Kokkos::ViewAllocateWithoutInitializing("y"), fem_length, stoch_length);
473 
474  block_vector_type x =
475  block_vector_type("x", fem_length, stoch_length);
476  block_vector_type y =
477  block_vector_type("y", fem_length, stoch_length);
478 
479  typename block_vector_type::HostMirror hx = Kokkos::create_mirror_view( x );
480  typename block_vector_type::HostMirror hy = Kokkos::create_mirror_view( y );
481 
482  // View the block vector as an array of the embedded intrinsic type.
483  typename block_vector_type::HostMirror::array_type hax = hx ;
484  typename block_vector_type::HostMirror::array_type hay = hy ;
485 
486  for (ordinal_type iRowFEM=0; iRowFEM<fem_length; ++iRowFEM) {
487  for (ordinal_type iRowStoch=0; iRowStoch<stoch_length; ++iRowStoch) {
488  hax(iRowFEM,iRowStoch) =
489  generate_vector_coefficient<scalar_type>(
490  fem_length, stoch_length, iRowFEM, iRowStoch );
491  hay(iRowFEM,iRowStoch) = 0.0;
492  }
493  }
494 
495  Kokkos::deep_copy( x, hx );
496  Kokkos::deep_copy( y, hy );
497 
498  //------------------------------
499  // Generate block matrix
500 
501  matrix_graph_type matrix_graph =
502  Kokkos::create_staticcrsgraph<matrix_graph_type>(
503  std::string("test crs graph"), fem_graph);
504  // FIXME:
505  // matrix_values_type matrix_values =
506  // matrix_values_type(
507  // Kokkos::ViewAllocateWithoutInitializing("matrix"), fem_graph_length, stoch_length);
508  matrix_values_type matrix_values =
509  matrix_values_type("matrix", fem_graph_length, stoch_length);
510  block_matrix_type matrix(
511  "block_matrix", fem_length, matrix_values, matrix_graph);
512  matrix.dev_config = dev_config;
513 
514  typename matrix_values_type::HostMirror hM =
515  Kokkos::create_mirror_view( matrix.values );
516 
517  typename matrix_values_type::HostMirror::array_type haM = hM ;
518 
519  for (ordinal_type iRowFEM=0, iEntryFEM=0; iRowFEM<fem_length; ++iRowFEM) {
520  const ordinal_type row_size = fem_graph[iRowFEM].size();
521  for (ordinal_type iRowEntryFEM=0; iRowEntryFEM<row_size;
522  ++iRowEntryFEM, ++iEntryFEM) {
523  const ordinal_type iColFEM = fem_graph[iRowFEM][iRowEntryFEM];
524 
525  for (ordinal_type k=0; k<stoch_length; ++k) {
526  haM(iEntryFEM,k) =
527  generate_matrix_coefficient<scalar_type>(
528  fem_length, stoch_length, iRowFEM, iColFEM, k);
529  }
530  }
531  }
532 
533  Kokkos::deep_copy( matrix.values, hM );
534 
535  //------------------------------
536  // multiply
537 
538  multiply_op( matrix, x, y );
539 
540  //------------------------------
541  // generate correct answer
542 
543  typedef typename block_vector_type::array_type array_type;
544 #ifdef KOKKOS_ENABLE_DEPRECATED_CODE
545  array_type ay_expected =
546  array_type("ay_expected", fem_length, stoch_length);
547 #else
548  array_type ay_expected =
550 #endif
551  typename array_type::HostMirror hay_expected =
552  Kokkos::create_mirror_view(ay_expected);
553  for (ordinal_type iRowFEM=0, iEntryFEM=0; iRowFEM<fem_length; ++iRowFEM) {
554  const ordinal_type row_size = fem_graph[iRowFEM].size();
555  for (ordinal_type iRowEntryFEM=0; iRowEntryFEM<row_size;
556  ++iRowEntryFEM, ++iEntryFEM) {
557  const ordinal_type iColFEM = fem_graph[iRowFEM][iRowEntryFEM];
558  for (ordinal_type k=0; k<stoch_length; ++k) {
559  hay_expected(iRowFEM, k) +=
560  generate_matrix_coefficient<scalar_type>(
561  fem_length, stoch_length, iRowFEM, iColFEM, k) *
562  generate_vector_coefficient<scalar_type>(
563  fem_length, stoch_length, iColFEM, k );
564  }
565  }
566  }
567  Kokkos::deep_copy( ay_expected, hay_expected );
568 
569  //------------------------------
570  // check
571 
572  typename block_vector_type::array_type ay = y;
573  scalar_type rel_tol = ScalarTol<scalar_type>::tol();
574  scalar_type abs_tol = ScalarTol<scalar_type>::tol();
575  bool success = compare_rank_2_views(ay, ay_expected, rel_tol, abs_tol, out);
576 
577  return success;
578 }
579 
581  template <typename Matrix, typename InputVector, typename OutputVector>
582  void operator() (const Matrix& A,
583  const InputVector& x,
584  OutputVector& y) const {
585  KokkosSparse::spmv("N", typename Matrix::value_type(1.0) , A, x, typename Matrix::value_type(0.0), y);
586  }
587 };
588 
589 template <typename Tag>
591  Tag tag;
592  Stokhos_MV_Multiply_Op(const Tag& tg = Tag()) : tag(tg) {}
593 
594  template <typename Matrix, typename InputVector, typename OutputVector>
595  void operator() (const Matrix& A,
596  const InputVector& x,
597  OutputVector& y) const {
598  Stokhos::multiply(A, x, y, tag);
599  }
600 };
601 
604 
605 #define CRSMATRIX_MP_VECTOR_TESTS_MATRIXSCALAR( SCALAR ) \
606  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( \
607  Kokkos_CrsMatrix_MP, ReplaceValues, SCALAR ) \
608  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( \
609  Kokkos_CrsMatrix_MP, SumIntoValues, SCALAR ) \
610  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( \
611  Kokkos_CrsMatrix_MP, SumIntoValuesAtomic, SCALAR )
612 
613 #define CRSMATRIX_MP_VECTOR_TESTS_STORAGE( STORAGE ) \
614  typedef Sacado::MP::Vector<STORAGE> MP_Vector_ ## STORAGE; \
615  CRSMATRIX_MP_VECTOR_TESTS_MATRIXSCALAR( MP_Vector_ ## STORAGE )
616 
617 #define CRSMATRIX_MP_VECTOR_TESTS_ORDINAL_SCALAR_DEVICE( ORDINAL, SCALAR, DEVICE ) \
618  typedef Stokhos::StaticFixedStorage<ORDINAL,SCALAR,VectorSize,DEVICE> SFS; \
619  typedef Stokhos::DynamicStorage<ORDINAL,SCALAR,DEVICE> DS; \
620  CRSMATRIX_MP_VECTOR_TESTS_STORAGE( SFS ) \
621  CRSMATRIX_MP_VECTOR_TESTS_STORAGE( DS )
622 
623 #define CRSMATRIX_MP_VECTOR_TESTS_DEVICE( DEVICE ) \
624  CRSMATRIX_MP_VECTOR_TESTS_ORDINAL_SCALAR_DEVICE( int, double, DEVICE )
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
Stokhos::StandardStorage< int, double > storage_type
bool compareVecs(const VectorType1 &a1, const std::string &a1_name, const VectorType2 &a2, const std::string &a2_name, const ValueType &rel_tol, const ValueType &abs_tol, Teuchos::FancyOStream &out)
ordinal generate_fem_graph(ordinal N, std::vector< std::vector< ordinal > > &graph)
Definition: TestEpetra.cpp:77
static bool check(const MatrixType matrix, Teuchos::FancyOStream &out)
Stokhos_MV_Multiply_Op< Stokhos::DefaultMultiply > DefaultMultiply
Kokkos::DefaultExecutionSpace execution_space
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
MatrixType buildDiagonalMatrix(typename MatrixType::ordinal_type nrow, typename MatrixType::ordinal_type mp_vector_size)
scalar generate_vector_coefficient(const ordinal nFEM, const ordinal nStoch, const ordinal iColFEM, const ordinal iStoch)
KOKKOS_INLINE_FUNCTION constexpr std::enable_if< is_view_uq_pce< View< T, P...> >::value, unsigned >::type dimension_scalar(const View< T, P...> &view)
IntType map_fem_graph_coord(const IntType &N, const IntType &i, const IntType &j, const IntType &k)
Definition: TestEpetra.cpp:67
void multiply(const CrsMatrix< MatrixValue, Device, Layout > &A, const InputMultiVectorType &x, OutputMultiVectorType &y, const std::vector< OrdinalType > &col_indices, SingleColumnMultivectorMultiply)
static bool check(const MatrixType matrix, Teuchos::FancyOStream &out)
Kokkos_MV_Multiply_Op KokkosMultiply
void operator()(const Matrix &A, const InputVector &x, OutputVector &y) const
static void apply(const MatrixType matrix)
static bool check(const MatrixType matrix, Teuchos::FancyOStream &out)
void deep_copy(const Stokhos::CrsMatrix< ValueType, DstDevice, Layout > &dst, const Stokhos::CrsMatrix< ValueType, SrcDevice, Layout > &src)
KOKKOS_INLINE_FUNCTION PCE< Storage > abs(const PCE< Storage > &a)
static ViewType create_view(const std::string &name, const OrdinalType &fem_length, const OrdinalType &stoch_length)
static ViewType create_view(const std::string &name, const OrdinalType &fem_length, const OrdinalType &stoch_length)
static ViewType create_view(const std::string &name, const OrdinalType &fem_length, const OrdinalType &stoch_length)
scalar generate_matrix_coefficient(const ordinal nFEM, const ordinal nStoch, const ordinal iRowFEM, const ordinal iColFEM, const ordinal iStoch)
expr val()
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL(Kokkos_CrsMatrix_MP, ReplaceValues, MatrixScalar)
void operator()(const Matrix &A, const InputVector &x, OutputVector &y) const
bool compare_rank_2_views(const array_type &y, const array_type &y_exp, const scalar_type rel_tol, const scalar_type abs_tol, Teuchos::FancyOStream &out)
static void apply(const MatrixType matrix)
bool test_embedded_vector(const typename VectorType::ordinal_type nGrid, const typename VectorType::ordinal_type stoch_length, KokkosSparse::DeviceConfig dev_config, Multiply multiply_op, Teuchos::FancyOStream &out)
std::enable_if< Kokkos::is_view_uq_pce< Kokkos::View< InputType, InputP... > >::value &&Kokkos::is_view_uq_pce< Kokkos::View< OutputType, OutputP... > >::value >::type spmv(const char mode[], const AlphaType &a, const MatrixType &A, const Kokkos::View< InputType, InputP... > &x, const BetaType &b, const Kokkos::View< OutputType, OutputP... > &y, const RANK_ONE)
Stokhos::CrsMatrix< ValueType, Device, Layout >::HostMirror create_mirror_view(const Stokhos::CrsMatrix< ValueType, Device, Layout > &A)