Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
Stokhos_KokkosCrsMatrixUQPCEUnitTest.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
42#include "Teuchos_UnitTestHelpers.hpp"
44
46#include "KokkosSparse_CrsMatrix.hpp"
47#include "KokkosSparse_spmv.hpp"
53
54// For computing DeviceConfig
55#include "Kokkos_Core.hpp"
56
57// Helper functions
58template< typename IntType >
59inline
60IntType map_fem_graph_coord( const IntType& N,
61 const IntType& i,
62 const IntType& j,
63 const IntType& k )
64{
65 return k + N * ( j + N * i );
66}
67
68template < typename ordinal >
69inline
70ordinal generate_fem_graph( ordinal N,
71 std::vector< std::vector<ordinal> >& graph )
72{
73 graph.resize( N * N * N, std::vector<ordinal>() );
74
75 ordinal total = 0;
76
77 for ( int i = 0; i < (int) N; ++i ) {
78 for ( int j = 0; j < (int) N; ++j ) {
79 for ( int k = 0; k < (int) N; ++k ) {
80
81 const ordinal row = map_fem_graph_coord((int)N,i,j,k);
82
83 graph[row].reserve(27);
84
85 for ( int ii = -1; ii < 2; ++ii ) {
86 for ( int jj = -1; jj < 2; ++jj ) {
87 for ( int kk = -1; kk < 2; ++kk ) {
88 if ( 0 <= i + ii && i + ii < (int) N &&
89 0 <= j + jj && j + jj < (int) N &&
90 0 <= k + kk && k + kk < (int) N ) {
91 ordinal col = map_fem_graph_coord((int)N,i+ii,j+jj,k+kk);
92
93 graph[row].push_back(col);
94 }
95 }}}
96 total += graph[row].size();
97 }}}
98
99 return total;
100}
101
102template <typename scalar, typename ordinal>
103inline
104scalar generate_matrix_coefficient( const ordinal nFEM,
105 const ordinal nStoch,
106 const ordinal iRowFEM,
107 const ordinal iColFEM,
108 const ordinal iStoch )
109{
110 const scalar A_fem = ( 10.0 + scalar(iRowFEM) / scalar(nFEM) ) +
111 ( 5.0 + scalar(iColFEM) / scalar(nFEM) );
112
113 const scalar A_stoch = ( 1.0 + scalar(iStoch) / scalar(nStoch) );
114
115 return A_fem + A_stoch;
116 //return 1.0;
117}
118
119template <typename scalar, typename ordinal>
120inline
121scalar generate_vector_coefficient( const ordinal nFEM,
122 const ordinal nStoch,
123 const ordinal iColFEM,
124 const ordinal iStoch )
125{
126 const scalar X_fem = 100.0 + scalar(iColFEM) / scalar(nFEM);
127 const scalar X_stoch = 1.0 + scalar(iStoch) / scalar(nStoch);
128 return X_fem + X_stoch;
129 //return 1.0;
130}
131
132template <typename kokkos_cijk_type, typename ordinal_type>
133kokkos_cijk_type build_cijk(ordinal_type stoch_dim,
134 ordinal_type poly_ord)
135{
136 using Teuchos::RCP;
137 using Teuchos::rcp;
138 using Teuchos::Array;
139
140 typedef typename kokkos_cijk_type::value_type value_type;
141 typedef typename kokkos_cijk_type::execution_space execution_space;
146
147 // Create product basis
148 Array< RCP<const one_d_basis> > bases(stoch_dim);
149 for (ordinal_type i=0; i<stoch_dim; i++)
150 bases[i] = rcp(new legendre_basis(poly_ord, true));
151 RCP<const product_basis> basis = rcp(new product_basis(bases));
152
153 // Triple product tensor
154 RCP<Cijk> cijk = basis->computeTripleProductTensor();
155
156 // Kokkos triple product tensor
157 kokkos_cijk_type kokkos_cijk =
158 Stokhos::create_product_tensor<execution_space>(*basis, *cijk);
159
160 return kokkos_cijk;
161}
162
163// Reasonable tolerances for common precisions
164template <typename Scalar> struct ScalarTol {};
165template <> struct ScalarTol<float> { static float tol() { return 1e-4; } };
166template <> struct ScalarTol<double> { static double tol() { return 1e-10; } };
167
168// Compare two rank-2 views for equality, to given precision
169template <typename array_type, typename scalar_type>
170bool compare_rank_2_views(const array_type& y,
171 const array_type& y_exp,
172 const scalar_type rel_tol,
173 const scalar_type abs_tol,
174 Teuchos::FancyOStream& out)
175{
176 typedef typename array_type::size_type size_type;
177 typename array_type::HostMirror hy = Kokkos::create_mirror_view(y);
178 typename array_type::HostMirror hy_exp = Kokkos::create_mirror_view(y_exp);
179 Kokkos::deep_copy(hy, y);
180 Kokkos::deep_copy(hy_exp, y_exp);
181
182 size_type num_rows = y.extent(0);
183 size_type num_cols = y.extent(1);
184 bool success = true;
185 for (size_type i=0; i<num_rows; ++i) {
186 for (size_type j=0; j<num_cols; ++j) {
187 scalar_type diff = std::abs( hy(i,j) - hy_exp(i,j) );
188 scalar_type tol = rel_tol*std::abs(hy_exp(i,j)) + abs_tol;
189 bool s = diff < tol;
190 out << "y_expected(" << i << "," << j << ") - "
191 << "y(" << i << "," << j << ") = " << hy_exp(i,j)
192 << " - " << hy(i,j) << " == "
193 << diff << " < " << tol << " : ";
194 if (s)
195 out << "passed";
196 else
197 out << "failed";
198 out << std::endl;
199 success = success && s;
200 }
201 }
202
203 return success;
204}
205
206template <typename vector_type, typename scalar_type>
207bool compareRank1(const vector_type& y,
208 const vector_type& y_exp,
209 const scalar_type rel_tol,
210 const scalar_type abs_tol,
211 Teuchos::FancyOStream& out)
212{
213 typedef typename vector_type::size_type size_type;
214 typename vector_type::HostMirror hy = Kokkos::create_mirror_view(y);
215 typename vector_type::HostMirror hy_exp = Kokkos::create_mirror_view(y_exp);
216 Kokkos::deep_copy(hy, y);
217 Kokkos::deep_copy(hy_exp, y_exp);
218
219 size_type num_rows = y.extent(0);
220 bool success = true;
221 for (size_type i=0; i<num_rows; ++i) {
222 for (size_type j=0; j<Kokkos::dimension_scalar(y); ++j) {
223 scalar_type diff = std::abs( hy(i).fastAccessCoeff(j) - hy_exp(i).fastAccessCoeff(j) );
224 scalar_type tol = rel_tol*std::abs(hy_exp(i).fastAccessCoeff(j)) + abs_tol;
225 bool s = diff < tol;
226 out << "y_expected(" << i << ").coeff(" << j << ") - "
227 << "y(" << i << ").coeff(" << j << ") = " << hy_exp(i).fastAccessCoeff(j)
228 << " - " << hy(i).fastAccessCoeff(j) << " == "
229 << diff << " < " << tol << " : ";
230 if (s)
231 out << "passed";
232 else
233 out << "failed";
234 out << std::endl;
235 success = success && s;
236 }
237 }
238 return success;
239}
240
241template <typename vector_type, typename scalar_type>
242bool compareRank2(const vector_type& y,
243 const vector_type& y_exp,
244 const scalar_type rel_tol,
245 const scalar_type abs_tol,
246 Teuchos::FancyOStream& out)
247{
248 typedef typename vector_type::size_type size_type;
249 typename vector_type::HostMirror hy = Kokkos::create_mirror_view(y);
250 typename vector_type::HostMirror hy_exp = Kokkos::create_mirror_view(y_exp);
251 Kokkos::deep_copy(hy, y);
252 Kokkos::deep_copy(hy_exp, y_exp);
253
254 size_type num_rows = y.extent(0);
255 size_type num_cols = y.extent(1);
256 bool success = true;
257
258 for (size_type col = 0; col < num_cols; ++col){
259 for (size_type i=0; i<num_rows; ++i) {
260 for (size_type j=0; j<Kokkos::dimension_scalar(y); ++j) {
261 scalar_type diff = std::abs( hy(i,col).fastAccessCoeff(j) - hy_exp(i,col).fastAccessCoeff(j) );
262 scalar_type tol = rel_tol*std::abs(hy_exp(i,col).fastAccessCoeff(j)) + abs_tol;
263 bool s = diff < tol;
264 out << "y_expected(" << i << ").coeff(" << j << ") - "
265 << "y(" << i << ").coeff(" << j << ") = " << hy_exp(i,col).fastAccessCoeff(j)
266 << " - " << hy(i,col).fastAccessCoeff(j) << " == "
267 << diff << " < " << tol << " : ";
268 if (s)
269 out << "passed";
270 else
271 out << "failed";
272 out << std::endl;
273 success = success && s;
274 }
275 }
276 }
277
278
279 return success;
280}
281
282
283// Helper function to build a diagonal matrix
284template <typename MatrixType, typename CijkType>
285MatrixType
286buildDiagonalMatrix(typename MatrixType::ordinal_type nrow,
287 typename MatrixType::ordinal_type pce_size,
288 const CijkType& cijk) {
289 typedef typename MatrixType::ordinal_type ordinal_type;
290 typedef typename MatrixType::StaticCrsGraphType matrix_graph_type;
291 typedef typename MatrixType::values_type matrix_values_type;
292
293 std::vector< std::vector<ordinal_type> > graph(nrow);
294 for (ordinal_type i=0; i<nrow; ++i)
295 graph[i] = std::vector<ordinal_type>(1, i);
296 ordinal_type graph_length = nrow;
297
298 matrix_graph_type matrix_graph =
299 Kokkos::create_staticcrsgraph<matrix_graph_type>("graph", graph);
300 matrix_values_type matrix_values =
301 Kokkos::make_view<matrix_values_type>("values", cijk, graph_length, pce_size);
302
303 MatrixType matrix("matrix", nrow, matrix_values, matrix_graph);
304 return matrix;
305}
306
307//
308// Tests
309//
310
311// Kernel to set diagonal of a matrix to prescribed values
312template <typename MatrixType>
314 typedef typename MatrixType::execution_space execution_space;
315 typedef typename MatrixType::size_type size_type;
316 typedef typename MatrixType::value_type value_type;
317 typedef typename MatrixType::ordinal_type ordinal_type;
318
319 const MatrixType m_matrix;
320 ReplaceDiagonalValuesKernel(const MatrixType matrix) : m_matrix(matrix) {};
321
322 // Replace diagonal entry for row 'i' with a value
323 KOKKOS_INLINE_FUNCTION
324 void operator() (const size_type i) const {
325 const ordinal_type row = i;
326 const ordinal_type col = i;
328 m_matrix.replaceValues(row, &col, 1, &val, false, true);
329 }
330
331 // Kernel launch
332 static void apply(const MatrixType matrix) {
333 const size_type nrow = matrix.numRows();
334 Kokkos::parallel_for( nrow, ReplaceDiagonalValuesKernel(matrix) );
335 }
336
337 // Check the result is as expected
338 static bool check(const MatrixType matrix,
339 Teuchos::FancyOStream& out) {
340 typedef typename MatrixType::values_type matrix_values_type;
341 typename matrix_values_type::HostMirror host_matrix_values =
342 Kokkos::create_mirror_view(matrix.values);
343 Kokkos::deep_copy(host_matrix_values, matrix.values);
344 const ordinal_type nrow = matrix.numRows();
345 bool success = true;
346 value_type val_expected(Kokkos::cijk(matrix.values));
347 for (ordinal_type row=0; row<nrow; ++row) {
348 val_expected = row;
349 bool s = compareVecs(host_matrix_values(row),
350 "matrix_values(row)",
351 val_expected,
352 "val_expected",
353 0.0, 0.0, out);
354 success = success && s;
355 }
356 return success;
357 }
358};
359
360// Kernel to add values to the diagonal of a matrix
361template <typename MatrixType>
363 typedef typename MatrixType::execution_space execution_space;
364 typedef typename MatrixType::size_type size_type;
365 typedef typename MatrixType::value_type value_type;
366 typedef typename MatrixType::ordinal_type ordinal_type;
367
368 const MatrixType m_matrix;
369 AddDiagonalValuesKernel(const MatrixType matrix) : m_matrix(matrix) {};
370
371 // Replace diagonal entry for row 'i' with a value
372 KOKKOS_INLINE_FUNCTION
373 void operator() (const size_type i) const {
374 const ordinal_type row = i;
375 const ordinal_type col = i;
377 m_matrix.sumIntoValues(row, &col, 1, &val, false, true);
378 }
379
380 // Kernel launch
381 static void apply(const MatrixType matrix) {
382 const size_type nrow = matrix.numRows();
383 Kokkos::parallel_for( nrow, AddDiagonalValuesKernel(matrix) );
384 }
385
386 // Check the result is as expected
387 static bool check(const MatrixType matrix,
388 Teuchos::FancyOStream& out) {
389 typedef typename MatrixType::values_type matrix_values_type;
390 typename matrix_values_type::HostMirror host_matrix_values =
391 Kokkos::create_mirror_view(matrix.values);
392 Kokkos::deep_copy(host_matrix_values, matrix.values);
393 const ordinal_type nrow = matrix.numRows();
394 bool success = true;
395 value_type val_expected(Kokkos::cijk(matrix.values));
396 for (ordinal_type row=0; row<nrow; ++row) {
397 val_expected = row;
398 bool s = compareVecs(host_matrix_values(row),
399 "matrix_values(row)",
400 val_expected,
401 "val_expected",
402 0.0, 0.0, out);
403 success = success && s;
404 }
405 return success;
406 }
407};
408
409// Kernel to add values to the diagonal of a matrix where each thread
410// adds to the same row (checks atomic really works)
411template <typename MatrixType>
413 typedef typename MatrixType::execution_space execution_space;
414 typedef typename MatrixType::size_type size_type;
415 typedef typename MatrixType::value_type value_type;
416 typedef typename MatrixType::ordinal_type ordinal_type;
417
418 const MatrixType m_matrix;
419 AddDiagonalValuesAtomicKernel(const MatrixType matrix) : m_matrix(matrix) {};
420
421 // Replace diagonal entry for row 'i' with a value
422 KOKKOS_INLINE_FUNCTION
423 void operator() (const size_type i) const {
424 const ordinal_type row = 0;
425 const ordinal_type col = 0;
427 m_matrix.sumIntoValues(row, &col, 1, &val, false, true);
428 }
429
430 // Kernel launch
431 static void apply(const MatrixType matrix) {
432 const size_type nrow = matrix.numRows();
433 Kokkos::parallel_for( nrow, AddDiagonalValuesAtomicKernel(matrix) );
434 }
435
436 // Check the result is as expected
437 static bool check(const MatrixType matrix,
438 Teuchos::FancyOStream& out) {
439 typedef typename MatrixType::values_type matrix_values_type;
440 typename matrix_values_type::HostMirror host_matrix_values =
441 Kokkos::create_mirror_view(matrix.values);
442 Kokkos::deep_copy(host_matrix_values, matrix.values);
443 const ordinal_type nrow = matrix.numRows();
444 bool success = true;
445 value_type val_expected(Kokkos::cijk(matrix.values));
446 for (ordinal_type row=0; row<nrow; ++row) {
448 if (row == 0)
449 val_expected = nrow*(nrow-1)/2;
450 else
451 val_expected = 0.0;
452 bool s = compareVecs(host_matrix_values(row),
453 "matrix_values(row)",
454 val_expected,
455 "val_expected",
456 0.0, 0.0, out);
457 success = success && s;
458 }
459 return success;
460 }
461};
462
464 Kokkos_CrsMatrix_PCE, ReplaceValues, MatrixScalar )
465{
466 typedef typename MatrixScalar::ordinal_type Ordinal;
467 typedef typename MatrixScalar::execution_space Device;
468 typedef typename MatrixScalar::cijk_type Cijk;
469 typedef KokkosSparse::CrsMatrix<MatrixScalar,Ordinal,Device> Matrix;
470
471 // Build Cijk tensor
472 const Ordinal stoch_dim = 2;
473 const Ordinal poly_ord = 3;
474 Cijk cijk = build_cijk<Cijk>(stoch_dim, poly_ord);
475
476 // Build diagonal matrix
477 const Ordinal nrow = 10;
478 const Ordinal pce_size = cijk.dimension();
479 Matrix matrix = buildDiagonalMatrix<Matrix>(nrow, pce_size, cijk);
480
481 // Launch our kernel
483 kernel::apply(matrix);
484
485 // Check the result
486 success = kernel::check(matrix, out);
487}
488
490 Kokkos_CrsMatrix_PCE, SumIntoValues, MatrixScalar )
491{
492 typedef typename MatrixScalar::ordinal_type Ordinal;
493 typedef typename MatrixScalar::execution_space Device;
494 typedef typename MatrixScalar::cijk_type Cijk;
495 typedef KokkosSparse::CrsMatrix<MatrixScalar,Ordinal,Device> Matrix;
496
497 // Build Cijk tensor
498 const Ordinal stoch_dim = 2;
499 const Ordinal poly_ord = 3;
500 Cijk cijk = build_cijk<Cijk>(stoch_dim, poly_ord);
501
502 // Build diagonal matrix
503 const Ordinal nrow = 10;
504 const Ordinal pce_size = cijk.dimension();
505 Matrix matrix = buildDiagonalMatrix<Matrix>(nrow, pce_size, cijk);
506
507 // Launch our kernel
508 typedef AddDiagonalValuesKernel<Matrix> kernel;
509 kernel::apply(matrix);
510
511 // Check the result
512 success = kernel::check(matrix, out);
513}
514
516 Kokkos_CrsMatrix_PCE, SumIntoValuesAtomic, MatrixScalar )
517{
518 typedef typename MatrixScalar::ordinal_type Ordinal;
519 typedef typename MatrixScalar::execution_space Device;
520 typedef typename MatrixScalar::cijk_type Cijk;
521 typedef KokkosSparse::CrsMatrix<MatrixScalar,Ordinal,Device> Matrix;
522
523 // Build Cijk tensor
524 const Ordinal stoch_dim = 2;
525 const Ordinal poly_ord = 3;
526 Cijk cijk = build_cijk<Cijk>(stoch_dim, poly_ord);
527
528 // Build diagonal matrix
529 const Ordinal nrow = 10;
530 const Ordinal pce_size = cijk.dimension();
531 Matrix matrix = buildDiagonalMatrix<Matrix>(nrow, pce_size, cijk);
532
533 // Launch our kernel
535 kernel::apply(matrix);
536
537 // Check the result
538 success = kernel::check(matrix, out);
539}
540
541template <typename PCEType, typename Multiply>
542bool test_embedded_pce(const typename PCEType::ordinal_type nGrid,
543 const typename PCEType::ordinal_type stoch_dim,
544 const typename PCEType::ordinal_type poly_ord,
545 KokkosSparse::DeviceConfig dev_config,
546 Multiply multiply_op,
547 Teuchos::FancyOStream& out)
548{
549 typedef typename PCEType::ordinal_type ordinal_type;
550 typedef typename PCEType::value_type scalar_type;
551 typedef typename PCEType::storage_type storage_type;
552 typedef typename PCEType::cijk_type cijk_type;
553 typedef typename storage_type::execution_space execution_space;
554 typedef Kokkos::LayoutLeft Layout;
555 typedef Kokkos::View< PCEType*, Layout, execution_space > block_vector_type;
556 typedef KokkosSparse::CrsMatrix< PCEType, ordinal_type, execution_space > block_matrix_type;
557 typedef typename block_matrix_type::StaticCrsGraphType matrix_graph_type;
558 typedef typename block_matrix_type::values_type matrix_values_type;
559
560 // Build Cijk tensor
561 cijk_type cijk = build_cijk<cijk_type>(stoch_dim, poly_ord);
562 const ordinal_type stoch_length = cijk.dimension();
563 // const ordinal_type align = 8;
564 // const ordinal_type stoch_length_aligned = (stoch_length+align-1) & ~(align-1);
565 const ordinal_type stoch_length_aligned = stoch_length;
566
567 // Check pce_length == storage_type::static_size for static storage
568 TEUCHOS_TEST_FOR_EXCEPTION(
569 storage_type::is_static && storage_type::static_size != stoch_length,
570 std::logic_error,
571 "Static storage size must equal pce size");
572
573 // Generate FEM graph:
574 const ordinal_type fem_length = nGrid * nGrid * nGrid;
575 std::vector< std::vector<ordinal_type> > fem_graph;
576 const ordinal_type fem_graph_length = generate_fem_graph( nGrid, fem_graph );
577
578 //------------------------------
579 // Generate input/output multivectors -- Sacado dimension is always last,
580 // regardless of LayoutLeft/Right
581
582 block_vector_type x =
583 Kokkos::make_view<block_vector_type>("x", cijk, fem_length, stoch_length_aligned);
584 block_vector_type y =
585 Kokkos::make_view<block_vector_type>("y", cijk, fem_length, stoch_length_aligned);
586
587 typename block_vector_type::HostMirror hx = Kokkos::create_mirror_view( x );
588 typename block_vector_type::HostMirror hy = Kokkos::create_mirror_view( y );
589
590 // View the block vector as an array of the embedded intrinsic type.
591 typename block_vector_type::HostMirror::array_type hax = hx ;
592 typename block_vector_type::HostMirror::array_type hay = hy ;
593
594 for (ordinal_type iRowFEM=0; iRowFEM<fem_length; ++iRowFEM) {
595 for (ordinal_type iRowStoch=0; iRowStoch<stoch_length; ++iRowStoch) {
596 hax(iRowStoch, iRowFEM) =
597 generate_vector_coefficient<scalar_type>(
598 fem_length, stoch_length, iRowFEM, iRowStoch );
599 hay(iRowStoch, iRowFEM) = 0.0;
600 }
601 }
602
603 Kokkos::deep_copy( x, hx );
604 Kokkos::deep_copy( y, hy );
605
606 //------------------------------
607 // Generate block matrix -- it is always LayoutRight (currently)
608
609 matrix_graph_type matrix_graph =
610 Kokkos::create_staticcrsgraph<matrix_graph_type>(
611 std::string("test crs graph"), fem_graph);
612 matrix_values_type matrix_values =
613 Kokkos::make_view<matrix_values_type>(
614 Kokkos::ViewAllocateWithoutInitializing("matrix"), cijk, fem_graph_length, stoch_length_aligned);
615 block_matrix_type matrix(
616 "block_matrix", fem_length, matrix_values, matrix_graph);
617 matrix.dev_config = dev_config;
618
619 typename matrix_values_type::HostMirror hM =
620 Kokkos::create_mirror_view( matrix.values );
621
622 typename matrix_values_type::HostMirror::array_type haM = hM ;
623
624 for (ordinal_type iRowFEM=0, iEntryFEM=0; iRowFEM<fem_length; ++iRowFEM) {
625 const ordinal_type row_size = fem_graph[iRowFEM].size();
626 for (ordinal_type iRowEntryFEM=0; iRowEntryFEM<row_size;
627 ++iRowEntryFEM, ++iEntryFEM) {
628 const ordinal_type iColFEM = fem_graph[iRowFEM][iRowEntryFEM];
629
630 for (ordinal_type k=0; k<stoch_length; ++k) {
631 haM(iEntryFEM, k) =
632 generate_matrix_coefficient<scalar_type>(
633 fem_length, stoch_length, iRowFEM, iColFEM, k);
634 }
635 }
636 }
637
638 Kokkos::deep_copy( matrix.values, hM );
639
640 //------------------------------
641 // multiply
642
643 multiply_op( matrix, x, y );
644
645 //------------------------------
646 // generate correct answer
647
648 typedef typename block_vector_type::array_type array_type;
649 array_type ay_expected =
650 array_type("ay_expected", stoch_length_aligned, fem_length);
651 typename array_type::HostMirror hay_expected =
652 Kokkos::create_mirror_view(ay_expected);
653 typename cijk_type::HostMirror host_cijk =
654 Kokkos::create_mirror_view(cijk);
655 Kokkos::deep_copy(host_cijk, cijk);
656 for (ordinal_type iRowFEM=0, iEntryFEM=0; iRowFEM<fem_length; ++iRowFEM) {
657 const ordinal_type row_size = fem_graph[iRowFEM].size();
658 for (ordinal_type iRowEntryFEM=0; iRowEntryFEM<row_size;
659 ++iRowEntryFEM, ++iEntryFEM) {
660 const ordinal_type iColFEM = fem_graph[iRowFEM][iRowEntryFEM];
661 for (ordinal_type i=0; i<stoch_length; ++i) {
662 const ordinal_type num_entry = host_cijk.num_entry(i);
663 const ordinal_type entry_beg = host_cijk.entry_begin(i);
664 const ordinal_type entry_end = entry_beg + num_entry;
665 scalar_type tmp = 0;
666 for (ordinal_type entry = entry_beg; entry < entry_end; ++entry) {
667 const ordinal_type j = host_cijk.coord(entry,0);
668 const ordinal_type k = host_cijk.coord(entry,1);
669 const scalar_type a_j =
670 generate_matrix_coefficient<scalar_type>(
671 fem_length, stoch_length, iRowFEM, iColFEM, j);
672 const scalar_type a_k =
673 generate_matrix_coefficient<scalar_type>(
674 fem_length, stoch_length, iRowFEM, iColFEM, k);
675 const scalar_type x_j =
676 generate_vector_coefficient<scalar_type>(
677 fem_length, stoch_length, iColFEM, j);
678 const scalar_type x_k =
679 generate_vector_coefficient<scalar_type>(
680 fem_length, stoch_length, iColFEM, k);
681 tmp += host_cijk.value(entry) * ( a_j * x_k + a_k * x_j );
682 }
683 hay_expected(i, iRowFEM) += tmp;
684 }
685 }
686 }
687 Kokkos::deep_copy( ay_expected, hay_expected );
688
689 //------------------------------
690 // check
691
692 typename block_vector_type::array_type ay = y;
693 scalar_type rel_tol = ScalarTol<scalar_type>::tol();
694 scalar_type abs_tol = ScalarTol<scalar_type>::tol();
695 bool success = compare_rank_2_views(ay, ay_expected, rel_tol, abs_tol, out);
696
697 return success;
698}
699
701 Kokkos_CrsMatrix_PCE, Multiply, Scalar, MultiplyOp )
702{
703 typedef typename Scalar::ordinal_type Ordinal;
704
705 const Ordinal nGrid = 5;
706 const Ordinal stoch_dim = 2;
707 const Ordinal poly_ord = 3;
708 KokkosSparse::DeviceConfig dev_config;
709
710 success = test_embedded_pce<Scalar>(
711 nGrid, stoch_dim, poly_ord, dev_config, MultiplyOp(), out);
712}
713
715 template <typename Matrix, typename InputVector, typename OutputVector>
716 void operator() (const Matrix& A,
717 const InputVector& x,
718 OutputVector& y) const {
719 KokkosSparse::spmv("N", typename Matrix::value_type(1.0) , A, x, typename Matrix::value_type(0.0), y);
720 }
721};
722
723template <typename Tag>
725 Tag tag;
726 Stokhos_MV_Multiply_Op(const Tag& tg = Tag()) : tag(tg) {}
727
728 template <typename Matrix, typename InputVector, typename OutputVector>
729 void operator() (const Matrix& A,
730 const InputVector& x,
731 OutputVector& y) const {
732 Stokhos::multiply(A, x, y, tag);
733 }
734};
735
737 Kokkos_CrsMatrix_PCE, MeanMultiplyRank1, Scalar )
738{
739 typedef typename Scalar::ordinal_type Ordinal;
740
741 const Ordinal nGrid = 5;
742 const Ordinal stoch_dim = 5;
743 const Ordinal poly_ord = 3;
744 KokkosSparse::DeviceConfig dev_config;
745
746 typedef typename Scalar::ordinal_type ordinal_type;
747 typedef typename Scalar::value_type scalar_type;
748 typedef typename Scalar::storage_type storage_type;
749 typedef typename Scalar::cijk_type cijk_type;
750 typedef typename storage_type::execution_space execution_space;
751 typedef Kokkos::LayoutLeft Layout;
752 typedef Kokkos::View< Scalar*, Layout, execution_space > block_vector_type;
753 typedef KokkosSparse::CrsMatrix< Scalar, ordinal_type, execution_space > block_matrix_type;
754 typedef typename block_matrix_type::StaticCrsGraphType matrix_graph_type;
755 typedef typename block_matrix_type::values_type matrix_values_type;
756
757 // Build Cijk tensor
758 cijk_type cijk = build_cijk<cijk_type>(stoch_dim, poly_ord);
759 cijk_type mean_cijk =
760 Stokhos::create_mean_based_product_tensor<execution_space, ordinal_type, scalar_type>();
761 const ordinal_type stoch_length = cijk.dimension();
762 const ordinal_type align = 8;
763 const ordinal_type stoch_length_aligned = (stoch_length+align-1) & ~(align-1);
764
765 // Generate FEM graph:
766 const ordinal_type fem_length = nGrid * nGrid * nGrid;
767 std::vector< std::vector<ordinal_type> > fem_graph;
768 const ordinal_type fem_graph_length = generate_fem_graph( nGrid, fem_graph );
769
770 block_vector_type x =
771 Kokkos::make_view<block_vector_type>("x", cijk, fem_length, stoch_length_aligned);
772 block_vector_type y =
773 Kokkos::make_view<block_vector_type>("y", cijk, fem_length, stoch_length_aligned);
774
775 block_vector_type y_expected =
776 Kokkos::make_view<block_vector_type>("y", cijk, fem_length, stoch_length_aligned);
777
778 typename block_vector_type::HostMirror hx = Kokkos::create_mirror_view( x );
779 typename block_vector_type::HostMirror hy = Kokkos::create_mirror_view( y );
780 typename block_vector_type::HostMirror hy_expected =
781 Kokkos::create_mirror_view( y_expected );
782
783 // View the block vector as an array of the embedded intrinsic type.
784 typename block_vector_type::HostMirror::array_type hax = hx ;
785 typename block_vector_type::HostMirror::array_type hay = hy ;
786 typename block_vector_type::HostMirror::array_type hay_expected =
787 hy_expected ;
788
789 for (ordinal_type iRowFEM=0; iRowFEM<fem_length; ++iRowFEM) {
790 for (ordinal_type iRowStoch=0; iRowStoch<stoch_length; ++iRowStoch) {
791 hax(iRowStoch,iRowFEM) =
792 generate_vector_coefficient<scalar_type>(
793 fem_length, stoch_length, iRowFEM, iRowStoch );
794 hay(iRowStoch,iRowFEM) = 0.0;
795 hay_expected(iRowStoch,iRowFEM) = 0.0;
796 }
797 }
798 Kokkos::deep_copy( x, hx );
799 Kokkos::deep_copy( y, hy );
800 Kokkos::deep_copy( y_expected, hy_expected );
801
802 //------------------------------
803 // Generate block matrix -- it is always LayoutRight (currently)
804 matrix_graph_type matrix_graph =
805 Kokkos::create_staticcrsgraph<matrix_graph_type>(
806 std::string("test crs graph"), fem_graph);
807 matrix_values_type matrix_values =
808 Kokkos::make_view<matrix_values_type>(
809 Kokkos::ViewAllocateWithoutInitializing("matrix"), mean_cijk, fem_graph_length, ordinal_type(1)); //instead of stoch_length
810 block_matrix_type matrix(
811 "block_matrix", fem_length, matrix_values, matrix_graph);
812 matrix.dev_config = dev_config;
813
814 typename matrix_values_type::HostMirror hM =
815 Kokkos::create_mirror_view( matrix.values );
816 typename matrix_values_type::HostMirror::array_type haM = hM ;
817
818 for (ordinal_type iRowFEM=0, iEntryFEM=0; iRowFEM<fem_length; ++iRowFEM) {
819 const ordinal_type row_size = fem_graph[iRowFEM].size();
820 for (ordinal_type iRowEntryFEM=0; iRowEntryFEM<row_size;
821 ++iRowEntryFEM, ++iEntryFEM) {
822 const ordinal_type iColFEM = fem_graph[iRowFEM][iRowEntryFEM];
823
824 haM(iEntryFEM, 0) =
825 generate_matrix_coefficient<scalar_type>(
826 fem_length, 1, iRowFEM, iColFEM, 0);
827 for (ordinal_type iRowStoch=0; iRowStoch<stoch_length; ++iRowStoch) {
828 hay_expected(iRowStoch,iRowFEM) +=
829 haM(iEntryFEM, 0) * hax(iRowStoch,iColFEM);
830 }
831 }
832 }
833 Kokkos::deep_copy( matrix.values, hM );
834 Kokkos::deep_copy( y_expected, hy_expected );
835
836 /*
837 //Generate same matrix with stochastic dim = Kokkos::dimension_scalar(x) (i.e. not = 1)
838 matrix_values_type full_matrix_values =
839 Kokkos::make_view<matrix_values_type>(
840 Kokkos::ViewAllocateWithoutInitializing("matrix"), cijk, fem_graph_length, stoch_length_aligned);
841 block_matrix_type full_matrix(
842 "block_matrix", fem_length, full_matrix_values, matrix_graph);
843 matrix.dev_config = dev_config;
844
845 typename matrix_values_type::HostMirror full_hM =
846 Kokkos::create_mirror_view( full_matrix.values );
847 typename matrix_values_type::HostMirror::array_type full_haM = full_hM ;
848
849 for (ordinal_type iRowFEM=0, iEntryFEM=0; iRowFEM<fem_length; ++iRowFEM) {
850 const ordinal_type row_size = fem_graph[iRowFEM].size();
851 for (ordinal_type iRowEntryFEM=0; iRowEntryFEM<row_size;
852 ++iRowEntryFEM, ++iEntryFEM) {
853 const ordinal_type iColFEM = fem_graph[iRowFEM][iRowEntryFEM];
854
855 for (ordinal_type k=0; k<stoch_length; ++k) {
856 if (k == 0)
857 full_haM(iEntryFEM, k) =
858 generate_matrix_coefficient<scalar_type>(
859 fem_length, 1, iRowFEM, iColFEM, k);
860 else
861 full_haM(iEntryFEM, k) = 0.0;
862 }
863 }
864 }
865
866 Kokkos::deep_copy( full_matrix.values, full_hM );
867 */
868
869 //------------------------------
870 // multiply
871
872 KokkosSparse::spmv("N", Scalar(1.0) , matrix, x, Scalar(0.0), y);
873
874 //------------------------------
875 // multiply with same matrix but with sacado_size = x.sacado_size
876
877 //Kokkos::MV_Multiply( y_expected, full_matrix, x );
878
879 //------------------------------
880 // check
881 scalar_type rel_tol = ScalarTol<scalar_type>::tol();
882 scalar_type abs_tol = ScalarTol<scalar_type>::tol();
883 success = compareRank1(y, y_expected, rel_tol, abs_tol, out);
884}
885
886
888 Kokkos_CrsMatrix_PCE, MeanMultiplyRank2, Scalar )
889{
890 typedef typename Scalar::ordinal_type ordinal_type;
891 typedef typename Scalar::value_type scalar_type;
892 typedef typename Scalar::storage_type storage_type;
893 typedef typename Scalar::cijk_type cijk_type;
894 typedef typename storage_type::execution_space execution_space;
895 typedef Kokkos::LayoutLeft Layout;
896 typedef Kokkos::View< Scalar**, Layout, execution_space > block_vector_type;
897 typedef KokkosSparse::CrsMatrix< Scalar, ordinal_type, execution_space > block_matrix_type;
898 typedef typename block_matrix_type::StaticCrsGraphType matrix_graph_type;
899 typedef typename block_matrix_type::values_type matrix_values_type;
900
901
902 const ordinal_type nGrid = 5;
903 const ordinal_type stoch_dim = 2;
904 const ordinal_type poly_ord = 3;
905 KokkosSparse::DeviceConfig dev_config;
906
907 // Build Cijk tensor
908 cijk_type cijk = build_cijk<cijk_type>(stoch_dim, poly_ord);
909 cijk_type mean_cijk =
910 Stokhos::create_mean_based_product_tensor<execution_space, ordinal_type, scalar_type>();
911 const ordinal_type stoch_length = cijk.dimension();
912 const ordinal_type align = 8;
913 const ordinal_type stoch_length_aligned = (stoch_length+align-1) & ~(align-1);
914 const ordinal_type num_cols = 2;
915 // Generate FEM graph:
916 const ordinal_type fem_length = nGrid * nGrid * nGrid;
917 std::vector< std::vector<ordinal_type> > fem_graph;
918 const ordinal_type fem_graph_length = generate_fem_graph( nGrid, fem_graph );
919
920 block_vector_type x =
921 Kokkos::make_view<block_vector_type>("x", cijk, fem_length, num_cols, stoch_length_aligned);
922 block_vector_type y =
923 Kokkos::make_view<block_vector_type>("y", cijk, fem_length, num_cols, stoch_length_aligned);
924
925 block_vector_type y_expected =
926 Kokkos::make_view<block_vector_type>("y_expected", cijk, fem_length, num_cols, stoch_length_aligned);
927
928 typename block_vector_type::HostMirror hx = Kokkos::create_mirror_view( x );
929 typename block_vector_type::HostMirror hy = Kokkos::create_mirror_view( y );
930 typename block_vector_type::HostMirror hy_expected =
931 Kokkos::create_mirror_view( y_expected );
932
933 for (ordinal_type i=0; i<num_cols; ++i){
934 for (ordinal_type iRowFEM=0; iRowFEM<fem_length; ++iRowFEM) {
935 for (ordinal_type iRowStoch=0; iRowStoch<stoch_length; ++iRowStoch) {
936 hx(iRowFEM,i).fastAccessCoeff(iRowStoch) =
937 generate_vector_coefficient<scalar_type>(
938 fem_length, stoch_length, iRowFEM, iRowStoch );
939 hy(iRowFEM,i).fastAccessCoeff(iRowStoch) = 0.0;
940 hy_expected(iRowFEM,i).fastAccessCoeff(iRowStoch) = 0.0;
941 }
942 }
943 }
944 Kokkos::deep_copy( x, hx );
945 Kokkos::deep_copy( y, hy );
946 Kokkos::deep_copy( y_expected, hy_expected );
947
948 //------------------------------
949 // Generate matrix with stochastic dimension 1
950 matrix_graph_type matrix_graph =
951 Kokkos::create_staticcrsgraph<matrix_graph_type>(
952 std::string("test crs graph"), fem_graph);
953 matrix_values_type matrix_values =
954 Kokkos::make_view<matrix_values_type>(
955 Kokkos::ViewAllocateWithoutInitializing("matrix"), mean_cijk, fem_graph_length, ordinal_type(1));
956 block_matrix_type matrix(
957 "block_matrix", fem_length, matrix_values, matrix_graph);
958 matrix.dev_config = dev_config;
959
960 typename matrix_values_type::HostMirror hM =
961 Kokkos::create_mirror_view( matrix.values );
962
963 typename matrix_values_type::HostMirror::array_type haM = hM ;
964
965 for (ordinal_type iRowFEM=0, iEntryFEM=0; iRowFEM<fem_length; ++iRowFEM) {
966 const ordinal_type row_size = fem_graph[iRowFEM].size();
967 for (ordinal_type iRowEntryFEM=0; iRowEntryFEM<row_size;
968 ++iRowEntryFEM, ++iEntryFEM) {
969 const ordinal_type iColFEM = fem_graph[iRowFEM][iRowEntryFEM];
970
971 haM(iEntryFEM, 0) =
972 generate_matrix_coefficient<scalar_type>(
973 fem_length, 1, iRowFEM, iColFEM, 0);
974 for (ordinal_type i=0; i<num_cols; ++i){
975 for (ordinal_type iRowStoch=0; iRowStoch<stoch_length; ++iRowStoch) {
976 hy_expected(iRowFEM,i).fastAccessCoeff(iRowStoch) +=
977 haM(iEntryFEM, 0) * hx(iColFEM,i).fastAccessCoeff(iRowStoch);
978 }
979 }
980 }
981 }
982
983 Kokkos::deep_copy( matrix.values, hM );
984 Kokkos::deep_copy( y_expected, hy_expected );
985
986 /*
987 //Generate same matrix with stochastic dim = Kokkos::dimension_scalar(x) (i.e. not = 1)
988 matrix_values_type full_matrix_values =
989 Kokkos::make_view<matrix_values_type>(
990 Kokkos::ViewAllocateWithoutInitializing("matrix"), cijk, fem_graph_length, stoch_length_aligned);
991 block_matrix_type full_matrix(
992 "block_matrix", fem_length, full_matrix_values, matrix_graph);
993 matrix.dev_config = dev_config;
994
995 typename matrix_values_type::HostMirror full_hM =
996 Kokkos::create_mirror_view( full_matrix.values );
997
998 typename matrix_values_type::HostMirror::array_type full_haM = full_hM ;
999
1000 for (ordinal_type iRowFEM=0, iEntryFEM=0; iRowFEM<fem_length; ++iRowFEM) {
1001 const ordinal_type row_size = fem_graph[iRowFEM].size();
1002 for (ordinal_type iRowEntryFEM=0; iRowEntryFEM<row_size;
1003 ++iRowEntryFEM, ++iEntryFEM) {
1004 const ordinal_type iColFEM = fem_graph[iRowFEM][iRowEntryFEM];
1005
1006 for (ordinal_type k=0; k<stoch_length; ++k) {
1007 if (k == 0)
1008 full_haM(iEntryFEM, k) =
1009 generate_matrix_coefficient<scalar_type>(
1010 fem_length, 1, iRowFEM, iColFEM, k);
1011 else
1012 full_haM(iEntryFEM, k) = 0.0;
1013 }
1014 }
1015 }
1016
1017 Kokkos::deep_copy( full_matrix.values, full_hM );
1018 */
1019
1020 //------------------------------
1021 // multiply
1022
1023 KokkosSparse::spmv("N", Scalar(1.0) , matrix, x, Scalar(0.0), y);
1024
1025 //------------------------------
1026 // multiply with full matrix
1027
1028 //Kokkos::MV_Multiply( y_expected, full_matrix, x );
1029
1030 //------------------------------
1031 // check
1032
1033 scalar_type rel_tol = ScalarTol<scalar_type>::tol();
1034 scalar_type abs_tol = ScalarTol<scalar_type>::tol();
1035 success = compareRank2(y, y_expected, rel_tol, abs_tol, out);
1036}
1037
1038
1040
1041#define CRSMATRIX_UQ_PCE_TESTS_MATRIXSCALAR( SCALAR ) \
1042 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( \
1043 Kokkos_CrsMatrix_PCE, ReplaceValues, SCALAR ) \
1044 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( \
1045 Kokkos_CrsMatrix_PCE, SumIntoValues, SCALAR ) \
1046 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( \
1047 Kokkos_CrsMatrix_PCE, SumIntoValuesAtomic, SCALAR )
1048#define CRSMATRIX_UQ_PCE_MEAN_MULTIPLY_TESTS( SCALAR ) \
1049 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( \
1050 Kokkos_CrsMatrix_PCE, MeanMultiplyRank1, SCALAR ) \
1051 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( \
1052 Kokkos_CrsMatrix_PCE, MeanMultiplyRank2, SCALAR )
1053#define CRS_MATRIX_UQ_PCE_MULTIPLY_TESTS_SCALAR_OP( SCALAR, OP ) \
1054 TEUCHOS_UNIT_TEST_TEMPLATE_2_INSTANT( \
1055 Kokkos_CrsMatrix_PCE, Multiply, SCALAR, OP )
1056
1057#define CRS_MATRIX_UQ_PCE_MULTIPLY_TESTS_SCALAR( SCALAR ) \
1058 CRS_MATRIX_UQ_PCE_MULTIPLY_TESTS_SCALAR_OP( SCALAR, KokkosMultiply )
1059
1060#define CRSMATRIX_UQ_PCE_TESTS_STORAGE( STORAGE ) \
1061 typedef Sacado::UQ::PCE<STORAGE> UQ_PCE_ ## STORAGE; \
1062 CRSMATRIX_UQ_PCE_TESTS_MATRIXSCALAR( UQ_PCE_ ## STORAGE ) \
1063 CRS_MATRIX_UQ_PCE_MULTIPLY_TESTS_SCALAR( UQ_PCE_ ## STORAGE ) \
1064 CRSMATRIX_UQ_PCE_MEAN_MULTIPLY_TESTS( UQ_PCE_ ## STORAGE )
1065
1066#define CRSMATRIX_UQ_PCE_TESTS_ORDINAL_SCALAR_DEVICE( ORDINAL, SCALAR, DEVICE ) \
1067 typedef Stokhos::DynamicStorage<ORDINAL,SCALAR,DEVICE> DS; \
1068 CRSMATRIX_UQ_PCE_TESTS_STORAGE( DS )
1069
1070#define CRSMATRIX_UQ_PCE_TESTS_DEVICE( DEVICE ) \
1071 CRSMATRIX_UQ_PCE_TESTS_ORDINAL_SCALAR_DEVICE( int, double, DEVICE )
expr1 expr1 expr1 expr2 expr1 expr1 c expr2 expr1 c fastAccessCoeff(j) - expr2.val(j)
expr val()
bool test_embedded_pce(const typename PCEType::ordinal_type nGrid, const typename PCEType::ordinal_type stoch_dim, const typename PCEType::ordinal_type poly_ord, KokkosSparse::DeviceConfig dev_config, Multiply multiply_op, Teuchos::FancyOStream &out)
ordinal generate_fem_graph(ordinal N, std::vector< std::vector< ordinal > > &graph)
kokkos_cijk_type build_cijk(ordinal_type stoch_dim, ordinal_type poly_ord)
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)
TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL(Kokkos_CrsMatrix_PCE, ReplaceValues, MatrixScalar)
Kokkos_MV_Multiply_Op KokkosMultiply
MatrixType buildDiagonalMatrix(typename MatrixType::ordinal_type nrow, typename MatrixType::ordinal_type pce_size, const CijkType &cijk)
TEUCHOS_UNIT_TEST_TEMPLATE_2_DECL(Kokkos_CrsMatrix_PCE, Multiply, Scalar, MultiplyOp)
IntType map_fem_graph_coord(const IntType &N, const IntType &i, const IntType &j, const IntType &k)
scalar generate_vector_coefficient(const ordinal nFEM, const ordinal nStoch, const ordinal iColFEM, const ordinal iStoch)
bool compareRank2(const vector_type &y, const vector_type &y_exp, const scalar_type rel_tol, const scalar_type abs_tol, Teuchos::FancyOStream &out)
bool compareRank1(const vector_type &y, const vector_type &y_exp, const scalar_type rel_tol, const scalar_type abs_tol, Teuchos::FancyOStream &out)
scalar generate_matrix_coefficient(const ordinal nFEM, const ordinal nStoch, const ordinal iRowFEM, const ordinal iColFEM, const ordinal iStoch)
Stokhos::StandardStorage< int, double > storage_type
Kokkos::DefaultExecutionSpace execution_space
Multivariate orthogonal polynomial basis generated from a total-order complete-polynomial tensor prod...
Legendre polynomial basis.
Abstract base class for 1-D orthogonal polynomials.
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
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)
KOKKOS_INLINE_FUNCTION constexpr std::enable_if< is_view_uq_pce< View< T, P... > >::value, unsigned >::type dimension_scalar(const View< T, P... > &view)
KOKKOS_INLINE_FUNCTION constexpr std::enable_if< is_view_uq_pce< view_type >::value, typenameCijkType< view_type >::type >::type cijk(const view_type &view)
void multiply(const CrsMatrix< MatrixValue, Device, Layout > &A, const InputMultiVectorType &x, OutputMultiVectorType &y, const std::vector< OrdinalType > &col_indices, SingleColumnMultivectorMultiply)
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
static bool check(const MatrixType matrix, Teuchos::FancyOStream &out)
static void apply(const MatrixType matrix)
static void apply(const MatrixType matrix)
static bool check(const MatrixType matrix, Teuchos::FancyOStream &out)
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
void operator()(const Matrix &A, const InputVector &x, OutputVector &y) const
static void apply(const MatrixType matrix)
KOKKOS_INLINE_FUNCTION void operator()(const size_type i) const
static bool check(const MatrixType matrix, Teuchos::FancyOStream &out)
void operator()(const Matrix &A, const InputVector &x, OutputVector &y) const