Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
MPAssembly/fenl_functors.hpp
Go to the documentation of this file.
1/*
2//@HEADER
3// ************************************************************************
4//
5// Kokkos: Manycore Performance-Portable Multidimensional Arrays
6// Copyright (2012) Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact H. Carter Edwards (hcedwar@sandia.gov)
39//
40// ************************************************************************
41//@HEADER
42*/
43
44#ifndef KOKKOS_EXAMPLE_FENLFUNCTORS_HPP
45#define KOKKOS_EXAMPLE_FENLFUNCTORS_HPP
46
47#include <stdio.h>
48
49#include <iostream>
50#include <fstream>
51#include <iomanip>
52#include <cstdlib>
53#include <cmath>
54#include <limits>
55
56#include <Kokkos_Pair.hpp>
57#include <Kokkos_UnorderedMap.hpp>
58
59#include <Kokkos_Timer.hpp>
60
61#include <BoxElemFixture.hpp>
62#include <HexElement.hpp>
63
65
66//----------------------------------------------------------------------------
67//----------------------------------------------------------------------------
68
69namespace Kokkos {
70namespace Example {
71namespace FENL {
72
73struct DeviceConfig {
74 struct Dim3 {
75 size_t x, y, z;
76 Dim3(const size_t x_, const size_t y_ = 1, const size_t z_ = 1) :
77 x(x_), y(y_), z(z_) {}
78 };
79
80 Dim3 block_dim;
81 size_t num_blocks;
83
84 DeviceConfig(const size_t num_blocks_ = 0,
85 const size_t threads_per_block_x_ = 0,
86 const size_t threads_per_block_y_ = 0,
87 const size_t threads_per_block_z_ = 1) :
88 block_dim(threads_per_block_x_,threads_per_block_y_,threads_per_block_z_),
89 num_blocks(num_blocks_),
91 {}
92};
93
94template< typename ValueType , class Space >
95struct CrsMatrix {
96#ifdef KOKKOS_ENABLE_DEPRECATED_CODE // Don't remove this until Kokkos has removed the deprecated code path probably around September 2018
97 typedef Kokkos::StaticCrsGraph< unsigned , Space , void , unsigned > StaticCrsGraphType ;
98#else
99 typedef Kokkos::StaticCrsGraph< unsigned , Space , void , void , unsigned > StaticCrsGraphType ;
100#endif
101 typedef View< ValueType * , Space > values_type ;
102
105
107
108 CrsMatrix( const StaticCrsGraphType & arg_graph )
109 : graph( arg_graph )
110 , values( "crs_matrix_values" , arg_graph.entries.extent(0) )
111 {}
112};
113
114// Traits class for creating strided local views for embedded ensemble-based,
115// specialized for ensemble UQ scalar type
116template <typename ViewType, typename Enabled = void>
117struct LocalViewTraits {
118 typedef ViewType view_type;
119 // typedef Kokkos::View<typename view_type::data_type,
120 // typename view_type::array_layout,
121 // typename view_type::execution_space,
122 // Kokkos::MemoryUnmanaged> local_view_type;
124 typedef typename view_type::value_type local_value_type;
125 static const bool use_team = false;
126 KOKKOS_INLINE_FUNCTION
128 const unsigned local_rank)
129 { return v; }
130};
131
132#if defined( KOKKOS_ENABLE_CUDA )
133
134template <typename ViewType>
135struct LocalViewTraits<
136 ViewType,
137 typename std::enable_if< std::is_same<typename ViewType::execution_space,
138 Kokkos::Cuda>::value &&
139 Kokkos::is_view_mp_vector<ViewType>::value
140 >::type > {
141 typedef ViewType view_type;
143 typedef typename local_view_type::value_type local_value_type;
144 static const bool use_team = true;
145
146 KOKKOS_INLINE_FUNCTION
148 const unsigned local_rank)
149 {
150 return Kokkos::partition<1>(v, local_rank);
151 }
152};
153
154#endif /* #if defined( KOKKOS_ENABLE_CUDA ) */
155
156// Compute DeviceConfig struct's based on scalar type
157template <typename ScalarType>
158struct CreateDeviceConfigs {
159 static void eval( Kokkos::Example::FENL::DeviceConfig& dev_config_elem,
160 Kokkos::Example::FENL::DeviceConfig& dev_config_bc ) {
161 dev_config_elem = Kokkos::Example::FENL::DeviceConfig( 0 , 1 , 1 );
162 dev_config_bc = Kokkos::Example::FENL::DeviceConfig( 0 , 1 , 1 );
163 }
164};
165
166// Compute DeviceConfig struct's based on scalar type
167template <typename StorageType>
168struct CreateDeviceConfigs< Sacado::MP::Vector<StorageType> > {
169 typedef typename StorageType::execution_space execution_space;
170 static void eval( Kokkos::Example::FENL::DeviceConfig& dev_config_elem,
171 Kokkos::Example::FENL::DeviceConfig& dev_config_bc ) {
172 static const unsigned VectorSize = StorageType::static_size;
173#if defined( KOKKOS_ENABLE_CUDA )
174 enum { is_cuda = std::is_same< execution_space, Kokkos::Cuda >::value };
175#else
176 enum { is_cuda = false };
177#endif /* #if defined( KOKKOS_ENABLE_CUDA ) */
178 if ( is_cuda ) {
179 dev_config_elem = Kokkos::Example::FENL::DeviceConfig( 0 , VectorSize , 64/VectorSize );
180 dev_config_bc = Kokkos::Example::FENL::DeviceConfig( 0 , VectorSize , 256/VectorSize );
181 }
182 else {
183 dev_config_elem = Kokkos::Example::FENL::DeviceConfig( 0 , 1 , 1 );
184 dev_config_bc = Kokkos::Example::FENL::DeviceConfig( 0 , 1 , 1 );
185 }
186 }
187};
188
189template< class ElemNodeIdView , class CrsGraphType , unsigned ElemNode >
190class NodeNodeGraph {
191public:
192
193 typedef typename ElemNodeIdView::execution_space execution_space ;
194 typedef pair<unsigned,unsigned> key_type ;
195
196 typedef Kokkos::UnorderedMap< key_type, void , execution_space > SetType ;
197 typedef typename CrsGraphType::row_map_type::non_const_type RowMapType ;
198 typedef Kokkos::View< unsigned , execution_space > UnsignedValue ;
199
200 // Static dimensions of 0 generate compiler warnings or errors.
201 typedef Kokkos::View< unsigned*[ElemNode][ElemNode] , execution_space >
203
204private:
205
211
212 const unsigned node_count ;
213 const ElemNodeIdView elem_node_id ;
219
220public:
221
222 CrsGraphType graph ;
224
225 struct Times
226 {
227 double ratio;
228 double fill_node_set;
229 double scan_node_count;
230 double fill_graph_entries;
231 double sort_graph_entries;
232 double fill_element_graph;
233 };
234
235 NodeNodeGraph( const ElemNodeIdView & arg_elem_node_id ,
236 const unsigned arg_node_count,
237 Times & results
238 )
239 : node_count(arg_node_count)
240 , elem_node_id( arg_elem_node_id )
241 , row_total( "row_total" )
242 , row_count(Kokkos::ViewAllocateWithoutInitializing("row_count") , node_count ) // will deep_copy to 0 inside loop
243 , row_map( "graph_row_map" , node_count + 1 )
244 , node_node_set()
246 , graph()
247 , elem_graph()
248 {
249 //--------------------------------
250 // Guess at capacity required for the map:
251
252 Kokkos::Timer wall_clock ;
253
254 wall_clock.reset();
256
257 // upper bound on the capacity
258 size_t set_capacity = (((28ull * node_count) / 2ull)*4ull)/3ull;
259
260
261 // Increase capacity until the (node,node) map is successfully filled.
262 {
263 // Zero the row count to restart the fill
265
266 node_node_set = SetType( set_capacity );
267
268 // May be larger that requested:
269 set_capacity = node_node_set.capacity();
270
271 Kokkos::parallel_for( elem_node_id.extent(0) , *this );
272 }
273
274 execution_space().fence();
275 results.ratio = (double)node_node_set.size() / (double)node_node_set.capacity();
276 results.fill_node_set = wall_clock.seconds();
277 //--------------------------------
278
279 wall_clock.reset();
281
282 // Exclusive scan of row_count into row_map
283 // including the final total in the 'node_count + 1' position.
284 // Zero the 'row_count' values.
285 Kokkos::parallel_scan( node_count , *this );
286
287 // Zero the row count for the fill:
289
290 unsigned graph_entry_count = 0 ;
291
292 Kokkos::deep_copy( graph_entry_count , row_total );
293
294 // Assign graph's row_map and allocate graph's entries
295 graph.row_map = row_map ;
296 graph.entries = typename CrsGraphType::entries_type( "graph_entries" , graph_entry_count );
297
298 //--------------------------------
299 // Fill graph's entries from the (node,node) set.
300
301 execution_space().fence();
302 results.scan_node_count = wall_clock.seconds();
303
304 wall_clock.reset();
306 Kokkos::parallel_for( node_node_set.capacity() , *this );
307
308 execution_space().fence();
309 results.fill_graph_entries = wall_clock.seconds();
310
311 //--------------------------------
312 // Done with the temporary sets and arrays
313 wall_clock.reset();
315
319 node_node_set.clear();
320
321 //--------------------------------
322
323 Kokkos::parallel_for( node_count , *this );
324
325 execution_space().fence();
326 results.sort_graph_entries = wall_clock.seconds();
327
328 //--------------------------------
329 // Element-to-graph mapping:
330 wall_clock.reset();
332 elem_graph = ElemGraphType("elem_graph", elem_node_id.extent(0) );
333 Kokkos::parallel_for( elem_node_id.extent(0) , *this );
334
335 execution_space().fence();
336 results.fill_element_graph = wall_clock.seconds();
337 }
338
339 //------------------------------------
340 // parallel_for: create map and count row length
341
342 KOKKOS_INLINE_FUNCTION
343 void fill_set( const unsigned ielem ) const
344 {
345 // Loop over element's (row_local_node,col_local_node) pairs:
346 for ( unsigned row_local_node = 0 ; row_local_node < elem_node_id.extent(1) ; ++row_local_node ) {
347
348 const unsigned row_node = elem_node_id( ielem , row_local_node );
349
350 for ( unsigned col_local_node = row_local_node ; col_local_node < elem_node_id.extent(1) ; ++col_local_node ) {
351
352 const unsigned col_node = elem_node_id( ielem , col_local_node );
353
354 // If either node is locally owned then insert the pair into the unordered map:
355
356 if ( row_node < row_count.extent(0) || col_node < row_count.extent(0) ) {
357
358 const key_type key = (row_node < col_node) ? make_pair( row_node, col_node ) : make_pair( col_node, row_node ) ;
359
360 const typename SetType::insert_result result = node_node_set.insert( key );
361
362 if ( result.success() ) {
363 if ( row_node < row_count.extent(0) ) { atomic_fetch_add( & row_count( row_node ) , (typename RowMapType::value_type)1 ); }
364 if ( col_node < row_count.extent(0) && col_node != row_node ) { atomic_fetch_add( & row_count( col_node ) , (typename RowMapType::value_type)1 ); }
365 }
366 }
367 }
368 }
369 }
370
371 KOKKOS_INLINE_FUNCTION
372 void fill_graph_entries( const unsigned iset ) const
373 {
374 if ( node_node_set.valid_at(iset) ) {
375 const key_type key = node_node_set.key_at(iset) ;
376 const unsigned row_node = key.first ;
377 const unsigned col_node = key.second ;
378
379 if ( row_node < row_count.extent(0) ) {
380 const unsigned offset = graph.row_map( row_node ) + atomic_fetch_add( & row_count( row_node ) , (typename RowMapType::value_type)1 );
381 graph.entries( offset ) = col_node ;
382 }
383
384 if ( col_node < row_count.extent(0) && col_node != row_node ) {
385 const unsigned offset = graph.row_map( col_node ) + atomic_fetch_add( & row_count( col_node ) , (typename RowMapType::value_type)1 );
386 graph.entries( offset ) = row_node ;
387 }
388 }
389 }
390
391 KOKKOS_INLINE_FUNCTION
392 void sort_graph_entries( const unsigned irow ) const
393 {
394 typedef typename CrsGraphType::size_type size_type;
395 const size_type row_beg = graph.row_map( irow );
396 const size_type row_end = graph.row_map( irow + 1 );
397 for ( size_type i = row_beg + 1 ; i < row_end ; ++i ) {
398 const typename CrsGraphType::data_type col = graph.entries(i);
399 size_type j = i ;
400 for ( ; row_beg < j && col < graph.entries(j-1) ; --j ) {
401 graph.entries(j) = graph.entries(j-1);
402 }
403 graph.entries(j) = col ;
404 }
405 }
406
407 KOKKOS_INLINE_FUNCTION
408 void fill_elem_graph_map( const unsigned ielem ) const
409 {
410 typedef typename CrsGraphType::data_type entry_type;
411 for ( unsigned row_local_node = 0 ; row_local_node < elem_node_id.extent(1) ; ++row_local_node ) {
412
413 const unsigned row_node = elem_node_id( ielem , row_local_node );
414
415 for ( unsigned col_local_node = 0 ; col_local_node < elem_node_id.extent(1) ; ++col_local_node ) {
416
417 const unsigned col_node = elem_node_id( ielem , col_local_node );
418
419 entry_type entry = 0 ;
420
421 if ( row_node + 1 < graph.row_map.extent(0) ) {
422
423 const entry_type entry_end = static_cast<entry_type> (graph.row_map( row_node + 1 ));
424
425 entry = graph.row_map( row_node );
426
427 for ( ; entry < entry_end && graph.entries(entry) != static_cast<entry_type> (col_node) ; ++entry );
428
429 if ( entry == entry_end ) entry = ~0u ;
430 }
431
432 elem_graph( ielem , row_local_node , col_local_node ) = entry ;
433 }
434 }
435 }
436
437 KOKKOS_INLINE_FUNCTION
438 void operator()( const unsigned iwork ) const
439 {
440 if ( phase == FILL_NODE_SET ) {
441 fill_set( iwork );
442 }
443 else if ( phase == FILL_GRAPH_ENTRIES ) {
444 fill_graph_entries( iwork );
445 }
446 else if ( phase == SORT_GRAPH_ENTRIES ) {
447 sort_graph_entries( iwork );
448 }
449 else if ( phase == FILL_ELEMENT_GRAPH ) {
450 fill_elem_graph_map( iwork );
451 }
452 }
453
454 //------------------------------------
455 // parallel_scan: row offsets
456
457 typedef unsigned value_type ;
458
459 KOKKOS_INLINE_FUNCTION
460 void operator()( const unsigned irow , unsigned & update , const bool final ) const
461 {
462 // exclusive scan
463 if ( final ) { row_map( irow ) = update ; }
464
465 update += row_count( irow );
466
467 if ( final ) {
468 if ( irow + 1 == row_count.extent(0) ) {
469 row_map( irow + 1 ) = update ;
470 row_total() = update ;
471 }
472 }
473 }
474
475 KOKKOS_INLINE_FUNCTION
476 void init( unsigned & update ) const { update = 0 ; }
477
478 KOKKOS_INLINE_FUNCTION
479 void join( unsigned & update , const unsigned & input ) const { update += input ; }
480
481 //------------------------------------
482};
483
484} /* namespace FENL */
485} /* namespace Example */
486} /* namespace Kokkos */
487
488//----------------------------------------------------------------------------
489//----------------------------------------------------------------------------
490
491namespace Kokkos {
492namespace Example {
493namespace FENL {
494
495struct ElementComputationConstantCoefficient {
496 enum { is_constant = false };
497
498 const double coeff_k ;
499
500 KOKKOS_INLINE_FUNCTION
501 double operator()( double pt[], unsigned ensemble_rank) const
502 { return coeff_k * std::sin(pt[0]) * std::sin(pt[1]) * std::sin(pt[2]); }
503
505 : coeff_k( val ) {}
506
508 : coeff_k( rhs.coeff_k ) {}
509};
510
511// Exponential KL from Stokhos
512template < typename Scalar, typename MeshScalar, typename Device >
513class ExponentialKLCoefficient {
514public:
515
516 // Turn into a meta-function class usable with Sacado::mpl::apply
517 template <typename T1, typename T2 = MeshScalar, typename T3 = Device>
518 struct apply {
520 };
521
522 enum { is_constant = false };
523 typedef Kokkos::View<Scalar*, Kokkos::LayoutLeft, Device> RandomVariableView;
524 typedef typename RandomVariableView::size_type size_type;
525
530
531 rf_type m_rf; // Exponential random field
532 const MeshScalar m_mean; // Mean of random field
533 const MeshScalar m_variance; // Variance of random field
534 const MeshScalar m_corr_len; // Correlation length of random field
535 const size_type m_num_rv; // Number of random variables
536 RandomVariableView m_rv; // KL random variables
537
538public:
539
541 const MeshScalar mean ,
542 const MeshScalar variance ,
543 const MeshScalar correlation_length ,
544 const size_type num_rv ) :
545 m_mean( mean ),
546 m_variance( variance ),
547 m_corr_len( correlation_length ),
548 m_num_rv( num_rv ),
549 m_rv( "KL Random Variables", m_num_rv )
550 {
551 Teuchos::ParameterList solverParams;
552 solverParams.set("Number of KL Terms", int(num_rv));
553 solverParams.set("Mean", mean);
554 solverParams.set("Standard Deviation", std::sqrt(variance));
555 int ndim = 3;
556 Teuchos::Array<double> domain_upper(ndim, 1.0), domain_lower(ndim, 0.0),
557 correlation_lengths(ndim, correlation_length);
558 solverParams.set("Domain Upper Bounds", domain_upper);
559 solverParams.set("Domain Lower Bounds", domain_lower);
560 solverParams.set("Correlation Lengths", correlation_lengths);
561
562 m_rf = rf_type(solverParams);
563 }
564
566 m_rf( rhs.m_rf ) ,
567 m_mean( rhs.m_mean ) ,
568 m_variance( rhs.m_variance ) ,
569 m_corr_len( rhs.m_corr_len ) ,
570 m_num_rv( rhs.m_num_rv ) ,
571 m_rv( rhs.m_rv ) {}
572
573 KOKKOS_INLINE_FUNCTION
574 void setRandomVariables( const RandomVariableView& rv) { m_rv = rv; }
575
576 KOKKOS_INLINE_FUNCTION
578
579 KOKKOS_INLINE_FUNCTION
580 local_scalar_type operator() ( const MeshScalar point[],
581 const size_type ensemble_rank ) const
582 {
583 local_rv_view_type local_rv =
585
586 local_scalar_type val = m_rf.evaluate(point, local_rv);
587
588 return val;
589 }
590};
591
592template< class FiniteElementMeshType , class SparseMatrixType
593 , class CoeffFunctionType = ElementComputationConstantCoefficient
594 >
595class ElementComputation ;
596
597
598template< class ExecutionSpace , BoxElemPart::ElemOrder Order , class CoordinateMap ,
599 typename ScalarType , class CoeffFunctionType >
601 < Kokkos::Example::BoxElemFixture< ExecutionSpace , Order , CoordinateMap >
602 , Kokkos::Example::FENL::CrsMatrix< ScalarType , ExecutionSpace >
603 , CoeffFunctionType >
604{
605public:
606
609
610 //------------------------------------
611
612 typedef ExecutionSpace execution_space ;
613 typedef ScalarType scalar_type ;
614
618 typedef Kokkos::View< scalar_type* , Kokkos::LayoutLeft, execution_space > vector_type ;
619
620 //------------------------------------
621
627 static const bool use_team = local_vector_view_traits::use_team;
628
629 static const unsigned SpatialDim = element_data_type::spatial_dimension ;
630 static const unsigned TensorDim = SpatialDim * SpatialDim ;
631 static const unsigned ElemNodeCount = element_data_type::element_node_count ;
632 static const unsigned FunctionCount = element_data_type::function_count ;
633 static const unsigned IntegrationCount = element_data_type::integration_count ;
634
635 //------------------------------------
636
639 typedef Kokkos::View< scalar_type*[FunctionCount][FunctionCount] , execution_space > elem_matrices_type ;
640 typedef Kokkos::View< scalar_type*[FunctionCount] , execution_space > elem_vectors_type ;
641
646
648
649 //------------------------------------
650
651
652 //------------------------------------
653 // Computational data:
654
664 const CoeffFunctionType coeff_function ;
666
668 : elem_data()
669 , elem_node_ids( rhs.elem_node_ids )
670 , node_coords( rhs.node_coords )
671 , elem_graph( rhs.elem_graph )
672 , elem_jacobians( rhs.elem_jacobians )
673 , elem_residuals( rhs.elem_residuals )
674 , solution( rhs.solution )
675 , residual( rhs.residual )
676 , jacobian( rhs.jacobian )
677 , coeff_function( rhs.coeff_function )
678 , dev_config( rhs.dev_config )
679 {}
680
681 // If the element->sparse_matrix graph is provided then perform atomic updates
682 // Otherwise fill per-element contributions for subequent gather-add into a residual and jacobian.
683 ElementComputation( const mesh_type & arg_mesh ,
684 const CoeffFunctionType & arg_coeff_function ,
685 const vector_type & arg_solution ,
686 const elem_graph_type & arg_elem_graph ,
687 const sparse_matrix_type & arg_jacobian ,
688 const vector_type & arg_residual ,
689 const Kokkos::Example::FENL::DeviceConfig arg_dev_config )
690 : elem_data()
691 , elem_node_ids( arg_mesh.elem_node() )
692 , node_coords( arg_mesh.node_coord() )
693 , elem_graph( arg_elem_graph )
694 , elem_jacobians()
695 , elem_residuals()
696 , solution( arg_solution )
697 , residual( arg_residual )
698 , jacobian( arg_jacobian )
699 , coeff_function( arg_coeff_function )
700 , dev_config( arg_dev_config )
701 {}
702
703 //------------------------------------
704
705 void apply() const
706 {
707 const size_t nelem = elem_node_ids.extent(0);
708 if ( use_team ) {
709 const size_t team_size = dev_config.block_dim.x * dev_config.block_dim.y;
710 const size_t league_size =
711 (nelem + dev_config.block_dim.y-1) / dev_config.block_dim.y;
712 Kokkos::TeamPolicy< execution_space > config( league_size, team_size );
713 parallel_for( config , *this );
714 }
715 else {
716 parallel_for( nelem , *this );
717 }
718 }
719
720 //------------------------------------
721
722 static const unsigned FLOPS_transform_gradients =
723 /* Jacobian */ FunctionCount * TensorDim * 2 +
724 /* Inverse jacobian */ TensorDim * 6 + 6 +
725 /* Gradient transform */ FunctionCount * 15 ;
726
727 KOKKOS_INLINE_FUNCTION
729 const double grad[][ FunctionCount ] , // Gradient of bases master element
730 const double x[] ,
731 const double y[] ,
732 const double z[] ,
733 double dpsidx[] ,
734 double dpsidy[] ,
735 double dpsidz[] ) const
736 {
737 enum { j11 = 0 , j12 = 1 , j13 = 2 ,
738 j21 = 3 , j22 = 4 , j23 = 5 ,
739 j31 = 6 , j32 = 7 , j33 = 8 };
740
741 // Jacobian accumulation:
742
743 double J[ TensorDim ] = { 0, 0, 0, 0, 0, 0, 0, 0, 0 };
744
745 for( unsigned i = 0; i < FunctionCount ; ++i ) {
746 const double x1 = x[i] ;
747 const double x2 = y[i] ;
748 const double x3 = z[i] ;
749
750 const double g1 = grad[0][i] ;
751 const double g2 = grad[1][i] ;
752 const double g3 = grad[2][i] ;
753
754 J[j11] += g1 * x1 ;
755 J[j12] += g1 * x2 ;
756 J[j13] += g1 * x3 ;
757
758 J[j21] += g2 * x1 ;
759 J[j22] += g2 * x2 ;
760 J[j23] += g2 * x3 ;
761
762 J[j31] += g3 * x1 ;
763 J[j32] += g3 * x2 ;
764 J[j33] += g3 * x3 ;
765 }
766
767 // Inverse jacobian:
768
769 double invJ[ TensorDim ] = {
770 static_cast<double>( J[j22] * J[j33] - J[j23] * J[j32] ) ,
771 static_cast<double>( J[j13] * J[j32] - J[j12] * J[j33] ) ,
772 static_cast<double>( J[j12] * J[j23] - J[j13] * J[j22] ) ,
773
774 static_cast<double>( J[j23] * J[j31] - J[j21] * J[j33] ) ,
775 static_cast<double>( J[j11] * J[j33] - J[j13] * J[j31] ) ,
776 static_cast<double>( J[j13] * J[j21] - J[j11] * J[j23] ) ,
777
778 static_cast<double>( J[j21] * J[j32] - J[j22] * J[j31] ) ,
779 static_cast<double>( J[j12] * J[j31] - J[j11] * J[j32] ) ,
780 static_cast<double>( J[j11] * J[j22] - J[j12] * J[j21] ) };
781
782 const double detJ = J[j11] * invJ[j11] +
783 J[j21] * invJ[j12] +
784 J[j31] * invJ[j13] ;
785
786 const double detJinv = 1.0 / detJ ;
787
788 for ( unsigned i = 0 ; i < TensorDim ; ++i ) { invJ[i] *= detJinv ; }
789
790 // Transform gradients:
791
792 for( unsigned i = 0; i < FunctionCount ; ++i ) {
793 const double g0 = grad[0][i];
794 const double g1 = grad[1][i];
795 const double g2 = grad[2][i];
796
797 dpsidx[i] = g0 * invJ[j11] + g1 * invJ[j12] + g2 * invJ[j13];
798 dpsidy[i] = g0 * invJ[j21] + g1 * invJ[j22] + g2 * invJ[j23];
799 dpsidz[i] = g0 * invJ[j31] + g1 * invJ[j32] + g2 * invJ[j33];
800 }
801
802 return detJ ;
803 }
804
805 KOKKOS_INLINE_FUNCTION
807 const local_scalar_type dof_values[] ,
808 const double dpsidx[] ,
809 const double dpsidy[] ,
810 const double dpsidz[] ,
811 const double detJ ,
812 const local_scalar_type coeff_k ,
813 const double integ_weight ,
814 const double bases_vals[] ,
815 local_scalar_type elem_res[] ,
816 local_scalar_type elem_mat[][ FunctionCount ] ) const
817 {
818 local_scalar_type value_at_pt = 0 ;
819 local_scalar_type gradx_at_pt = 0 ;
820 local_scalar_type grady_at_pt = 0 ;
821 local_scalar_type gradz_at_pt = 0 ;
822
823 for ( unsigned m = 0 ; m < FunctionCount ; m++ ) {
824 value_at_pt += dof_values[m] * bases_vals[m] ;
825 gradx_at_pt += dof_values[m] * dpsidx[m] ;
826 grady_at_pt += dof_values[m] * dpsidy[m] ;
827 gradz_at_pt += dof_values[m] * dpsidz[m] ;
828 }
829
830 const local_scalar_type k_detJ_weight = coeff_k * detJ * integ_weight ;
831 const local_scalar_type res_val = value_at_pt * value_at_pt * detJ * integ_weight ;
832 const local_scalar_type mat_val = 2.0 * value_at_pt * detJ * integ_weight ;
833
834 // $$ R_i = \int_{\Omega} \nabla \phi_i \cdot (k \nabla T) + \phi_i T^2 d \Omega $$
835 // $$ J_{i,j} = \frac{\partial R_i}{\partial T_j} = \int_{\Omega} k \nabla \phi_i \cdot \nabla \phi_j + 2 \phi_i \phi_j T d \Omega $$
836
837 for ( unsigned m = 0; m < FunctionCount; ++m) {
838 local_scalar_type * const mat = elem_mat[m] ;
839 const double bases_val_m = bases_vals[m];
840 const double dpsidx_m = dpsidx[m] ;
841 const double dpsidy_m = dpsidy[m] ;
842 const double dpsidz_m = dpsidz[m] ;
843
844 elem_res[m] += k_detJ_weight * ( dpsidx_m * gradx_at_pt +
845 dpsidy_m * grady_at_pt +
846 dpsidz_m * gradz_at_pt ) +
847 res_val * bases_val_m ;
848
849 for( unsigned n = 0; n < FunctionCount; n++) {
850
851 mat[n] += k_detJ_weight * ( dpsidx_m * dpsidx[n] +
852 dpsidy_m * dpsidy[n] +
853 dpsidz_m * dpsidz[n] ) +
854 mat_val * bases_val_m * bases_vals[n];
855 }
856 }
857 }
858
859 KOKKOS_INLINE_FUNCTION
860 void operator()( const typename TeamPolicy< execution_space >::member_type & dev ) const
861 {
862
863 const unsigned num_ensemble_threads = dev_config.block_dim.x ;
864 const unsigned num_element_threads = dev_config.block_dim.y ;
865 const unsigned element_rank = dev.team_rank() / num_ensemble_threads ;
866 const unsigned ensemble_rank = dev.team_rank() % num_ensemble_threads ;
867
868 const unsigned ielem =
869 dev.league_rank() * num_element_threads + element_rank;
870
871 if (ielem >= elem_node_ids.extent(0))
872 return;
873
874 (*this)( ielem, ensemble_rank );
875
876 }
877
878 KOKKOS_INLINE_FUNCTION
879 void operator()( const unsigned ielem ,
880 const unsigned ensemble_rank = 0 ) const
881 {
882 local_vector_type local_solution =
883 local_vector_view_traits::create_local_view(solution,
884 ensemble_rank);
885 local_vector_type local_residual =
886 local_vector_view_traits::create_local_view(residual,
887 ensemble_rank);
888 local_matrix_type local_jacobian_values =
889 local_matrix_view_traits::create_local_view(jacobian.values,
890 ensemble_rank);
891
892 // Gather nodal coordinates and solution vector:
893
894 double x[ FunctionCount ] ;
895 double y[ FunctionCount ] ;
896 double z[ FunctionCount ] ;
897 local_scalar_type val[ FunctionCount ] ;
898 unsigned node_index[ ElemNodeCount ];
899
900 for ( unsigned i = 0 ; i < ElemNodeCount ; ++i ) {
901 const unsigned ni = elem_node_ids( ielem , i );
902
903 node_index[i] = ni ;
904
905 x[i] = node_coords( ni , 0 );
906 y[i] = node_coords( ni , 1 );
907 z[i] = node_coords( ni , 2 );
908
909 val[i] = local_solution( ni );
910 }
911
912
913 local_scalar_type elem_vec[ FunctionCount ] ;
914 local_scalar_type elem_mat[ FunctionCount ][ FunctionCount ] ;
915
916 for( unsigned i = 0; i < FunctionCount ; i++ ) {
917 elem_vec[i] = 0 ;
918 for( unsigned j = 0; j < FunctionCount ; j++){
919 elem_mat[i][j] = 0 ;
920 }
921 }
922
923
924 for ( unsigned i = 0 ; i < IntegrationCount ; ++i ) {
925 double dpsidx[ FunctionCount ] ;
926 double dpsidy[ FunctionCount ] ;
927 double dpsidz[ FunctionCount ] ;
928
929 local_scalar_type coeff_k = 0 ;
930
931 {
932 double pt[] = {0.0, 0.0, 0.0};
933
934 // If function is not constant
935 // then compute physical coordinates of integration point
936 if ( ! coeff_function.is_constant ) {
937 for ( unsigned j = 0 ; j < FunctionCount ; ++j ) {
938 pt[0] += x[j] * elem_data.values[i][j] ;
939 pt[1] += y[j] * elem_data.values[i][j] ;
940 pt[2] += z[j] * elem_data.values[i][j] ;
941 }
942 }
943
944 // Need to fix this for local_scalar_type!!!!!!
945 coeff_k = coeff_function(pt, ensemble_rank);
946 }
947
948 const double detJ =
949 transform_gradients( elem_data.gradients[i] , x , y , z ,
950 dpsidx , dpsidy , dpsidz );
951
952 contributeResidualJacobian( val , dpsidx , dpsidy , dpsidz ,
953 detJ , coeff_k ,
954 elem_data.weights[i] ,
955 elem_data.values[i] ,
956 elem_vec , elem_mat );
957 }
958
959 for( unsigned i = 0 ; i < FunctionCount ; i++ ) {
960 const unsigned row = node_index[i] ;
961 if ( row < residual.extent(0) ) {
962 atomic_add( & local_residual( row ) , elem_vec[i] );
963
964 for( unsigned j = 0 ; j < FunctionCount ; j++ ) {
965 const unsigned entry = elem_graph( ielem , i , j );
966 if ( entry != ~0u ) {
967 atomic_add( & local_jacobian_values( entry ) , elem_mat[i][j] );
968 }
969 }
970 }
971 }
972 }
973}; /* ElementComputation */
974
975//----------------------------------------------------------------------------
976
977template< class FixtureType , class SparseMatrixType >
978class DirichletComputation ;
979
980template< class ExecutionSpace , BoxElemPart::ElemOrder Order , class CoordinateMap ,
981 typename ScalarType >
982class DirichletComputation<
983 Kokkos::Example::BoxElemFixture< ExecutionSpace , Order , CoordinateMap > ,
984 Kokkos::Example::FENL::CrsMatrix< ScalarType , ExecutionSpace > >
985{
986public:
987
990 typedef typename node_coord_type::value_type scalar_coord_type ;
991
992 typedef ExecutionSpace execution_space ;
993 typedef ScalarType scalar_type ;
994
998 typedef Kokkos::View< scalar_type* , execution_space > vector_type ;
999
1000 //------------------------------------
1001
1006 static const bool use_team = local_vector_view_traits::use_team;
1007
1008 typedef double bc_scalar_type ;
1009
1010 //------------------------------------
1011 // Computational data:
1012
1013 const node_coord_type node_coords ;
1014 const vector_type solution ;
1015 const sparse_matrix_type jacobian ;
1016 const vector_type residual ;
1017 const bc_scalar_type bc_lower_value ;
1018 const bc_scalar_type bc_upper_value ;
1019 const scalar_coord_type bc_lower_limit ;
1020 const scalar_coord_type bc_upper_limit ;
1021 const unsigned bc_plane ;
1022 const unsigned node_count ;
1023 bool init ;
1024 const Kokkos::Example::FENL::DeviceConfig dev_config ;
1025
1026
1028 const vector_type & arg_solution ,
1029 const sparse_matrix_type & arg_jacobian ,
1030 const vector_type & arg_residual ,
1031 const unsigned arg_bc_plane ,
1032 const bc_scalar_type arg_bc_lower_value ,
1033 const bc_scalar_type arg_bc_upper_value ,
1034 const Kokkos::Example::FENL::DeviceConfig arg_dev_config )
1035 : node_coords( arg_mesh.node_coord() )
1036 , solution( arg_solution )
1037 , jacobian( arg_jacobian )
1038 , residual( arg_residual )
1039 , bc_lower_value( arg_bc_lower_value )
1040 , bc_upper_value( arg_bc_upper_value )
1041 , bc_lower_limit( std::numeric_limits<scalar_coord_type>::epsilon() )
1042 , bc_upper_limit( scalar_coord_type(1) - std::numeric_limits<scalar_coord_type>::epsilon() )
1043 , bc_plane( arg_bc_plane )
1044 , node_count( arg_mesh.node_count_owned() )
1045 , init( false )
1046 , dev_config( arg_dev_config )
1047 {
1048 parallel_for( node_count , *this );
1049 init = true ;
1050 }
1051
1052 void apply() const
1053 {
1054 if ( use_team ) {
1055 const size_t team_size = dev_config.block_dim.x * dev_config.block_dim.y;
1056 const size_t league_size =
1057 (node_count + dev_config.block_dim.y-1) / dev_config.block_dim.y;
1058 Kokkos::TeamPolicy< execution_space > config( league_size, team_size );
1059 parallel_for( config , *this );
1060 }
1061 else
1062 parallel_for( node_count , *this );
1063 }
1064
1065 //------------------------------------
1066
1067 KOKKOS_INLINE_FUNCTION
1068 void operator()( const typename TeamPolicy< execution_space >::member_type & dev ) const
1069 {
1070
1071 const unsigned num_ensemble_threads = dev_config.block_dim.x ;
1072 const unsigned num_node_threads = dev_config.block_dim.y ;
1073 const unsigned node_rank = dev.team_rank() / num_ensemble_threads ;
1074 const unsigned ensemble_rank = dev.team_rank() % num_ensemble_threads ;
1075
1076 const unsigned inode = dev.league_rank() * num_node_threads + node_rank;
1077
1078 if (inode >= node_count)
1079 return;
1080
1081 (*this)( inode, ensemble_rank );
1082
1083 }
1084
1085 KOKKOS_INLINE_FUNCTION
1086 void operator()( const unsigned inode ,
1087 const unsigned ensemble_rank = 0) const
1088 {
1089 local_vector_type local_residual =
1090 local_vector_view_traits::create_local_view(residual,
1091 ensemble_rank);
1092 local_matrix_type local_jacobian_values =
1093 local_matrix_view_traits::create_local_view(jacobian.values,
1094 ensemble_rank);
1095
1096 // Apply dirichlet boundary condition on the Solution and Residual vectors.
1097 // To maintain the symmetry of the original global stiffness matrix,
1098 // zero out the columns that correspond to boundary conditions, and
1099 // update the residual vector accordingly
1100
1101 const unsigned iBeg = jacobian.graph.row_map[inode];
1102 const unsigned iEnd = jacobian.graph.row_map[inode+1];
1103
1104 const scalar_coord_type c = node_coords(inode,bc_plane);
1105 const bool bc_lower = c <= bc_lower_limit ;
1106 const bool bc_upper = bc_upper_limit <= c ;
1107
1108 if ( ! init ) {
1109 solution(inode) = bc_lower ? bc_lower_value : (
1110 bc_upper ? bc_upper_value : 0 );
1111 }
1112 else {
1113 if ( bc_lower || bc_upper ) {
1114
1115 local_residual(inode) = 0 ;
1116
1117 // zero each value on the row, and leave a one
1118 // on the diagonal
1119
1120 for( unsigned i = iBeg ; i < iEnd ; ++i ) {
1121 local_jacobian_values(i) = int(inode) == int(jacobian.graph.entries(i)) ? 1 : 0 ;
1122 }
1123 }
1124 else {
1125
1126 // Find any columns that are boundary conditions.
1127 // Clear them and adjust the residual vector
1128
1129 for( unsigned i = iBeg ; i < iEnd ; ++i ) {
1130 const unsigned cnode = jacobian.graph.entries(i) ;
1131 const scalar_coord_type cc = node_coords(cnode,bc_plane);
1132
1133 if ( ( cc <= bc_lower_limit ) || ( bc_upper_limit <= cc ) ) {
1134 local_jacobian_values(i) = 0 ;
1135 }
1136 }
1137 }
1138 }
1139 }
1140};
1141
1142template< typename FixtureType , typename VectorType >
1144{
1145public:
1146
1147 typedef FixtureType fixture_type ;
1148 typedef VectorType vector_type ;
1149 typedef typename vector_type::execution_space execution_space ;
1150 typedef typename vector_type::value_type value_type ;
1151
1154 static const unsigned TensorDim = SpatialDim * SpatialDim ;
1157
1158 //------------------------------------
1159 // Computational data:
1160
1164
1166 : elem_data()
1167 , fixture( rhs.fixture )
1168 , solution( rhs.solution )
1169 {}
1170
1171 ResponseComputation( const fixture_type& arg_fixture ,
1172 const vector_type & arg_solution )
1173 : elem_data()
1174 , fixture( arg_fixture )
1175 , solution( arg_solution )
1176 {}
1177
1178 //------------------------------------
1179
1181 {
1182 value_type response = 0;
1183 //Kokkos::parallel_reduce( fixture.elem_count() , *this , response );
1184 Kokkos::parallel_reduce( solution.extent(0) , *this , response );
1185 return response;
1186 }
1187
1188 //------------------------------------
1189
1190 KOKKOS_INLINE_FUNCTION
1192 const double grad[][ ElemNodeCount ] , // Gradient of bases master element
1193 const double x[] ,
1194 const double y[] ,
1195 const double z[] ) const
1196 {
1197 enum { j11 = 0 , j12 = 1 , j13 = 2 ,
1198 j21 = 3 , j22 = 4 , j23 = 5 ,
1199 j31 = 6 , j32 = 7 , j33 = 8 };
1200
1201 // Jacobian accumulation:
1202
1203 double J[ TensorDim ] = { 0, 0, 0, 0, 0, 0, 0, 0, 0 };
1204
1205 for( unsigned i = 0; i < ElemNodeCount ; ++i ) {
1206 const double x1 = x[i] ;
1207 const double x2 = y[i] ;
1208 const double x3 = z[i] ;
1209
1210 const double g1 = grad[0][i] ;
1211 const double g2 = grad[1][i] ;
1212 const double g3 = grad[2][i] ;
1213
1214 J[j11] += g1 * x1 ;
1215 J[j12] += g1 * x2 ;
1216 J[j13] += g1 * x3 ;
1217
1218 J[j21] += g2 * x1 ;
1219 J[j22] += g2 * x2 ;
1220 J[j23] += g2 * x3 ;
1221
1222 J[j31] += g3 * x1 ;
1223 J[j32] += g3 * x2 ;
1224 J[j33] += g3 * x3 ;
1225 }
1226
1227 // Inverse jacobian:
1228
1229 double invJ[ TensorDim ] = {
1230 static_cast<double>( J[j22] * J[j33] - J[j23] * J[j32] ) ,
1231 static_cast<double>( J[j13] * J[j32] - J[j12] * J[j33] ) ,
1232 static_cast<double>( J[j12] * J[j23] - J[j13] * J[j22] ) ,
1233
1234 static_cast<double>( J[j23] * J[j31] - J[j21] * J[j33] ) ,
1235 static_cast<double>( J[j11] * J[j33] - J[j13] * J[j31] ) ,
1236 static_cast<double>( J[j13] * J[j21] - J[j11] * J[j23] ) ,
1237
1238 static_cast<double>( J[j21] * J[j32] - J[j22] * J[j31] ) ,
1239 static_cast<double>( J[j12] * J[j31] - J[j11] * J[j32] ) ,
1240 static_cast<double>( J[j11] * J[j22] - J[j12] * J[j21] ) };
1241
1242 const double detJ = J[j11] * invJ[j11] +
1243 J[j21] * invJ[j12] +
1244 J[j31] * invJ[j13] ;
1245
1246 return detJ ;
1247 }
1248
1249 KOKKOS_INLINE_FUNCTION
1251 const value_type dof_values[] ,
1252 const double detJ ,
1253 const double integ_weight ,
1254 const double bases_vals[] ) const
1255 {
1256 // $$ g_i = \int_{\Omega} T^2 d \Omega $$
1257
1258 value_type value_at_pt = 0 ;
1259 for ( unsigned m = 0 ; m < ElemNodeCount ; m++ ) {
1260 value_at_pt += dof_values[m] * bases_vals[m] ;
1261 }
1262
1263 value_type elem_response =
1264 value_at_pt * value_at_pt * detJ * integ_weight ;
1265
1266 return elem_response;
1267 }
1268
1269 /*
1270 KOKKOS_INLINE_FUNCTION
1271 void operator()( const unsigned ielem , value_type& response ) const
1272 {
1273 // Gather nodal coordinates and solution vector:
1274
1275 double x[ ElemNodeCount ] ;
1276 double y[ ElemNodeCount ] ;
1277 double z[ ElemNodeCount ] ;
1278 value_type val[ ElemNodeCount ] ;
1279
1280 for ( unsigned i = 0 ; i < ElemNodeCount ; ++i ) {
1281 const unsigned ni = fixture.elem_node( ielem , i );
1282
1283 x[i] = fixture.node_coord( ni , 0 );
1284 y[i] = fixture.node_coord( ni , 1 );
1285 z[i] = fixture.node_coord( ni , 2 );
1286
1287 val[i] = solution( ni );
1288 }
1289
1290 for ( unsigned i = 0 ; i < IntegrationCount ; ++i ) {
1291
1292 const double detJ = compute_detJ( elem_data.gradients[i] , x , y , z );
1293
1294 response += contributeResponse( val , detJ , elem_data.weights[i] ,
1295 elem_data.values[i] );
1296 }
1297 }
1298 */
1299
1300 KOKKOS_INLINE_FUNCTION
1301 void operator()( const unsigned i , value_type& response ) const
1302 {
1303 value_type u = solution(i);
1304 response += (u * u) / fixture.node_count_global();
1305 }
1306
1307 KOKKOS_INLINE_FUNCTION
1308 void init( value_type & response ) const
1309 { response = 0 ; }
1310
1311 KOKKOS_INLINE_FUNCTION
1312 void join( value_type & response ,
1313 const value_type & input ) const
1314 { response += input ; }
1315
1316}; /* ResponseComputation */
1317
1318} /* namespace FENL */
1319} /* namespace Example */
1320} /* namespace Kokkos */
1321
1322//----------------------------------------------------------------------------
1323
1324/* A Cuda-specific specialization for the element computation functor. */
1325#if defined( __CUDACC__ )
1326// #include <NonlinearElement_Cuda.hpp>
1327#endif
1328
1329//----------------------------------------------------------------------------
1330
1331#endif /* #ifndef KOKKOS_EXAMPLE_FENLFUNCTORS_HPP */
expr val()
Generate a distributed unstructured finite element mesh from a partitioned NX*NY*NZ box of elements.
Kokkos::View< const double *[SpaceDim], Device > node_coord_type
Kokkos::View< const unsigned *[ElemNode], Device > elem_node_type
DirichletComputation(const mesh_type &arg_mesh, const vector_type &arg_solution, const sparse_matrix_type &arg_jacobian, const vector_type &arg_residual, const unsigned arg_bc_plane, const bc_scalar_type arg_bc_lower_value, const bc_scalar_type arg_bc_upper_value, const Kokkos::Example::FENL::DeviceConfig arg_dev_config)
KOKKOS_INLINE_FUNCTION void contributeResidualJacobian(const local_scalar_type dof_values[], const double dpsidx[], const double dpsidy[], const double dpsidz[], const double detJ, const local_scalar_type coeff_k, const double integ_weight, const double bases_vals[], local_scalar_type elem_res[], local_scalar_type elem_mat[][FunctionCount]) const
ElementComputation(const mesh_type &arg_mesh, const CoeffFunctionType &arg_coeff_function, const vector_type &arg_solution, const elem_graph_type &arg_elem_graph, const sparse_matrix_type &arg_jacobian, const vector_type &arg_residual, const Kokkos::Example::FENL::DeviceConfig arg_dev_config)
KOKKOS_INLINE_FUNCTION double transform_gradients(const double grad[][FunctionCount], const double x[], const double y[], const double z[], double dpsidx[], double dpsidy[], double dpsidz[]) const
local_rv_view_traits::local_view_type local_rv_view_type
KOKKOS_INLINE_FUNCTION void setRandomVariables(const RandomVariableView &rv)
LocalViewTraits< RandomVariableView > local_rv_view_traits
Stokhos::KL::ExponentialRandomField< MeshScalar, Device > rf_type
Kokkos::View< Scalar *, Kokkos::LayoutLeft, Device > RandomVariableView
KOKKOS_INLINE_FUNCTION RandomVariableView getRandomVariables() const
ExponentialKLCoefficient(const MeshScalar mean, const MeshScalar variance, const MeshScalar correlation_length, const size_type num_rv)
KOKKOS_INLINE_FUNCTION local_scalar_type operator()(const MeshScalar point[], const size_type ensemble_rank) const
ExponentialKLCoefficient(const ExponentialKLCoefficient &rhs)
local_rv_view_traits::local_value_type local_scalar_type
KOKKOS_INLINE_FUNCTION void operator()(const unsigned irow, unsigned &update, const bool final) const
KOKKOS_INLINE_FUNCTION void sort_graph_entries(const unsigned irow) const
KOKKOS_INLINE_FUNCTION void fill_elem_graph_map(const unsigned ielem) const
KOKKOS_INLINE_FUNCTION void fill_graph_entries(const unsigned iset) const
KOKKOS_INLINE_FUNCTION void init(unsigned &update) const
KOKKOS_INLINE_FUNCTION void operator()(const unsigned iwork) const
Kokkos::UnorderedMap< key_type, void, execution_space > SetType
Kokkos::View< unsigned *[ElemNode][ElemNode], execution_space > ElemGraphType
KOKKOS_INLINE_FUNCTION void fill_set(const unsigned ielem) const
ElemNodeIdView::execution_space execution_space
Kokkos::View< unsigned, execution_space > UnsignedValue
KOKKOS_INLINE_FUNCTION void join(unsigned &update, const unsigned &input) const
NodeNodeGraph(const ElemNodeIdView &arg_elem_node_id, const unsigned arg_node_count, Times &results)
CrsGraphType::row_map_type::non_const_type RowMapType
KOKKOS_INLINE_FUNCTION value_type contributeResponse(const value_type dof_values[], const double detJ, const double integ_weight, const double bases_vals[]) const
KOKKOS_INLINE_FUNCTION void join(value_type &response, const value_type &input) const
ResponseComputation(const fixture_type &arg_fixture, const vector_type &arg_solution)
KOKKOS_INLINE_FUNCTION void init(value_type &response) const
Kokkos::Example::HexElement_Data< fixture_type::ElemNode > element_data_type
KOKKOS_INLINE_FUNCTION double compute_detJ(const double grad[][ElemNodeCount], const double x[], const double y[], const double z[]) const
KOKKOS_INLINE_FUNCTION void operator()(const unsigned i, value_type &response) const
double values[integration_count][function_count]
double gradients[integration_count][spatial_dimension][function_count]
Class representing a KL expansion of an exponential random field.
KOKKOS_INLINE_FUNCTION Teuchos::PromotionTraits< typenamerv_type::value_type, value_type >::promote evaluate(const point_type &point, const rv_type &random_variables) const
Evaluate random field at a point.
KOKKOS_INLINE_FUNCTION void atomic_add(volatile Sacado::UQ::PCE< Storage > *const dest, const Sacado::UQ::PCE< Storage > &src)
void deep_copy(const Stokhos::CrsMatrix< ValueType, DstDevice, Layout > &dst, const Stokhos::CrsMatrix< ValueType, SrcDevice, Layout > &src)
static void eval(Kokkos::Example::FENL::DeviceConfig &dev_config_elem, Kokkos::Example::FENL::DeviceConfig &dev_config_bc)
static void eval(Kokkos::Example::FENL::DeviceConfig &dev_config_elem, Kokkos::Example::FENL::DeviceConfig &dev_config_bc)
Kokkos::StaticCrsGraph< unsigned, Space, void, void, unsigned > StaticCrsGraphType
CrsMatrix(const StaticCrsGraphType &arg_graph)
Dim3(const size_t x_, const size_t y_=1, const size_t z_=1)
DeviceConfig(const size_t num_blocks_=0, const size_t threads_per_block_x_=0, const size_t threads_per_block_y_=0, const size_t threads_per_block_z_=1)
ElementComputationConstantCoefficient(const ElementComputationConstantCoefficient &rhs)
KOKKOS_INLINE_FUNCTION double operator()(double pt[], unsigned ensemble_rank) const
static KOKKOS_INLINE_FUNCTION local_view_type create_local_view(const view_type &v, const unsigned local_rank)