Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Stokhos_KokkosArrayKernelsUnitTest.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 #ifndef STOKHOS_KOKKOS_CORE_KERNELS_UNIT_TEST_HPP
43 #define STOKHOS_KOKKOS_CORE_KERNELS_UNIT_TEST_HPP
44 
45 #include "Teuchos_TestingHelpers.hpp"
46 #include "Teuchos_UnitTestHelpers.hpp"
47 #include "Teuchos_TestForException.hpp"
48 
49 #include "Stokhos_Epetra.hpp"
50 #include "Stokhos_Sparse3TensorUtilities.hpp"
51 #include "EpetraExt_BlockUtility.h"
52 #include "Stokhos_UnitTestHelpers.hpp"
53 
54 #ifdef HAVE_MPI
55 #include "Epetra_MpiComm.h"
56 #else
57 #include "Epetra_SerialComm.h"
58 #endif
59 
60 #include "Kokkos_Macros.hpp"
61 
62 #include "Stokhos_Update.hpp"
63 #include "Stokhos_CrsMatrix.hpp"
64 #include "Stokhos_BlockCrsMatrix.hpp"
65 #include "Stokhos_StochasticProductTensor.hpp"
66 #include "Stokhos_CrsProductTensor.hpp"
67 #include "Stokhos_TiledCrsProductTensor.hpp"
68 #include "Stokhos_SimpleTiledCrsProductTensor.hpp"
69 #include "Stokhos_CooProductTensor.hpp"
70 #include "Stokhos_SymmetricDiagonalSpec.hpp"
71 #include "Stokhos_FlatSparse3Tensor.hpp"
72 #include "Stokhos_FlatSparse3Tensor_kji.hpp"
73 #include "Stokhos_LinearSparse3Tensor.hpp"
74 #include "Stokhos_LexicographicBlockSparse3Tensor.hpp"
75 
76 #ifdef HAVE_STOKHOS_KOKKOSLINALG
77 #include "KokkosSparse_CrsMatrix.hpp"
78 #include "KokkosSparse_spmv.hpp"
79 #include "KokkosBlas1_update.hpp"
80 #endif
81 
82 namespace KokkosKernelsUnitTest {
83 
84 using Teuchos::rcp;
85 using Teuchos::RCP;
86 using Teuchos::Array;
87 using Teuchos::ParameterList;
88 
89 template< typename IntType >
90 inline
91 IntType map_fem_graph_coord( const IntType & N ,
92  const IntType & i ,
93  const IntType & j ,
94  const IntType & k )
95 {
96  return k + N * ( j + N * i );
97 }
98 
99 template < typename ordinal >
100 inline
101 ordinal generate_fem_graph( ordinal N ,
102  std::vector< std::vector<ordinal> > & graph )
103 {
104  graph.resize( N * N * N , std::vector<ordinal>() );
105 
106  ordinal total = 0 ;
107 
108  for ( int i = 0 ; i < (int) N ; ++i ) {
109  for ( int j = 0 ; j < (int) N ; ++j ) {
110  for ( int k = 0 ; k < (int) N ; ++k ) {
111 
112  const ordinal row = map_fem_graph_coord((int)N,i,j,k);
113 
114  graph[row].reserve(27);
115 
116  for ( int ii = -1 ; ii < 2 ; ++ii ) {
117  for ( int jj = -1 ; jj < 2 ; ++jj ) {
118  for ( int kk = -1 ; kk < 2 ; ++kk ) {
119  if ( 0 <= i + ii && i + ii < (int) N &&
120  0 <= j + jj && j + jj < (int) N &&
121  0 <= k + kk && k + kk < (int) N ) {
122  ordinal col = map_fem_graph_coord((int)N,i+ii,j+jj,k+kk);
123 
124  graph[row].push_back(col);
125  }
126  }}}
127  total += graph[row].size();
128  }}}
129 
130  return total ;
131 }
132 
133 template <typename scalar, typename ordinal>
134 inline
135 scalar generate_matrix_coefficient( const ordinal nFEM ,
136  const ordinal nStoch ,
137  const ordinal iRowFEM ,
138  const ordinal iColFEM ,
139  const ordinal iStoch )
140 {
141  const scalar A_fem = ( 10.0 + scalar(iRowFEM) / scalar(nFEM) ) +
142  ( 5.0 + scalar(iColFEM) / scalar(nFEM) );
143 
144  const scalar A_stoch = ( 1.0 + scalar(iStoch) / scalar(nStoch) );
145 
146  return A_fem + A_stoch ;
147  //return 1.0;
148 }
149 
150 template <typename scalar, typename ordinal>
151 inline
152 scalar generate_vector_coefficient( const ordinal nFEM ,
153  const ordinal nStoch ,
154  const ordinal iColFEM ,
155  const ordinal iStoch )
156 {
157  const scalar X_fem = 100.0 + scalar(iColFEM) / scalar(nFEM);
158  const scalar X_stoch = 1.0 + scalar(iStoch) / scalar(nStoch);
159  return X_fem + X_stoch ;
160  //return 1.0;
161 }
162 
163 template <typename Device>
164 struct UnitTestSetup {
165  typedef double value_type ;
166  typedef Stokhos::OneDOrthogPolyBasis<int,value_type> abstract_basis_type;
167  typedef Stokhos::LegendreBasis<int,value_type> basis_type;
168  //typedef Stokhos::CompletePolynomialBasis<int,value_type> product_basis_type;
169  typedef Stokhos::LexographicLess< Stokhos::MultiIndex<int> > less_type;
170  typedef Stokhos::TotalOrderBasis<int,value_type,less_type> product_basis_type;
171  typedef Stokhos::Sparse3Tensor<int,value_type> Cijk_type;
172 
173  int p, d, nGrid, fem_length, stoch_length, stoch_length_aligned, fem_graph_length;
174  double rel_tol, abs_tol;
175  std::vector< std::vector<int> > fem_graph ;
176  RCP< product_basis_type> basis;
177  RCP<Cijk_type> Cijk;
178  RCP<const Epetra_Comm> globalComm;
179  RCP<Stokhos::EpetraVectorOrthogPoly> sg_x, sg_y;
180  RCP<Stokhos::ProductEpetraVector> sg_y_commuted;
181  Teuchos::Array<int> perm, inv_perm;
182 
183  // Can't be a constructor because MPI will not be initialized
184  void setup(int p_ = 5, int d_ = 2) {
185 
186  p = p_;
187  d = d_;
188  nGrid = 3;
189  rel_tol = 1e-12;
190  abs_tol = 1e-12;
191 
192  // Create a communicator for Epetra objects
193 #ifdef HAVE_MPI
194  globalComm = rcp(new Epetra_MpiComm(MPI_COMM_WORLD));
195 #else
196  globalComm = rcp(new Epetra_SerialComm);
197 #endif
198  //int MyPID = globalComm->MyPID();
199 
200  //------------------------------
201  // Generate FEM graph:
202 
203  fem_length = nGrid * nGrid * nGrid ;
204  fem_graph_length = generate_fem_graph( nGrid , fem_graph );
205 
206  // Create Stochastic Galerkin basis and expansion
207  Array< RCP<const abstract_basis_type> > bases(d);
208  for (int i=0; i<d; i++)
209  bases[i] = rcp(new basis_type(p,true));
210  basis = rcp(new product_basis_type(bases, 1e-12));
211  stoch_length = basis->size();
212  Cijk = basis->computeTripleProductTensor();
213 
214  // Create stochastic parallel distribution
215  ParameterList parallelParams;
216  RCP<Stokhos::ParallelData> sg_parallel_data =
217  rcp(new Stokhos::ParallelData(basis, Cijk, globalComm, parallelParams));
218  RCP<const EpetraExt::MultiComm> sg_comm =
219  sg_parallel_data->getMultiComm();
220  RCP<const Epetra_Comm> app_comm =
221  sg_parallel_data->getSpatialComm();
222  RCP<const Stokhos::EpetraSparse3Tensor> epetraCijk =
223  sg_parallel_data->getEpetraCijk();
224  RCP<const Epetra_BlockMap> stoch_row_map =
225  epetraCijk->getStochasticRowMap();
226 
227  //------------------------------
228 
229  // Generate Epetra objects
230  RCP<const Epetra_Map> x_map =
231  rcp(new Epetra_Map(fem_length, 0, *app_comm));
232  RCP<const Epetra_Map> sg_x_map =
233  rcp(EpetraExt::BlockUtility::GenerateBlockMap(
234  *x_map, *stoch_row_map, *sg_comm));
235 
236  RCP<Epetra_CrsGraph> graph = rcp(new Epetra_CrsGraph(Copy, *x_map, 27));
237  int *my_GIDs = x_map->MyGlobalElements();
238  int num_my_GIDs = x_map->NumMyElements();
239  for (int i=0; i<num_my_GIDs; ++i) {
240  int row = my_GIDs[i];
241  int num_indices = fem_graph[row].size();
242  int *indices = &(fem_graph[row][0]);
243  graph->InsertGlobalIndices(row, num_indices, indices);
244  }
245  graph->FillComplete();
246 
247  RCP<ParameterList> sg_op_params = rcp(new ParameterList);
248  RCP<Stokhos::MatrixFreeOperator> sg_A =
249  rcp(new Stokhos::MatrixFreeOperator(sg_comm, basis, epetraCijk,
250  x_map, x_map, sg_x_map, sg_x_map,
251  sg_op_params));
252  RCP<Epetra_BlockMap> sg_A_overlap_map =
253  rcp(new Epetra_LocalMap(
254  stoch_length, 0, *(sg_parallel_data->getStochasticComm())));
255  RCP< Stokhos::EpetraOperatorOrthogPoly > A_sg_blocks =
256  rcp(new Stokhos::EpetraOperatorOrthogPoly(
257  basis, sg_A_overlap_map, x_map, x_map, sg_x_map, sg_comm));
258  for (int i=0; i<stoch_length; i++) {
259  RCP<Epetra_CrsMatrix> A = rcp(new Epetra_CrsMatrix(Copy, *graph));
260 
261  for ( int iRowFEM = 0 , iEntryFEM = 0 ; iRowFEM < fem_length ; ++iRowFEM ) {
262  for ( size_t iRowEntryFEM = 0 ; iRowEntryFEM < fem_graph[iRowFEM].size() ; ++iRowEntryFEM , ++iEntryFEM ) {
263  const int iColFEM = fem_graph[iRowFEM][iRowEntryFEM] ;
264 
265  double v = generate_matrix_coefficient<double>(
266  fem_length , stoch_length , iRowFEM , iColFEM , i );
267  A->ReplaceGlobalValues(iRowFEM, 1, &v, &iColFEM);
268  }
269  }
270  A->FillComplete();
271  A_sg_blocks->setCoeffPtr(i, A);
272  }
273  sg_A->setupOperator(A_sg_blocks);
274 
275  sg_x = rcp(new Stokhos::EpetraVectorOrthogPoly(
276  basis, stoch_row_map, x_map, sg_x_map, sg_comm));
277  sg_y = rcp(new Stokhos::EpetraVectorOrthogPoly(
278  basis, stoch_row_map, x_map, sg_x_map, sg_comm));
279 
280  // Initialize vectors
281  for (int iColFEM=0; iColFEM < fem_length; ++iColFEM ) {
282  for (int iColStoch=0 ; iColStoch < stoch_length; ++iColStoch ) {
283  (*sg_x)[iColStoch][iColFEM] = generate_vector_coefficient<double>(
284  fem_length , stoch_length , iColFEM , iColStoch );
285  }
286  }
287  sg_y->init(0.0);
288 
289  // Apply operator
290  sg_A->Apply( *(sg_x->getBlockVector()), *(sg_y->getBlockVector()) );
291 
292  // Transpose y to commuted layout
293  sg_y_commuted = rcp(new Stokhos::ProductEpetraVector(
294  x_map, stoch_row_map, sg_comm));
295  for (int block=0; block<stoch_length; ++block) {
296  for (int i=0; i<fem_length; ++i)
297  (*sg_y_commuted)[i][block] = (*sg_y)[block][i];
298  }
299 
300  typedef typename Kokkos::ViewTraits<value_type,Device,void,void>::host_mirror_space host_device;
301  typedef Stokhos::CrsProductTensor< value_type , host_device > tensor_type;
302  typedef Stokhos::StochasticProductTensor< value_type , tensor_type ,
303  host_device > stoch_tensor_type ;
304 
305  stoch_tensor_type tensor =
306  Stokhos::create_stochastic_product_tensor< tensor_type >( *basis, *Cijk );
307  stoch_length_aligned = tensor.aligned_dimension();
308 
309  perm.resize(stoch_length);
310  inv_perm.resize(stoch_length);
311  for (int i=0; i<stoch_length; ++i) {
312  Stokhos::MultiIndex<int> term(d);
313  for (int j=0; j<d; ++j)
314  term[j] = tensor.bases_degree(i,j);
315  int idx = basis->index(term);
316  perm[idx] = i;
317  inv_perm[i] = idx;
318  }
319  }
320 
321  template <typename vec_type>
322  bool test_original(const std::vector<vec_type>& y,
323  Teuchos::FancyOStream& out) const {
324  bool success = true;
325  for (int block=0; block<stoch_length; ++block) {
326  for (int i=0; i<fem_length; ++i) {
327  double diff = std::abs( (*sg_y)[block][i] - y[block][i] );
328  double tol = rel_tol*std::abs((*sg_y)[block][i]) + abs_tol;
329  bool s = diff < tol;
330  if (!s) {
331  out << "y_expected[" << block << "][" << i << "] - "
332  << "y[" << block << "][" << i << "] = " << (*sg_y)[block][i]
333  << " - " << y[block][i] << " == "
334  << diff << " < " << tol << " : failed"
335  << std::endl;
336  }
337  success = success && s;
338  }
339  }
340 
341  return success;
342  }
343 
344  template <typename multi_vec_type>
345  bool test_original(const multi_vec_type& y,
346  Teuchos::FancyOStream& out) const {
347  bool success = true;
348  for (int block=0; block<stoch_length; ++block) {
349  for (int i=0; i<fem_length; ++i) {
350  double diff = std::abs( (*sg_y)[block][i] - y(i,block) );
351  double tol = rel_tol*std::abs((*sg_y)[block][i]) + abs_tol;
352  bool s = diff < tol;
353  if (!s) {
354  out << "y_expected[" << block << "][" << i << "] - "
355  << "y(" << i << "," << block << ") = " << (*sg_y)[block][i]
356  << " - " << y(i,block) << " == "
357  << diff << " < " << tol << " : failed"
358  << std::endl;
359  }
360  success = success && s;
361  }
362  }
363 
364  return success;
365  }
366 
367  template <typename vec_type>
368  bool test_commuted_perm(const vec_type& y,
369  Teuchos::FancyOStream& out) const {
370  bool success = true;
371  for (int block=0; block<stoch_length; ++block) {
372  int b = perm[block];
373  for (int i=0; i<fem_length; ++i) {
374  double diff = std::abs( (*sg_y)[block][i] - y(b,i) );
375  double tol = rel_tol*std::abs((*sg_y)[block][i]) + abs_tol;
376  bool s = diff < tol;
377  if (!s) {
378  out << "y_expected[" << block << "][" << i << "] - "
379  << "y(" << b << "," << i << ") = " << (*sg_y)[block][i]
380  << " - " << y(b,i) << " == "
381  << diff << " < " << tol << " : failed"
382  << std::endl;
383  }
384  success = success && s;
385  }
386  }
387 
388  return success;
389  }
390 
391  template <typename vec_type>
392  bool test_commuted(const vec_type& y,
393  Teuchos::FancyOStream& out) const {
394  bool success = true;
395  for (int block=0; block<stoch_length; ++block) {
396  int b = block;
397  for (int i=0; i<fem_length; ++i) {
398  double diff = std::abs( (*sg_y)[block][i] - y(b,i) );
399  double tol = rel_tol*std::abs((*sg_y)[block][i]) + abs_tol;
400  bool s = diff < tol;
401  if (!s) {
402  out << "y_expected[" << block << "][" << i << "] - "
403  << "y(" << b << "," << i << ") = " << (*sg_y)[block][i] << " - "
404  << y(b,i) << " == "
405  << diff << " < " << tol << " : failed"
406  << std::endl;
407  }
408  success = success && s;
409  }
410  }
411 
412  return success;
413  }
414 
415  template <typename vec_type>
416  bool test_commuted_flat(const vec_type& y,
417  Teuchos::FancyOStream& out) const {
418  bool success = true;
419  for (int block=0; block<stoch_length; ++block) {
420  int b = block;
421  for (int i=0; i<fem_length; ++i) {
422  double diff = std::abs( (*sg_y)[block][i] - y(i*stoch_length+b) );
423  double tol = rel_tol*std::abs((*sg_y)[block][i]) + abs_tol;
424  bool s = diff < tol;
425  if (!s) {
426  out << "y_expected[" << block << "][" << i << "] - "
427  << "y(" << b << "," << i << ") == "
428  << diff << " < " << tol << " : failed"
429  << std::endl;
430  }
431  success = success && s;
432  }
433  }
434 
435  return success;
436  }
437 
438  template <typename vec_type>
439  bool test_original_flat(const vec_type& y,
440  Teuchos::FancyOStream& out) const {
441  bool success = true;
442  for (int block=0; block<stoch_length; ++block) {
443  int b = block;
444  for (int i=0; i<fem_length; ++i) {
445  double diff = std::abs( (*sg_y)[block][i] - y(b*fem_length+i) );
446  double tol = rel_tol*std::abs((*sg_y)[block][i]) + abs_tol;
447  bool s = diff < tol;
448  if (!s) {
449  out << "y_expected[" << block << "][" << i << "] - "
450  << "y(" << b << "," << i << ") == "
451  << diff << " < " << tol << " : failed"
452  << std::endl;
453  }
454  success = success && s;
455  }
456  }
457 
458  return success;
459  }
460 
461 };
462 
463 template <typename value_type, typename Device, typename SparseMatOps>
464 bool test_crs_matrix_free(const UnitTestSetup<Device>& setup,
465  Teuchos::FancyOStream& out) {
466  typedef Stokhos::CrsMatrix<value_type,Device> matrix_type;
467  typedef typename matrix_type::values_type matrix_values_type;
468  typedef typename matrix_type::graph_type matrix_graph_type;
469  typedef Kokkos::View<value_type*,Device> vec_type;
470 
471  //------------------------------
472 
473  std::vector<matrix_type> matrix( setup.stoch_length ) ;
474  std::vector<vec_type> x( setup.stoch_length ) ;
475  std::vector<vec_type> y( setup.stoch_length ) ;
476  std::vector<vec_type> tmp( setup.stoch_length ) ;
477 
478  for (int block=0; block<setup.stoch_length; ++block) {
479  matrix[block].graph = Kokkos::create_staticcrsgraph<matrix_graph_type>(
480  std::string("testing") , setup.fem_graph );
481 
482  matrix[block].values =
483  matrix_values_type( "matrix" , setup.fem_graph_length );
484 
485  x[block] = vec_type( "x" , setup.fem_length );
486  y[block] = vec_type( "y" , setup.fem_length );
487  tmp[block] = vec_type( "tmp" , setup.fem_length );
488 
489  typename matrix_values_type::HostMirror hM =
490  Kokkos::create_mirror( matrix[block].values );
491 
492  for ( int iRowFEM = 0 , iEntryFEM = 0 ; iRowFEM < setup.fem_length ; ++iRowFEM ) {
493  for ( size_t iRowEntryFEM = 0 ; iRowEntryFEM < setup.fem_graph[iRowFEM].size() ; ++iRowEntryFEM , ++iEntryFEM ) {
494  const int iColFEM = setup.fem_graph[iRowFEM][iRowEntryFEM] ;
495 
496  hM(iEntryFEM) = generate_matrix_coefficient<value_type>(
497  setup.fem_length , setup.stoch_length , iRowFEM , iColFEM , block );
498  }
499  }
500 
501  Kokkos::deep_copy( matrix[block].values , hM );
502 
503  typename vec_type::HostMirror hx =
504  Kokkos::create_mirror( x[block] );
505  typename vec_type::HostMirror hy =
506  Kokkos::create_mirror( y[block] );
507 
508  for ( int i = 0 ; i < setup.fem_length ; ++i ) {
509  hx(i) = generate_vector_coefficient<value_type>(
510  setup.fem_length , setup.stoch_length , i , block );
511  hy(i) = 0.0;
512  }
513 
514  Kokkos::deep_copy( x[block] , hx );
515  Kokkos::deep_copy( y[block] , hy );
516  }
517 
518  // Original matrix-free multiply algorithm using a block apply
519  SparseMatOps smo;
520  typename UnitTestSetup<Device>::Cijk_type::k_iterator k_begin =
521  setup.Cijk->k_begin();
522  typename UnitTestSetup<Device>::Cijk_type::k_iterator k_end =
523  setup.Cijk->k_end();
524  for (typename UnitTestSetup<Device>::Cijk_type::k_iterator k_it=k_begin;
525  k_it!=k_end; ++k_it) {
526  int nj = setup.Cijk->num_j(k_it);
527  if (nj > 0) {
528  int k = index(k_it);
529  typename UnitTestSetup<Device>::Cijk_type::kj_iterator j_begin =
530  setup.Cijk->j_begin(k_it);
531  typename UnitTestSetup<Device>::Cijk_type::kj_iterator j_end =
532  setup.Cijk->j_end(k_it);
533  std::vector<vec_type> xx(nj), yy(nj);
534  int jdx = 0;
535  for (typename UnitTestSetup<Device>::Cijk_type::kj_iterator j_it = j_begin;
536  j_it != j_end;
537  ++j_it) {
538  int j = index(j_it);
539  xx[jdx] = x[j];
540  yy[jdx] = tmp[j];
541  jdx++;
542  }
543  Stokhos::multiply( matrix[k] , xx , yy, smo );
544  jdx = 0;
545  for (typename UnitTestSetup<Device>::Cijk_type::kj_iterator j_it = j_begin;
546  j_it != j_end; ++j_it) {
547  typename UnitTestSetup<Device>::Cijk_type::kji_iterator i_begin =
548  setup.Cijk->i_begin(j_it);
549  typename UnitTestSetup<Device>::Cijk_type::kji_iterator i_end =
550  setup.Cijk->i_end(j_it);
551  for (typename UnitTestSetup<Device>::Cijk_type::kji_iterator i_it = i_begin;
552  i_it != i_end;
553  ++i_it) {
554  int i = index(i_it);
555  value_type c = value(i_it);
556  Stokhos::update( value_type(1.0) , y[i] , c , yy[jdx] );
557  }
558  jdx++;
559  }
560  }
561  }
562 
563  std::vector<typename vec_type::HostMirror> hy(setup.stoch_length);
564  for (int i=0; i<setup.stoch_length; ++i) {
565  hy[i] = Kokkos::create_mirror( y[i] );
566  Kokkos::deep_copy( hy[i] , y[i] );
567  }
568  bool success = setup.test_original(hy, out);
569 
570  return success;
571 }
572 
573 template <typename value_type, typename Device, typename SparseMatOps>
574 bool test_crs_matrix_free_view(const UnitTestSetup<Device>& setup,
575  Teuchos::FancyOStream& out) {
576  typedef Stokhos::CrsMatrix<value_type,Device> matrix_type;
577  typedef typename matrix_type::values_type matrix_values_type;
578  typedef typename matrix_type::graph_type matrix_graph_type;
579  typedef Kokkos::View<value_type*, Kokkos::LayoutLeft, Device, Kokkos::MemoryUnmanaged> vec_type;
580  typedef Kokkos::View<value_type**, Kokkos::LayoutLeft, Device> multi_vec_type;
581 
582  //------------------------------
583 
584  std::vector<matrix_type> matrix( setup.stoch_length ) ;
585  multi_vec_type x( "x", setup.fem_length, setup.stoch_length ) ;
586  multi_vec_type y( "y", setup.fem_length, setup.stoch_length ) ;
587  multi_vec_type tmp_x( "tmp_x", setup.fem_length, setup.stoch_length ) ;
588  multi_vec_type tmp_y( "tmp_y", setup.fem_length, setup.stoch_length ) ;
589 
590  typename multi_vec_type::HostMirror hx = Kokkos::create_mirror( x );
591  typename multi_vec_type::HostMirror hy = Kokkos::create_mirror( y );
592 
593  for (int block=0; block<setup.stoch_length; ++block) {
594  matrix[block].graph = Kokkos::create_staticcrsgraph<matrix_graph_type>(
595  std::string("testing") , setup.fem_graph );
596 
597  matrix[block].values =
598  matrix_values_type( "matrix" , setup.fem_graph_length );
599 
600  typename matrix_values_type::HostMirror hM =
601  Kokkos::create_mirror( matrix[block].values );
602 
603  for ( int iRowFEM = 0 , iEntryFEM = 0 ; iRowFEM < setup.fem_length ; ++iRowFEM ) {
604  for ( size_t iRowEntryFEM = 0 ; iRowEntryFEM < setup.fem_graph[iRowFEM].size() ; ++iRowEntryFEM , ++iEntryFEM ) {
605  const int iColFEM = setup.fem_graph[iRowFEM][iRowEntryFEM] ;
606 
607  hM(iEntryFEM) = generate_matrix_coefficient<value_type>(
608  setup.fem_length , setup.stoch_length , iRowFEM , iColFEM , block );
609  }
610  }
611 
612  Kokkos::deep_copy( matrix[block].values , hM );
613 
614  for ( int i = 0 ; i < setup.fem_length ; ++i ) {
615  hx(i, block) = generate_vector_coefficient<value_type>(
616  setup.fem_length , setup.stoch_length , i , block );
617  hy(i, block) = 0.0;
618  }
619 
620  }
621 
622  Kokkos::deep_copy( x , hx );
623  Kokkos::deep_copy( y , hy );
624 
625  // Original matrix-free multiply algorithm using a block apply
626  typedef typename UnitTestSetup<Device>::Cijk_type::k_iterator k_iterator;
627  typedef typename UnitTestSetup<Device>::Cijk_type::kj_iterator kj_iterator;
628  typedef typename UnitTestSetup<Device>::Cijk_type::kji_iterator kji_iterator;
629  SparseMatOps smo;
630  k_iterator k_begin = setup.Cijk->k_begin();
631  k_iterator k_end = setup.Cijk->k_end();
632  for (k_iterator k_it=k_begin; k_it!=k_end; ++k_it) {
633  unsigned nj = setup.Cijk->num_j(k_it);
634  if (nj > 0) {
635  int k = index(k_it);
636  kj_iterator j_begin = setup.Cijk->j_begin(k_it);
637  kj_iterator j_end = setup.Cijk->j_end(k_it);
638  std::vector<int> j_indices(nj);
639  unsigned jdx = 0;
640  for (kj_iterator j_it = j_begin; j_it != j_end; ++j_it) {
641  int j = index(j_it);
642  vec_type xx = Kokkos::subview( x, Kokkos::ALL(), j );
643  vec_type tt = Kokkos::subview( tmp_x, Kokkos::ALL(), jdx++ );
644  Kokkos::deep_copy(tt, xx);
645  }
646  multi_vec_type tmp_x_view =
647  Kokkos::subview( tmp_x, Kokkos::ALL(),
648  std::make_pair(0u,nj));
649  multi_vec_type tmp_y_view =
650  Kokkos::subview( tmp_y, Kokkos::ALL(),
651  std::make_pair(0u,nj));
652  Stokhos::multiply( matrix[k] , tmp_x_view , tmp_y_view, smo );
653  jdx = 0;
654  for (kj_iterator j_it = j_begin; j_it != j_end; ++j_it) {
655  vec_type tmp_y_view =
656  Kokkos::subview( tmp_y, Kokkos::ALL(), jdx++ );
657  kji_iterator i_begin = setup.Cijk->i_begin(j_it);
658  kji_iterator i_end = setup.Cijk->i_end(j_it);
659  for (kji_iterator i_it = i_begin; i_it != i_end; ++i_it) {
660  int i = index(i_it);
661  value_type c = value(i_it);
662  vec_type y_view = Kokkos::subview( y, Kokkos::ALL(), i );
663  Stokhos::update( value_type(1.0) , y_view , c , tmp_y_view );
664  }
665  }
666  }
667  }
668 
669  Kokkos::deep_copy( hy , y );
670  bool success = setup.test_original(hy, out);
671 
672  return success;
673 }
674 
675 #ifdef HAVE_STOKHOS_KOKKOSLINALG
676 
677 template <typename value_type, typename Device>
678 bool test_crs_matrix_free_kokkos(const UnitTestSetup<Device>& setup,
679  Teuchos::FancyOStream& out) {
680  typedef int ordinal_type;
681  typedef KokkosSparse::CrsMatrix<value_type,ordinal_type,Device> matrix_type;
682  typedef typename matrix_type::values_type matrix_values_type;
683  typedef typename matrix_type::StaticCrsGraphType matrix_graph_type;
684  typedef Kokkos::View<value_type*, Kokkos::LayoutLeft, Device, Kokkos::MemoryUnmanaged> vec_type;
685  typedef Kokkos::View<value_type**, Kokkos::LayoutLeft, Device> multi_vec_type;
686 
687  //------------------------------
688 
689  std::vector<matrix_type> matrix( setup.stoch_length ) ;
690  multi_vec_type x( "x", setup.fem_length, setup.stoch_length ) ;
691  multi_vec_type y( "y", setup.fem_length, setup.stoch_length ) ;
692  multi_vec_type tmp_x( "tmp_x", setup.fem_length, setup.stoch_length ) ;
693  multi_vec_type tmp_y( "tmp_y", setup.fem_length, setup.stoch_length ) ;
694 
695  typename multi_vec_type::HostMirror hx = Kokkos::create_mirror( x );
696  typename multi_vec_type::HostMirror hy = Kokkos::create_mirror( y );
697 
698  for (int block=0; block<setup.stoch_length; ++block) {
699  matrix_graph_type matrix_graph =
700  Kokkos::create_staticcrsgraph<matrix_graph_type>(
701  std::string("test crs graph"), setup.fem_graph);
702  matrix_values_type matrix_values =
703  matrix_values_type( "matrix" , setup.fem_graph_length );
704  typename matrix_values_type::HostMirror hM =
705  Kokkos::create_mirror( matrix_values );
706 
707  for ( int iRowFEM = 0 , iEntryFEM = 0 ; iRowFEM < setup.fem_length ; ++iRowFEM ) {
708  for ( size_t iRowEntryFEM = 0 ; iRowEntryFEM < setup.fem_graph[iRowFEM].size() ; ++iRowEntryFEM , ++iEntryFEM ) {
709  const int iColFEM = setup.fem_graph[iRowFEM][iRowEntryFEM] ;
710 
711  hM(iEntryFEM) = generate_matrix_coefficient<value_type>(
712  setup.fem_length , setup.stoch_length , iRowFEM , iColFEM , block );
713  }
714  }
715 
716  Kokkos::deep_copy( matrix_values , hM );
717  matrix[block] = matrix_type("matrix", setup.fem_length, matrix_values,
718  matrix_graph);
719 
720  for ( int i = 0 ; i < setup.fem_length ; ++i ) {
721  hx(i, block) = generate_vector_coefficient<value_type>(
722  setup.fem_length , setup.stoch_length , i , block );
723  hy(i, block) = 0.0;
724  }
725 
726  }
727 
728  Kokkos::deep_copy( x , hx );
729  Kokkos::deep_copy( y , hy );
730 
731  // Original matrix-free multiply algorithm using a block apply
732  typedef typename UnitTestSetup<Device>::Cijk_type::k_iterator k_iterator;
733  typedef typename UnitTestSetup<Device>::Cijk_type::kj_iterator kj_iterator;
734  typedef typename UnitTestSetup<Device>::Cijk_type::kji_iterator kji_iterator;
735  k_iterator k_begin = setup.Cijk->k_begin();
736  k_iterator k_end = setup.Cijk->k_end();
737  for (k_iterator k_it=k_begin; k_it!=k_end; ++k_it) {
738  int nj = setup.Cijk->num_j(k_it);
739  if (nj > 0) {
740  int k = index(k_it);
741  kj_iterator j_begin = setup.Cijk->j_begin(k_it);
742  kj_iterator j_end = setup.Cijk->j_end(k_it);
743  unsigned jdx = 0;
744  for (kj_iterator j_it = j_begin; j_it != j_end; ++j_it) {
745  int j = index(j_it);
746  vec_type xx = Kokkos::subview( x, Kokkos::ALL(), j );
747  vec_type tt = Kokkos::subview( tmp_x, Kokkos::ALL(), jdx++ );
748  Kokkos::deep_copy(tt, xx);
749  }
750  multi_vec_type tmp_x_view =
751  Kokkos::subview( tmp_x, Kokkos::ALL(),
752  std::make_pair(0u,jdx));
753  multi_vec_type tmp_y_view =
754  Kokkos::subview( tmp_y, Kokkos::ALL(),
755  std::make_pair(0u,jdx));
756  KokkosSparse::spmv( "N", value_type(1.0), matrix[k] , tmp_x_view , value_type(0.0) , tmp_y_view );
757  jdx = 0;
758  for (kj_iterator j_it = j_begin; j_it != j_end; ++j_it) {
759  vec_type tmp_y_view =
760  Kokkos::subview( tmp_y, Kokkos::ALL(), jdx++ );
761  kji_iterator i_begin = setup.Cijk->i_begin(j_it);
762  kji_iterator i_end = setup.Cijk->i_end(j_it);
763  for (kji_iterator i_it = i_begin; i_it != i_end; ++i_it) {
764  int i = index(i_it);
765  value_type c = value(i_it);
766  vec_type y_view = Kokkos::subview( y, Kokkos::ALL(), i );
767  //Stokhos::update( value_type(1.0) , y_view , c , tmp_y_view );
768  KokkosBlas::update(c, tmp_y_view, value_type(1.0), y_view, value_type(0.0), y_view);
769  }
770  }
771  }
772  }
773 
774  Kokkos::deep_copy( hy , y );
775  bool success = setup.test_original(hy, out);
776 
777  return success;
778 }
779 
780 #endif
781 
782 template< typename ScalarType , class Device >
783 bool
784 test_crs_dense_block(const UnitTestSetup<Device>& setup,
785  Teuchos::FancyOStream& out)
786 {
787  typedef ScalarType value_type ;
788  typedef Kokkos::View< value_type**, Kokkos::LayoutLeft , Device > block_vector_type ;
789 
790  //------------------------------
791 
792  typedef Stokhos::BlockCrsMatrix< Stokhos::SymmetricDiagonalSpec< Device > , value_type , Device > matrix_type ;
793 
794  typedef typename matrix_type::graph_type graph_type ;
795 
796  // Convert sparse Cijk to dense for faster assembly
797  typedef typename UnitTestSetup<Device>::Cijk_type::k_iterator k_iterator;
798  typedef typename UnitTestSetup<Device>::Cijk_type::kj_iterator kj_iterator;
799  typedef typename UnitTestSetup<Device>::Cijk_type::kji_iterator kji_iterator;
800  Stokhos::Dense3Tensor<int,double> dense_cijk(setup.stoch_length);
801  for (k_iterator k_it=setup.Cijk->k_begin();
802  k_it!=setup.Cijk->k_end(); ++k_it) {
803  int k = index(k_it);
804  for (kj_iterator j_it = setup.Cijk->j_begin(k_it);
805  j_it != setup.Cijk->j_end(k_it); ++j_it) {
806  int j = index(j_it);
807  for (kji_iterator i_it = setup.Cijk->i_begin(j_it);
808  i_it != setup.Cijk->i_end(j_it); ++i_it) {
809  int i = index(i_it);
810  double c = value(i_it);
811  dense_cijk(i,j,k) = c;
812  }
813  }
814  }
815 
816  //------------------------------
817  // Generate input multivector:
818 
819  block_vector_type x =
820  block_vector_type( "x" , setup.stoch_length , setup.fem_length );
821  block_vector_type y =
822  block_vector_type( "y" , setup.stoch_length , setup.fem_length );
823 
824  typename block_vector_type::HostMirror hx = Kokkos::create_mirror( x );
825 
826  for ( int iColFEM = 0 ; iColFEM < setup.fem_length ; ++iColFEM ) {
827  for ( int iColStoch = 0 ; iColStoch < setup.stoch_length ; ++iColStoch ) {
828  hx(iColStoch,iColFEM) =
829  generate_vector_coefficient<ScalarType>(
830  setup.fem_length , setup.stoch_length , iColFEM , iColStoch );
831  }
832  }
833 
834  Kokkos::deep_copy( x , hx );
835 
836  //------------------------------
837  // Generate CRS matrix of blocks with symmetric diagonal storage
838 
839  matrix_type matrix ;
840 
841  matrix.block =
842  Stokhos::SymmetricDiagonalSpec< Device >( setup.stoch_length );
843  matrix.graph = Kokkos::create_staticcrsgraph<graph_type>(
844  std::string("test crs graph") , setup.fem_graph );
845  matrix.values = block_vector_type(
846  "matrix" , matrix.block.matrix_size() , setup.fem_graph_length );
847 
848  {
849  typename block_vector_type::HostMirror hM =
850  Kokkos::create_mirror( matrix.values );
851 
852  for ( int iRowStoch = 0 ; iRowStoch < setup.stoch_length ; ++iRowStoch ) {
853  for ( int iRowFEM = 0 , iEntryFEM = 0 ; iRowFEM < setup.fem_length ; ++iRowFEM ) {
854 
855  for ( size_t iRowEntryFEM = 0 ; iRowEntryFEM < setup.fem_graph[iRowFEM].size() ; ++iRowEntryFEM , ++iEntryFEM ) {
856  const int iColFEM = setup.fem_graph[iRowFEM][iRowEntryFEM];
857 
858  for ( int iColStoch = 0 ; iColStoch < setup.stoch_length ; ++iColStoch ) {
859 
860  const size_t offset =
861  matrix.block.matrix_offset( iRowStoch , iColStoch );
862 
863  ScalarType value = 0 ;
864 
865  for ( int k = 0 ; k < setup.stoch_length ; ++k ) {
866  value += dense_cijk( iRowStoch , iColStoch , k ) *
867  generate_matrix_coefficient<ScalarType>(
868  setup.fem_length , setup.stoch_length , iRowFEM , iColFEM , k );
869  }
870 
871  hM( offset , iEntryFEM ) = value ;
872  }
873  }
874 
875  }
876  }
877 
878  Kokkos::deep_copy( matrix.values , hM );
879  }
880 
881  //------------------------------
882 
883  Stokhos::multiply( matrix , x , y );
884 
885  typename block_vector_type::HostMirror hy = Kokkos::create_mirror( y );
886  Kokkos::deep_copy( hy , y );
887 
888  bool success = setup.test_commuted(hy, out);
889 
890  return success;
891 }
892 
893 template< typename ScalarType , class Device >
894 bool
895 test_crs_flat_commuted(const UnitTestSetup<Device>& setup,
896  Teuchos::FancyOStream& out)
897 {
898  typedef ScalarType value_type ;
899  typedef Kokkos::View< value_type* , Device > vector_type ;
900 
901  //------------------------------
902 
903  typedef Stokhos::CrsMatrix<value_type,Device> matrix_type;
904  typedef typename matrix_type::values_type matrix_values_type;
905  typedef typename matrix_type::graph_type matrix_graph_type;
906 
907  // Build stochastic graph
908  std::vector< std::vector< int > > stoch_graph( setup.stoch_length );
909  Teuchos::RCP<Epetra_CrsGraph> cijk_graph = Stokhos::sparse3Tensor2CrsGraph(
910  *setup.basis, *setup.Cijk, *setup.globalComm);
911  for ( int i = 0 ; i < setup.stoch_length ; ++i ) {
912  int len = cijk_graph->NumGlobalIndices(i);
913  stoch_graph[i].resize(len);
914  int len2;
915  cijk_graph->ExtractGlobalRowCopy(i, len, len2, &stoch_graph[i][0]);
916  }
917 
918  // Convert sparse Cijk to dense for faster assembly in debug builds
919  typedef typename UnitTestSetup<Device>::Cijk_type::k_iterator k_iterator;
920  typedef typename UnitTestSetup<Device>::Cijk_type::kj_iterator kj_iterator;
921  typedef typename UnitTestSetup<Device>::Cijk_type::kji_iterator kji_iterator;
922  Stokhos::Dense3Tensor<int,double> dense_cijk(setup.stoch_length);
923  for (k_iterator k_it=setup.Cijk->k_begin();
924  k_it!=setup.Cijk->k_end(); ++k_it) {
925  int k = index(k_it);
926  for (kj_iterator j_it = setup.Cijk->j_begin(k_it);
927  j_it != setup.Cijk->j_end(k_it); ++j_it) {
928  int j = index(j_it);
929  for (kji_iterator i_it = setup.Cijk->i_begin(j_it);
930  i_it != setup.Cijk->i_end(j_it); ++i_it) {
931  int i = index(i_it);
932  double c = value(i_it);
933  dense_cijk(i,j,k) = c;
934  }
935  }
936  }
937 
938  //------------------------------
939  // Generate flattened graph with FEM outer and stochastic inner
940 
941  const int flat_length = setup.fem_length * setup.stoch_length ;
942 
943  std::vector< std::vector<int> > flat_graph( flat_length );
944 
945  for ( int iOuterRow = 0 ; iOuterRow < setup.fem_length ; ++iOuterRow ) {
946 
947  const size_t iOuterRowNZ = setup.fem_graph[iOuterRow].size();
948 
949  for ( int iInnerRow = 0 ; iInnerRow < setup.stoch_length ; ++iInnerRow ) {
950 
951  const size_t iInnerRowNZ = stoch_graph[ iInnerRow ].size(); ;
952  const int iFlatRowNZ = iOuterRowNZ * iInnerRowNZ ;
953  const int iFlatRow = iInnerRow + iOuterRow * setup.stoch_length ;
954 
955  flat_graph[iFlatRow].resize( iFlatRowNZ );
956 
957  int iFlatEntry = 0 ;
958 
959  for ( size_t iOuterEntry = 0 ; iOuterEntry < iOuterRowNZ ; ++iOuterEntry ) {
960 
961  const int iOuterCol = setup.fem_graph[iOuterRow][iOuterEntry];
962 
963  for ( size_t iInnerEntry = 0 ; iInnerEntry < iInnerRowNZ ; ++iInnerEntry ) {
964 
965  const int iInnerCol = stoch_graph[iInnerRow][iInnerEntry] ;
966  const int iFlatColumn = iInnerCol + iOuterCol * setup.stoch_length ;
967 
968  flat_graph[iFlatRow][iFlatEntry] = iFlatColumn ;
969 
970  ++iFlatEntry ;
971  }
972  }
973  }
974  }
975 
976  //------------------------------
977 
978  vector_type x = vector_type( "x" , flat_length );
979  vector_type y = vector_type( "y" , flat_length );
980 
981  typename vector_type::HostMirror hx = Kokkos::create_mirror( x );
982 
983  for ( int iColFEM = 0 ; iColFEM < setup.fem_length ; ++iColFEM ) {
984  for ( int iColStoch = 0 ; iColStoch < setup.stoch_length ; ++iColStoch ) {
985  hx(iColStoch + iColFEM*setup.stoch_length) =
986  generate_vector_coefficient<ScalarType>(
987  setup.fem_length , setup.stoch_length , iColFEM , iColStoch );
988  }
989  }
990 
991  Kokkos::deep_copy( x , hx );
992 
993  //------------------------------
994 
995  matrix_type matrix ;
996 
997  matrix.graph = Kokkos::create_staticcrsgraph<matrix_graph_type>(
998  std::string("testing") , flat_graph );
999 
1000  const size_t flat_graph_length = matrix.graph.entries.extent(0);
1001 
1002  matrix.values = matrix_values_type( "matrix" , flat_graph_length );
1003  {
1004  typename matrix_values_type::HostMirror hM =
1005  Kokkos::create_mirror( matrix.values );
1006 
1007  for ( int iRow = 0 , iEntry = 0 ; iRow < flat_length ; ++iRow ) {
1008  const int iRowFEM = iRow / setup.stoch_length ;
1009  const int iRowStoch = iRow % setup.stoch_length ;
1010 
1011  for ( size_t iRowEntry = 0 ; iRowEntry < flat_graph[ iRow ].size() ; ++iRowEntry , ++iEntry ) {
1012  const int iCol = flat_graph[ iRow ][ iRowEntry ];
1013  const int iColFEM = iCol / setup.stoch_length ;
1014  const int iColStoch = iCol % setup.stoch_length ;
1015 
1016  double value = 0 ;
1017  for ( int k = 0 ; k < setup.stoch_length ; ++k ) {
1018  const double A_fem_k =
1019  generate_matrix_coefficient<ScalarType>(
1020  setup.fem_length , setup.stoch_length , iRowFEM, iColFEM, k );
1021  value += dense_cijk(iRowStoch,iColStoch,k) * A_fem_k ;
1022  }
1023  hM( iEntry ) = value ;
1024  }
1025  }
1026 
1027  Kokkos::deep_copy( matrix.values , hM );
1028  }
1029 
1030  //------------------------------
1031 
1032 
1033  Stokhos::multiply( matrix , x , y );
1034 
1035  typename vector_type::HostMirror hy = Kokkos::create_mirror( y );
1036  Kokkos::deep_copy( hy , y );
1037 
1038  bool success = setup.test_commuted_flat(hy, out);
1039  return success;
1040 }
1041 
1042 template< typename ScalarType , class Device >
1043 bool
1044 test_crs_flat_original(const UnitTestSetup<Device>& setup,
1045  Teuchos::FancyOStream& out)
1046 {
1047  typedef ScalarType value_type ;
1048  typedef Kokkos::View< value_type* , Device > vector_type ;
1049 
1050  //------------------------------
1051 
1052  typedef Stokhos::CrsMatrix<value_type,Device> matrix_type;
1053  typedef typename matrix_type::values_type matrix_values_type;
1054  typedef typename matrix_type::graph_type matrix_graph_type;
1055 
1056  // Build stochastic graph
1057  std::vector< std::vector< int > > stoch_graph( setup.stoch_length );
1058  Teuchos::RCP<Epetra_CrsGraph> cijk_graph = Stokhos::sparse3Tensor2CrsGraph(
1059  *setup.basis, *setup.Cijk, *setup.globalComm);
1060  for ( int i = 0 ; i < setup.stoch_length ; ++i ) {
1061  int len = cijk_graph->NumGlobalIndices(i);
1062  stoch_graph[i].resize(len);
1063  int len2;
1064  cijk_graph->ExtractGlobalRowCopy(i, len, len2, &stoch_graph[i][0]);
1065  }
1066 
1067  // Convert sparse Cijk to dense for faster assembly in debug builds
1068  typedef typename UnitTestSetup<Device>::Cijk_type::k_iterator k_iterator;
1069  typedef typename UnitTestSetup<Device>::Cijk_type::kj_iterator kj_iterator;
1070  typedef typename UnitTestSetup<Device>::Cijk_type::kji_iterator kji_iterator;
1071  Stokhos::Dense3Tensor<int,double> dense_cijk(setup.stoch_length);
1072  for (k_iterator k_it=setup.Cijk->k_begin();
1073  k_it!=setup.Cijk->k_end(); ++k_it) {
1074  int k = index(k_it);
1075  for (kj_iterator j_it = setup.Cijk->j_begin(k_it);
1076  j_it != setup.Cijk->j_end(k_it); ++j_it) {
1077  int j = index(j_it);
1078  for (kji_iterator i_it = setup.Cijk->i_begin(j_it);
1079  i_it != setup.Cijk->i_end(j_it); ++i_it) {
1080  int i = index(i_it);
1081  double c = value(i_it);
1082  dense_cijk(i,j,k) = c;
1083  }
1084  }
1085  }
1086 
1087  //------------------------------
1088  // Generate flattened graph with stochastic outer and FEM inner
1089 
1090  const size_t flat_length = setup.fem_length * setup.stoch_length ;
1091 
1092  std::vector< std::vector<int> > flat_graph( flat_length );
1093 
1094  for ( int iOuterRow = 0 ; iOuterRow < setup.stoch_length ; ++iOuterRow ) {
1095 
1096  const size_t iOuterRowNZ = stoch_graph[iOuterRow].size();
1097 
1098  for ( int iInnerRow = 0 ; iInnerRow < setup.fem_length ; ++iInnerRow ) {
1099 
1100  const size_t iInnerRowNZ = setup.fem_graph[iInnerRow].size();
1101  const int iFlatRowNZ = iOuterRowNZ * iInnerRowNZ ;
1102  const int iFlatRow = iInnerRow + iOuterRow * setup.fem_length ;
1103 
1104  flat_graph[iFlatRow].resize( iFlatRowNZ );
1105 
1106  int iFlatEntry = 0 ;
1107 
1108  for ( size_t iOuterEntry = 0 ; iOuterEntry < iOuterRowNZ ; ++iOuterEntry ) {
1109 
1110  const int iOuterCol = stoch_graph[ iOuterRow ][ iOuterEntry ];
1111 
1112  for ( size_t iInnerEntry = 0 ; iInnerEntry < iInnerRowNZ ; ++iInnerEntry ) {
1113 
1114  const int iInnerCol = setup.fem_graph[ iInnerRow][iInnerEntry];
1115  const int iFlatColumn = iInnerCol + iOuterCol * setup.fem_length ;
1116 
1117  flat_graph[iFlatRow][iFlatEntry] = iFlatColumn ;
1118  ++iFlatEntry ;
1119  }
1120  }
1121  }
1122  }
1123 
1124  //------------------------------
1125 
1126  vector_type x = vector_type( "x" , flat_length );
1127  vector_type y = vector_type( "y" , flat_length );
1128 
1129  typename vector_type::HostMirror hx = Kokkos::create_mirror( x );
1130 
1131  for ( size_t iCol = 0 ; iCol < flat_length ; ++iCol ) {
1132  const int iColStoch = iCol / setup.fem_length ;
1133  const int iColFEM = iCol % setup.fem_length ;
1134 
1135  hx(iCol) = generate_vector_coefficient<ScalarType>(
1136  setup.fem_length , setup.stoch_length , iColFEM , iColStoch );
1137  }
1138 
1139  Kokkos::deep_copy( x , hx );
1140 
1141  //------------------------------
1142 
1143  matrix_type matrix ;
1144 
1145  matrix.graph = Kokkos::create_staticcrsgraph<matrix_graph_type>( std::string("testing") , flat_graph );
1146 
1147  const size_t flat_graph_length = matrix.graph.entries.extent(0);
1148 
1149  matrix.values = matrix_values_type( "matrix" , flat_graph_length );
1150  {
1151  typename matrix_values_type::HostMirror hM =
1152  Kokkos::create_mirror( matrix.values );
1153 
1154  for ( size_t iRow = 0 , iEntry = 0 ; iRow < flat_length ; ++iRow ) {
1155  const int iRowStoch = iRow / setup.fem_length ;
1156  const int iRowFEM = iRow % setup.fem_length ;
1157 
1158  for ( size_t iRowEntry = 0 ; iRowEntry < flat_graph[ iRow ].size() ; ++iRowEntry , ++iEntry ) {
1159  const int iCol = flat_graph[ iRow ][ iRowEntry ];
1160  const int iColStoch = iCol / setup.fem_length ;
1161  const int iColFEM = iCol % setup.fem_length ;
1162 
1163  double value = 0 ;
1164  for ( int k = 0 ; k < setup.stoch_length ; ++k ) {
1165  const double A_fem_k =
1166  generate_matrix_coefficient<ScalarType>(
1167  setup.fem_length , setup.stoch_length ,
1168  iRowFEM , iColFEM , k );
1169  value += dense_cijk(iRowStoch,iColStoch,k) * A_fem_k ;
1170  }
1171  hM( iEntry ) = value ;
1172 
1173  }
1174 
1175  }
1176 
1177  Kokkos::deep_copy( matrix.values , hM );
1178  }
1179 
1180  Stokhos::multiply( matrix , x , y );
1181 
1182  typename vector_type::HostMirror hy = Kokkos::create_mirror( y );
1183  Kokkos::deep_copy( hy , y );
1184 
1185  bool success = setup.test_original_flat(hy, out);
1186  return success;
1187 }
1188 
1189 template< typename ScalarType , typename TensorType, class Device >
1190 bool test_crs_product_tensor(
1191  const UnitTestSetup<Device>& setup,
1192  Teuchos::FancyOStream& out,
1193  const Teuchos::ParameterList& params = Teuchos::ParameterList()) {
1194  typedef ScalarType value_type ;
1195  typedef Kokkos::View< value_type** , Kokkos::LayoutLeft ,
1196  Device > block_vector_type ;
1197 
1198  typedef Stokhos::StochasticProductTensor< value_type , TensorType , Device > tensor_type ;
1199 
1200  typedef Stokhos::BlockCrsMatrix< tensor_type , value_type , Device > matrix_type ;
1201  typedef typename matrix_type::graph_type graph_type ;
1202 
1203  //------------------------------
1204  // Generate input multivector:
1205 
1206  block_vector_type x =
1207  block_vector_type( "x" , setup.stoch_length_aligned , setup.fem_length );
1208  block_vector_type y =
1209  block_vector_type( "y" , setup.stoch_length_aligned , setup.fem_length );
1210 
1211  typename block_vector_type::HostMirror hx =
1212  Kokkos::create_mirror( x );
1213 
1214  for ( int iColFEM = 0 ; iColFEM < setup.fem_length ; ++iColFEM ) {
1215  for ( int iColStoch = 0 ; iColStoch < setup.stoch_length ; ++iColStoch ) {
1216  hx(setup.perm[iColStoch],iColFEM) =
1217  generate_vector_coefficient<ScalarType>(
1218  setup.fem_length , setup.stoch_length , iColFEM , iColStoch );
1219  }
1220  }
1221 
1222  Kokkos::deep_copy( x , hx );
1223 
1224  //------------------------------
1225 
1226  matrix_type matrix ;
1227 
1228  matrix.block =
1229  Stokhos::create_stochastic_product_tensor< TensorType >( *setup.basis,
1230  *setup.Cijk,
1231  params);
1232 
1233  matrix.graph = Kokkos::create_staticcrsgraph<graph_type>(
1234  std::string("test crs graph") , setup.fem_graph );
1235 
1236  matrix.values = block_vector_type(
1237  "matrix" , setup.stoch_length_aligned , setup.fem_graph_length );
1238 
1239  typename block_vector_type::HostMirror hM =
1240  Kokkos::create_mirror( matrix.values );
1241 
1242  for ( int iRowFEM = 0 , iEntryFEM = 0 ; iRowFEM < setup.fem_length ; ++iRowFEM ) {
1243  for ( size_t iRowEntryFEM = 0 ; iRowEntryFEM < setup.fem_graph[iRowFEM].size() ; ++iRowEntryFEM , ++iEntryFEM ) {
1244  const int iColFEM = setup.fem_graph[iRowFEM][iRowEntryFEM] ;
1245 
1246  for ( int k = 0 ; k < setup.stoch_length ; ++k ) {
1247  hM(setup.perm[k],iEntryFEM) =
1248  generate_matrix_coefficient<ScalarType>(
1249  setup.fem_length , setup.stoch_length , iRowFEM , iColFEM , k );
1250  //hM(k,iEntryFEM) = 1.0;
1251  }
1252  }
1253  }
1254 
1255  Kokkos::deep_copy( matrix.values , hM );
1256 
1257  //------------------------------
1258 
1259  Stokhos::multiply( matrix , x , y );
1260 
1261  typename block_vector_type::HostMirror hy = Kokkos::create_mirror( y );
1262  Kokkos::deep_copy( hy , y );
1263 
1264  bool success = setup.test_commuted_perm(hy, out);
1265  //bool success = true;
1266  return success;
1267 }
1268 
1269 template< typename ScalarType , class Device , int BlockSize >
1270 bool test_linear_tensor(const UnitTestSetup<Device>& setup,
1271  Teuchos::FancyOStream& out,
1272  const bool symmetric) {
1273  typedef ScalarType value_type ;
1274  typedef Kokkos::View< value_type** , Kokkos::LayoutLeft ,
1275  Device > block_vector_type ;
1276  typedef Stokhos::LinearSparse3Tensor<value_type,Device,BlockSize> TensorType;
1277  typedef Stokhos::StochasticProductTensor< value_type , TensorType , Device > tensor_type ;
1278 
1279  typedef Stokhos::BlockCrsMatrix< tensor_type , value_type , Device > matrix_type ;
1280  typedef typename matrix_type::graph_type graph_type ;
1281 
1282  // Build tensor
1283  matrix_type matrix ;
1284 
1285  Teuchos::ParameterList params;
1286  params.set("Symmetric",symmetric);
1287  matrix.block =
1288  Stokhos::create_stochastic_product_tensor< TensorType >( *setup.basis,
1289  *setup.Cijk,
1290  params );
1291  int aligned_stoch_length = matrix.block.tensor().aligned_dimension();
1292 
1293  //------------------------------
1294  // Generate input multivector:
1295 
1296  block_vector_type x =
1297  block_vector_type( "x" , aligned_stoch_length , setup.fem_length );
1298  block_vector_type y =
1299  block_vector_type( "y" , aligned_stoch_length , setup.fem_length );
1300 
1301  typename block_vector_type::HostMirror hx =
1302  Kokkos::create_mirror( x );
1303 
1304  for ( int iColFEM = 0 ; iColFEM < setup.fem_length ; ++iColFEM ) {
1305  for ( int iColStoch = 0 ; iColStoch < setup.stoch_length ; ++iColStoch ) {
1306  hx(iColStoch,iColFEM) =
1307  generate_vector_coefficient<ScalarType>(
1308  setup.fem_length , setup.stoch_length , iColFEM , iColStoch );
1309  //hx(iColStoch,iColFEM) = 1.0;
1310  }
1311  }
1312 
1313  Kokkos::deep_copy( x , hx );
1314 
1315  //------------------------------
1316 
1317  matrix.graph = Kokkos::create_staticcrsgraph<graph_type>(
1318  std::string("test crs graph") , setup.fem_graph );
1319 
1320  matrix.values = block_vector_type(
1321  "matrix" , aligned_stoch_length , setup.fem_graph_length );
1322 
1323  typename block_vector_type::HostMirror hM =
1324  Kokkos::create_mirror( matrix.values );
1325 
1326  for ( int iRowFEM = 0 , iEntryFEM = 0 ; iRowFEM < setup.fem_length ; ++iRowFEM ) {
1327  for ( size_t iRowEntryFEM = 0 ; iRowEntryFEM < setup.fem_graph[iRowFEM].size() ; ++iRowEntryFEM , ++iEntryFEM ) {
1328  const int iColFEM = setup.fem_graph[iRowFEM][iRowEntryFEM] ;
1329 
1330  for ( int k = 0 ; k < setup.stoch_length ; ++k ) {
1331  hM(k,iEntryFEM) =
1332  generate_matrix_coefficient<ScalarType>(
1333  setup.fem_length , setup.stoch_length , iRowFEM , iColFEM , k );
1334  //hM(k,iEntryFEM) = 1.0;
1335  }
1336  }
1337  }
1338 
1339  Kokkos::deep_copy( matrix.values , hM );
1340 
1341  //------------------------------
1342 
1343  Stokhos::multiply( matrix , x , y );
1344 
1345  typename block_vector_type::HostMirror hy = Kokkos::create_mirror( y );
1346  Kokkos::deep_copy( hy , y );
1347 
1348  bool success = setup.test_commuted(hy, out);
1349  //bool success = true;
1350  return success;
1351 }
1352 
1353 template< typename ScalarType , class Device >
1354 bool test_lexo_block_tensor(const UnitTestSetup<Device>& setup,
1355  Teuchos::FancyOStream& out) {
1356  typedef ScalarType value_type ;
1357  typedef int ordinal_type;
1358  typedef Kokkos::View< value_type** , Kokkos::LayoutLeft ,
1359  Device > block_vector_type ;
1360  typedef Stokhos::LexicographicBlockSparse3Tensor<value_type,Device> TensorType;
1361  typedef Stokhos::StochasticProductTensor< value_type , TensorType , Device > tensor_type ;
1362 
1363  typedef Stokhos::BlockCrsMatrix< tensor_type , value_type , Device > matrix_type ;
1364  typedef typename matrix_type::graph_type graph_type ;
1365 
1366  //------------------------------
1367  // Generate input multivector:
1368 
1369  block_vector_type x =
1370  block_vector_type( "x" , setup.stoch_length , setup.fem_length );
1371  block_vector_type y =
1372  block_vector_type( "y" , setup.stoch_length , setup.fem_length );
1373 
1374  typename block_vector_type::HostMirror hx =
1375  Kokkos::create_mirror( x );
1376 
1377  for ( int iColFEM = 0 ; iColFEM < setup.fem_length ; ++iColFEM ) {
1378  for ( int iColStoch = 0 ; iColStoch < setup.stoch_length ; ++iColStoch ) {
1379  hx(iColStoch,iColFEM) =
1380  generate_vector_coefficient<ScalarType>(
1381  setup.fem_length , setup.stoch_length , iColFEM , iColStoch );
1382  }
1383  }
1384 
1385  Kokkos::deep_copy( x , hx );
1386 
1387  //------------------------------
1388 
1389  matrix_type matrix ;
1390 
1391  /*
1392  typedef UnitTestSetup<Device>::abstract_basis_type abstract_basis_type;
1393  typedef UnitTestSetup<Device>::basis_type basis_type;
1394  typedef Stokhos::LexographicLess< Stokhos::MultiIndex<int> > less_type;
1395  typedef Stokhos::TotalOrderBasis<ordinal_type,value_type,less_type> product_basis_type;
1396  Teuchos::Array< Teuchos::RCP<const abstract_basis_type> > bases(setup.d);
1397  for (int i=0; i<setup.d; i++)
1398  bases[i] = rcp(new basis_type(setup.p,true));
1399  product_basis_type basis(bases, 1e-12);
1400  */
1401  const bool symmetric = true;
1402  Teuchos::RCP< Stokhos::LTBSparse3Tensor<ordinal_type, value_type> > Cijk =
1403  Stokhos::computeTripleProductTensorLTBBlockLeaf(*setup.basis, symmetric);
1404 
1405  matrix.block =
1406  Stokhos::create_stochastic_product_tensor< TensorType >( *setup.basis,
1407  *Cijk );
1408 
1409  matrix.graph = Kokkos::create_staticcrsgraph<graph_type>(
1410  std::string("test crs graph") , setup.fem_graph );
1411 
1412  matrix.values = block_vector_type(
1413  "matrix" , setup.stoch_length , setup.fem_graph_length );
1414 
1415  typename block_vector_type::HostMirror hM =
1416  Kokkos::create_mirror( matrix.values );
1417 
1418  for ( int iRowFEM = 0 , iEntryFEM = 0 ; iRowFEM < setup.fem_length ; ++iRowFEM ) {
1419  for ( size_t iRowEntryFEM = 0 ; iRowEntryFEM < setup.fem_graph[iRowFEM].size() ; ++iRowEntryFEM , ++iEntryFEM ) {
1420  const int iColFEM = setup.fem_graph[iRowFEM][iRowEntryFEM] ;
1421 
1422  for ( int k = 0 ; k < setup.stoch_length ; ++k ) {
1423  hM(k,iEntryFEM) =
1424  generate_matrix_coefficient<ScalarType>(
1425  setup.fem_length , setup.stoch_length , iRowFEM , iColFEM , k );
1426  }
1427  }
1428  }
1429 
1430  Kokkos::deep_copy( matrix.values , hM );
1431 
1432  //------------------------------
1433 
1434  Stokhos::multiply( matrix , x , y );
1435 
1436  typename block_vector_type::HostMirror hy = Kokkos::create_mirror( y );
1437  Kokkos::deep_copy( hy , y );
1438 
1439  bool success = setup.test_commuted(hy, out);
1440  return success;
1441 }
1442 
1443 }
1444 
1445 #endif
KokkosKernelsUnitTest::map_fem_graph_coord
IntType map_fem_graph_coord(const IntType &N, const IntType &i, const IntType &j, const IntType &k)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:91
KokkosKernelsUnitTest::UnitTestSetup::Cijk_type
Stokhos::Sparse3Tensor< int, value_type > Cijk_type
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:171
Epetra_Map
Stokhos::MatrixFreeOperator
An Epetra operator representing the block stochastic Galerkin operator.
Definition: Stokhos_MatrixFreeOperator.hpp:60
Stokhos::StochasticProductTensor
Bases defined by combinatorial product of polynomial bases.
Definition: Stokhos_StochasticProductTensor.hpp:79
Stokhos::EpetraOperatorOrthogPoly
A container class storing an orthogonal polynomial whose coefficients are vectors, operators, or in general any type that would have an expensive copy constructor.
Definition: Stokhos_EpetraOperatorOrthogPoly.hpp:55
Epetra_CrsGraph::FillComplete
int FillComplete()
Teuchos_UnitTestHelpers.hpp
Stokhos::Sparse3Tensor::k_begin
k_iterator k_begin() const
Iterator pointing to first k entry.
Definition: Stokhos_Sparse3TensorImp.hpp:287
Stokhos_FlatSparse3Tensor.hpp
Stokhos::TotalOrderBasis
Multivariate orthogonal polynomial basis generated from a total order tensor product of univariate po...
Definition: Stokhos_TotalOrderBasis.hpp:68
Epetra_CrsGraph::NumGlobalIndices
int NumGlobalIndices(long long Row) const
KokkosKernelsUnitTest::test_lexo_block_tensor
bool test_lexo_block_tensor(const UnitTestSetup< Device > &setup, Teuchos::FancyOStream &out)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:1354
Stokhos::ParallelData::getMultiComm
Teuchos::RCP< const EpetraExt::MultiComm > getMultiComm() const
Get global comm.
Definition: Stokhos_ParallelData.hpp:76
Stokhos::Sparse3Tensor::num_j
ordinal_type num_j(const k_iterator &k) const
Number of j entries in C(i,j,k) for given k.
Definition: Stokhos_Sparse3TensorImp.hpp:226
Stokhos::Sparse3Tensor
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
Definition: Stokhos_Sparse3Tensor.hpp:56
Epetra_BlockMap::MyGlobalElements
int MyGlobalElements(int *MyGlobalElementList) const
Teuchos::ParameterList::set
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
Stokhos_BlockCrsMatrix.hpp
Stokhos::ProductContainer::init
void init(const value_type &val)
Initialize coefficients.
Definition: Stokhos_ProductContainerImp.hpp:219
Stokhos::ProductEpetraVector::getBlockVector
Teuchos::RCP< EpetraExt::BlockVector > getBlockVector()
Get block vector.
Definition: Stokhos_ProductEpetraVector.cpp:240
KokkosKernelsUnitTest::test_crs_product_tensor
bool test_crs_product_tensor(const UnitTestSetup< Device > &setup, Teuchos::FancyOStream &out, const Teuchos::ParameterList &params=Teuchos::ParameterList())
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:1190
Stokhos::Sparse3Tensor::j_begin
kj_iterator j_begin(const k_iterator &k) const
Iterator pointing to first j entry for given k.
Definition: Stokhos_Sparse3TensorImp.hpp:335
KokkosKernelsUnitTest::test_crs_flat_commuted
bool test_crs_flat_commuted(const UnitTestSetup< Device > &setup, Teuchos::FancyOStream &out)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:895
KokkosKernelsUnitTest::UnitTestSetup::test_original_flat
bool test_original_flat(const vec_type &y, Teuchos::FancyOStream &out) const
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:439
KokkosKernelsUnitTest::UnitTestSetup::test_commuted
bool test_commuted(const vec_type &y, Teuchos::FancyOStream &out) const
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:392
KokkosKernelsUnitTest::UnitTestSetup::abstract_basis_type
Stokhos::OneDOrthogPolyBasis< int, value_type > abstract_basis_type
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:166
Stokhos_Epetra.hpp
Stokhos::SymmetricDiagonalSpec
Symmetric diagonal storage for a dense matrix.
Definition: Stokhos_SymmetricDiagonalSpec.hpp:63
KokkosKernelsUnitTest::generate_vector_coefficient
scalar generate_vector_coefficient(const ordinal nFEM, const ordinal nStoch, const ordinal iColFEM, const ordinal iStoch)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:152
Stokhos::MultiIndex::index
Teuchos::Array< ordinal_type > index
index terms
Definition: Stokhos_ProductBasisUtils.hpp:194
Stokhos_Sparse3TensorUtilities.hpp
Epetra_MpiComm
KokkosKernelsUnitTest::UnitTestSetup::product_basis_type
Stokhos::TotalOrderBasis< int, value_type, less_type > product_basis_type
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:170
Stokhos::computeTripleProductTensorLTBBlockLeaf
Teuchos::RCP< LTBSparse3Tensor< ordinal_type, value_type > > computeTripleProductTensorLTBBlockLeaf(const TotalOrderBasis< ordinal_type, value_type, LexographicLess< MultiIndex< ordinal_type > > > &product_basis, bool symmetric=false)
Definition: Stokhos_LTBSparse3Tensor.hpp:368
Epetra_CrsGraph::InsertGlobalIndices
int InsertGlobalIndices(int_type GlobalRow, int NumIndices, int_type *Indices)
Epetra_SerialComm.h
Stokhos::ParallelData::getStochasticComm
Teuchos::RCP< const Epetra_Comm > getStochasticComm() const
Get stochastic comm.
Definition: Stokhos_ParallelData.hpp:80
Epetra_CrsMatrix::FillComplete
int FillComplete(bool OptimizeDataStorage=true)
Stokhos_SymmetricDiagonalSpec.hpp
Stokhos::ParallelData::getEpetraCijk
Teuchos::RCP< const Stokhos::EpetraSparse3Tensor > getEpetraCijk() const
Get Epetra Cijk.
Definition: Stokhos_ParallelData.hpp:88
Stokhos::Sparse3Tensor::j_end
kj_iterator j_end(const k_iterator &k) const
Iterator pointing to last j entry for given k.
Definition: Stokhos_Sparse3TensorImp.hpp:347
Stokhos_StochasticProductTensor.hpp
Stokhos::SparseArrayIterator
Bi-directional iterator for traversing a sparse array.
Definition: Stokhos_SparseArray.hpp:56
j
j
Definition: Sacado_Fad_Exp_MP_Vector.hpp:527
KokkosKernelsUnitTest::test_linear_tensor
bool test_linear_tensor(const UnitTestSetup< Device > &setup, Teuchos::FancyOStream &out, const bool symmetric)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:1270
Stokhos::multiply
void multiply(const CrsMatrix< MatrixValue, Device, Layout > &A, const InputMultiVectorType &x, OutputMultiVectorType &y, const std::vector< OrdinalType > &col_indices, SingleColumnMultivectorMultiply)
Definition: Stokhos_CrsMatrix.hpp:621
Stokhos::CrsProductTensor
Sparse product tensor with replicated entries to provide subsets with a given coordinate.
Definition: Stokhos_CrsProductTensor.hpp:80
Epetra_BlockMap::NumMyElements
int NumMyElements() const
Stokhos::MultiIndex
A multidimensional index.
Definition: Stokhos_ProductBasisUtils.hpp:79
Teuchos::rcp
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Stokhos::EpetraSparse3Tensor::getStochasticRowMap
Teuchos::RCP< const Epetra_BlockMap > getStochasticRowMap() const
Get stochastic row map.
Definition: Stokhos_EpetraSparse3Tensor.hpp:124
Stokhos_CrsProductTensor.hpp
Stokhos_LinearSparse3Tensor.hpp
Stokhos::EpetraVectorOrthogPoly
A container class storing an orthogonal polynomial whose coefficients are vectors, operators, or in general any type that would have an expensive copy constructor.
Definition: Stokhos_EpetraVectorOrthogPoly.hpp:55
Teuchos::basic_FancyOStream
Stokhos::LinearSparse3Tensor
Sparse product tensor with replicated entries to provide subsets with a given coordinate.
Definition: Stokhos_LinearSparse3Tensor.hpp:62
Stokhos::ProductEpetraVector
A container class for products of Epetra_Vector&#39;s.
Definition: Stokhos_ProductEpetraVector.hpp:56
Stokhos_Update.hpp
KokkosKernelsUnitTest::test_crs_flat_original
bool test_crs_flat_original(const UnitTestSetup< Device > &setup, Teuchos::FancyOStream &out)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:1044
Stokhos::LexicographicBlockSparse3Tensor
Sparse product tensor with replicated entries to provide subsets with a given coordinate.
Definition: Stokhos_LexicographicBlockSparse3Tensor.hpp:64
Kokkos::deep_copy
void deep_copy(const Stokhos::CrsMatrix< ValueType, DstDevice, Layout > &dst, const Stokhos::CrsMatrix< ValueType, SrcDevice, Layout > &src)
Definition: Stokhos_CrsMatrix.hpp:688
Sacado::UQ::abs
KOKKOS_INLINE_FUNCTION PCE< Storage > abs(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:1181
Stokhos::MatrixFreeOperator::Apply
virtual int Apply(const Epetra_MultiVector &Input, Epetra_MultiVector &Result) const
Returns the result of a Epetra_Operator applied to a Epetra_MultiVector Input in Result as described ...
Definition: Stokhos_MatrixFreeOperator.cpp:209
Teuchos::Array::resize
void resize(size_type new_size, const value_type &x=value_type())
KokkosKernelsUnitTest::UnitTestSetup::test_original
bool test_original(const std::vector< vec_type > &y, Teuchos::FancyOStream &out) const
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:322
KokkosKernelsUnitTest::UnitTestSetup::less_type
Stokhos::LexographicLess< Stokhos::MultiIndex< int > > less_type
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:169
Stokhos_UnitTestHelpers.hpp
Epetra_SerialComm
Stokhos::Sparse3Tensor::k_end
k_iterator k_end() const
Iterator pointing to last k entry.
Definition: Stokhos_Sparse3TensorImp.hpp:299
Epetra_MpiComm.h
Stokhos::LegendreBasis
Legendre polynomial basis.
Definition: Stokhos_LegendreBasis.hpp:67
KokkosKernelsUnitTest::generate_matrix_coefficient
scalar generate_matrix_coefficient(const ordinal nFEM, const ordinal nStoch, const ordinal iRowFEM, const ordinal iColFEM, const ordinal iStoch)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:135
KokkosKernelsUnitTest::UnitTestSetup::test_original
bool test_original(const multi_vec_type &y, Teuchos::FancyOStream &out) const
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:345
Teuchos_TestingHelpers.hpp
Stokhos::OneDOrthogPolyBasis
Abstract base class for 1-D orthogonal polynomials.
Definition: Stokhos_OneDOrthogPolyBasis.hpp:81
Stokhos::sparse3Tensor2CrsGraph
Teuchos::RCP< Epetra_CrsGraph > sparse3Tensor2CrsGraph(const Stokhos::OrthogPolyBasis< ordinal_type, value_type > &basis, const Stokhos::Sparse3Tensor< ordinal_type, value_type > &Cijk, const Epetra_Comm &comm)
Build an Epetra_CrsGraph from a sparse 3 tensor.
Definition: Stokhos_Sparse3TensorUtilities.hpp:67
KokkosKernelsUnitTest::generate_fem_graph
ordinal generate_fem_graph(ordinal N, std::vector< std::vector< ordinal > > &graph)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:101
KokkosBlas::update
std::enable_if< Kokkos::is_view_uq_pce< Kokkos::View< XD, XP...> >::value &&Kokkos::is_view_uq_pce< Kokkos::View< YD, YP...> >::value &&Kokkos::is_view_uq_pce< Kokkos::View< ZD, ZP...> >::value >::type update(const typename Kokkos::View< XD, XP...>::array_type::non_const_value_type &alpha, const Kokkos::View< XD, XP...> &x, const typename Kokkos::View< YD, YP...>::array_type::non_const_value_type &beta, const Kokkos::View< YD, YP...> &y, const typename Kokkos::View< ZD, ZP...>::array_type::non_const_value_type &gamma, const Kokkos::View< ZD, ZP...> &z)
Definition: Kokkos_Blas1_UQ_PCE.hpp:218
Stokhos::BlockCrsMatrix
CRS matrix of dense blocks.
Definition: Stokhos_BlockCrsMatrix.hpp:64
Epetra_CrsMatrix
Teuchos_TestForException.hpp
Stokhos::LexographicLess
A comparison functor implementing a strict weak ordering based lexographic ordering.
Definition: Stokhos_ProductBasisUtils.hpp:799
Stokhos_SimpleTiledCrsProductTensor.hpp
int
int
Definition: tpetra_mat_vec.cpp:243
Stokhos::ParallelData::getSpatialComm
Teuchos::RCP< const Epetra_Comm > getSpatialComm() const
Get spatial comm.
Definition: Stokhos_ParallelData.hpp:84
Stokhos_CrsMatrix.hpp
KokkosKernelsUnitTest::UnitTestSetup::basis_type
Stokhos::LegendreBasis< int, value_type > basis_type
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:167
DynamicVecTest::vec_type
Sacado::MP::Vector< storage_type > vec_type
Definition: Stokhos_SacadoMPVectorCommTests.cpp:685
Epetra_CrsGraph::ExtractGlobalRowCopy
int ExtractGlobalRowCopy(int_type Row, int LenOfIndices, int &NumIndices, int_type *Indices) const
Stokhos::MatrixFreeOperator::setupOperator
virtual void setupOperator(const Teuchos::RCP< Stokhos::EpetraOperatorOrthogPoly > &poly)
Setup operator.
Definition: Stokhos_MatrixFreeOperator.cpp:173
KokkosKernelsUnitTest::UnitTestSetup::test_commuted_perm
bool test_commuted_perm(const vec_type &y, Teuchos::FancyOStream &out) const
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:368
Stokhos::Sparse3Tensor::i_begin
kji_iterator i_begin(const kj_iterator &j) const
Iterator pointing to first i entry for given j and k.
Definition: Stokhos_Sparse3TensorImp.hpp:383
Stokhos_CooProductTensor.hpp
Stokhos::Sparse3Tensor::i_end
kji_iterator i_end(const kj_iterator &j) const
Iterator pointing to last i entry for given j and k.
Definition: Stokhos_Sparse3TensorImp.hpp:395
Stokhos_LexicographicBlockSparse3Tensor.hpp
Stokhos_TiledCrsProductTensor.hpp
Epetra_CrsMatrix::ReplaceGlobalValues
virtual int ReplaceGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
Kokkos::create_mirror
Stokhos::CrsMatrix< ValueType, Device, Layout >::HostMirror create_mirror(const Stokhos::CrsMatrix< ValueType, Device, Layout > &A)
Definition: Stokhos_CrsMatrix.hpp:667
Stokhos::update
void update(const ValueType &alpha, VectorType &x, const ValueType &beta, const VectorType &y)
Definition: Stokhos_Update.hpp:87
KokkosKernelsUnitTest::UnitTestSetup::test_commuted_flat
bool test_commuted_flat(const vec_type &y, Teuchos::FancyOStream &out) const
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:416
Stokhos_FlatSparse3Tensor_kji.hpp
KokkosSparse::spmv
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)
Definition: Kokkos_CrsMatrix_UQ_PCE.hpp:1235
Epetra_CrsGraph
Epetra_LocalMap
KokkosKernelsUnitTest::test_crs_matrix_free
bool test_crs_matrix_free(const UnitTestSetup< Device > &setup, Teuchos::FancyOStream &out)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:464
UnitTestSetup::Cijk
RCP< Stokhos::Sparse3Tensor< int, double > > Cijk
Definition: Stokhos_SacadoPCECommTests.cpp:63
KokkosKernelsUnitTest::test_crs_matrix_free_view
bool test_crs_matrix_free_view(const UnitTestSetup< Device > &setup, Teuchos::FancyOStream &out)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:574
KokkosKernelsUnitTest::test_crs_dense_block
bool test_crs_dense_block(const UnitTestSetup< Device > &setup, Teuchos::FancyOStream &out)
Definition: Stokhos_KokkosArrayKernelsUnitTest.hpp:784
Generated by   doxygen 1.8.5