44 #ifndef KOKKOS_EXAMPLE_FENLFUNCTORS_HPP 45 #define KOKKOS_EXAMPLE_FENLFUNCTORS_HPP 56 #include <Kokkos_Core.hpp> 57 #include <Kokkos_Pair.hpp> 58 #include <Kokkos_UnorderedMap.hpp> 59 #include <Kokkos_StaticCrsGraph.hpp> 61 #include <impl/Kokkos_Timer.hpp> 75 template<
typename ValueType ,
class Space >
87 ,
coeff(
"crs_matrix_coeff" , arg_graph.entries.dimension_0() )
91 template<
class ElemNodeIdView ,
class CrsGraphType ,
unsigned ElemNode >
98 typedef Kokkos::UnorderedMap< key_type, void , execution_space >
SetType ;
99 typedef typename CrsGraphType::row_map_type::non_const_type
RowMapType ;
103 typedef Kokkos::View< unsigned*[ElemNode][ElemNode] , execution_space >
144 const unsigned arg_node_count,
160 Kokkos::Impl::Timer wall_clock ;
166 size_t set_capacity = (28ull *
node_count) / 2;
167 unsigned failed_insert_count = 0 ;
178 Kokkos::parallel_reduce( Kokkos::RangePolicy<execution_space,TagFillNodeSet>(0,
elem_node_id.dimension_0())
180 , failed_insert_count );
182 }
while ( failed_insert_count );
184 execution_space::fence();
200 unsigned graph_entry_count = 0 ;
202 Kokkos::deep_copy( graph_entry_count ,
row_total );
206 graph.entries =
typename CrsGraphType::entries_type(
"graph_entries" , graph_entry_count );
211 execution_space::fence();
218 execution_space::fence();
235 execution_space::fence();
243 Kokkos::parallel_for(
elem_node_id.dimension_0() , *this );
245 execution_space::fence();
256 for (
unsigned row_local_node = 0 ; row_local_node <
elem_node_id.dimension_1() ; ++row_local_node ) {
258 const unsigned row_node =
elem_node_id( ielem , row_local_node );
260 for (
unsigned col_local_node = row_local_node ; col_local_node <
elem_node_id.dimension_1() ; ++col_local_node ) {
262 const unsigned col_node =
elem_node_id( ielem , col_local_node );
268 const key_type key = (row_node < col_node) ? make_pair( row_node, col_node ) : make_pair( col_node, row_node ) ;
270 const typename SetType::insert_result result =
node_node_set.insert( key );
273 if ( result.success() ) {
276 if ( row_node <
row_count.dimension_0() ) { atomic_fetch_add( &
row_count( row_node ) , 1 ); }
279 if ( col_node <
row_count.dimension_0() && col_node != row_node ) { atomic_fetch_add( &
row_count( col_node ) , 1 ); }
281 else if ( result.failed() ) {
296 const unsigned row_node = key.first ;
297 const unsigned col_node = key.second ;
299 if ( row_node <
row_count.dimension_0() ) {
300 const unsigned offset =
graph.row_map( row_node ) + atomic_fetch_add( &
row_count( row_node ) , 1 );
301 graph.entries( offset ) = col_node ;
304 if ( col_node <
row_count.dimension_0() && col_node != row_node ) {
305 const unsigned offset =
graph.row_map( col_node ) + atomic_fetch_add( &
row_count( col_node ) , 1 );
306 graph.entries( offset ) = row_node ;
314 const unsigned row_beg =
graph.row_map( irow );
315 const unsigned row_end =
graph.row_map( irow + 1 );
316 for (
unsigned i = row_beg + 1 ; i < row_end ; ++i ) {
317 const unsigned col =
graph.entries(i);
319 for ( ; row_beg < j && col <
graph.entries(j-1) ; --j ) {
322 graph.entries(j) = col ;
329 for (
unsigned row_local_node = 0 ; row_local_node <
elem_node_id.dimension_1() ; ++row_local_node ) {
331 const unsigned row_node =
elem_node_id( ielem , row_local_node );
333 for (
unsigned col_local_node = 0 ; col_local_node <
elem_node_id.dimension_1() ; ++col_local_node ) {
335 const unsigned col_node =
elem_node_id( ielem , col_local_node );
337 unsigned entry = ~0u ;
339 if ( row_node + 1 <
graph.row_map.dimension_0() ) {
341 const unsigned entry_end =
graph.row_map( row_node + 1 );
343 entry =
graph.row_map( row_node );
345 for ( ; entry < entry_end &&
graph.entries(entry) != col_node ; ++entry );
347 if ( entry == entry_end ) entry = ~0u ;
350 elem_graph( ielem , row_local_node , col_local_node ) = entry ;
380 void operator()(
const unsigned irow ,
unsigned & update ,
const bool final )
const 383 if (
final ) {
row_map( irow ) = update ; }
388 if ( irow + 1 ==
row_count.dimension_0() ) {
401 ,
volatile unsigned & update
402 ,
volatile const unsigned & input )
const { update += input ; }
406 void init(
unsigned & update )
const { update = 0 ; }
409 void join(
volatile unsigned & update
410 ,
volatile const unsigned & input )
const { update += input ; }
427 class CoordinateMap ,
typename ScalarType >
442 typedef Kokkos::View< scalar_type* , Kokkos::LayoutLeft, execution_space >
vector_type ;
456 typedef Kokkos::View< scalar_type*[FunctionCount][FunctionCount] , execution_space >
elem_matrices_type ;
515 double dpsidz[] )
const 517 enum { j11 = 0 , j12 = 1 , j13 = 2 ,
518 j21 = 3 , j22 = 4 , j23 = 5 ,
519 j31 = 6 , j32 = 7 , j33 = 8 };
523 double J[
TensorDim ] = { 0, 0, 0, 0, 0, 0, 0, 0, 0 };
526 const double x1 = x[i] ;
527 const double x2 = y[i] ;
528 const double x3 = z[i] ;
530 const double g1 = grad[0][i] ;
531 const double g2 = grad[1][i] ;
532 const double g3 = grad[2][i] ;
550 static_cast<double>( J[j22] * J[j33] - J[j23] * J[j32] ) ,
551 static_cast<double>( J[j13] * J[j32] - J[j12] * J[j33] ) ,
552 static_cast<double>( J[j12] * J[j23] - J[j13] * J[j22] ) ,
554 static_cast<double>( J[j23] * J[j31] - J[j21] * J[j33] ) ,
555 static_cast<double>( J[j11] * J[j33] - J[j13] * J[j31] ) ,
556 static_cast<double>( J[j13] * J[j21] - J[j11] * J[j23] ) ,
558 static_cast<double>( J[j21] * J[j32] - J[j22] * J[j31] ) ,
559 static_cast<double>( J[j12] * J[j31] - J[j11] * J[j32] ) ,
560 static_cast<double>( J[j11] * J[j22] - J[j12] * J[j21] ) };
562 const double detJ = J[j11] * invJ[j11] +
566 const double detJinv = 1.0 / detJ ;
568 for (
unsigned i = 0 ; i <
TensorDim ; ++i ) { invJ[i] *= detJinv ; }
573 const double g0 = grad[0][i];
574 const double g1 = grad[1][i];
575 const double g2 = grad[2][i];
577 dpsidx[i] = g0 * invJ[j11] + g1 * invJ[j12] + g2 * invJ[j13];
578 dpsidy[i] = g0 * invJ[j21] + g1 * invJ[j22] + g2 * invJ[j23];
579 dpsidz[i] = g0 * invJ[j31] + g1 * invJ[j32] + g2 * invJ[j33];
594 template<
class FiniteElementMeshType ,
595 class SparseMatrixType ,
601 class CoordinateMap ,
typename ScalarType >
604 CrsMatrix< ScalarType , ExecutionSpace > ,
616 static const unsigned FunctionCount = base_type::FunctionCount;
617 static const unsigned IntegrationCount = base_type::IntegrationCount;
618 static const unsigned ElemNodeCount = base_type::ElemNodeCount;
628 base_type(arg_mesh, arg_solution, arg_elem_graph,
629 arg_jacobian, arg_residual) {}
635 const size_t nelem = this->elem_node_ids.dimension_0();
636 parallel_for( nelem , *
this );
642 unsigned node_index[],
643 double x[],
double y[],
double z[],
647 for (
unsigned i = 0 ; i < ElemNodeCount ; ++i ) {
648 const unsigned ni = this->elem_node_ids( ielem , i );
652 x[i] = this->node_coords( ni , 0 );
653 y[i] = this->node_coords( ni , 1 );
654 z[i] = this->node_coords( ni , 2 );
656 val[i] = this->solution( ni ) ;
659 for(
unsigned j = 0; j < FunctionCount ; j++){
667 const unsigned node_index[],
671 for(
unsigned i = 0 ; i < FunctionCount ; i++ ) {
672 const unsigned row = node_index[i] ;
673 if ( row < this->residual.dimension_0() ) {
674 atomic_add( & this->residual( row ) , res[i] );
676 for(
unsigned j = 0 ; j < FunctionCount ; j++ ) {
677 const unsigned entry = this->elem_graph( ielem , i , j );
678 if ( entry != ~0u ) {
679 atomic_add( & this->jacobian.coeff( entry ) , mat[i][j] );
695 double coeff_k = 3.456;
696 double coeff_src = 1.234;
697 double advection[] = { 1.1, 1.2, 1.3 };
698 double dpsidx[ FunctionCount ] ;
699 double dpsidy[ FunctionCount ] ;
700 double dpsidz[ FunctionCount ] ;
701 for (
unsigned i = 0 ; i < IntegrationCount ; ++i ) {
703 const double integ_weight = this->elem_data.weights[i];
704 const double* bases_vals = this->elem_data.values[i];
706 this->transform_gradients( this->elem_data.gradients[i] ,
708 dpsidx , dpsidy , dpsidz );
709 const double detJ_weight = detJ * integ_weight;
710 const double detJ_weight_coeff_k = detJ_weight * coeff_k;
716 for (
unsigned m = 0 ; m < FunctionCount ; m++ ) {
717 value_at_pt += dof_values[m] * bases_vals[m] ;
718 gradx_at_pt += dof_values[m] * dpsidx[m] ;
719 grady_at_pt += dof_values[m] * dpsidy[m] ;
720 gradz_at_pt += dof_values[m] * dpsidz[m] ;
724 coeff_src * value_at_pt * value_at_pt ;
726 2.0 * coeff_src * value_at_pt ;
733 advection_x*gradx_at_pt +
734 advection_y*grady_at_pt +
735 advection_z*gradz_at_pt ;
737 for (
unsigned m = 0; m < FunctionCount; ++m) {
739 const double bases_val_m = bases_vals[m] * detJ_weight ;
740 const double dpsidx_m = dpsidx[m] ;
741 const double dpsidy_m = dpsidy[m] ;
742 const double dpsidz_m = dpsidz[m] ;
745 detJ_weight_coeff_k * ( dpsidx_m * gradx_at_pt +
746 dpsidy_m * grady_at_pt +
747 dpsidz_m * gradz_at_pt ) +
748 bases_val_m * ( advection_term + source_term ) ;
750 for(
unsigned n = 0;
n < FunctionCount;
n++) {
751 const double dpsidx_n = dpsidx[
n] ;
752 const double dpsidy_n = dpsidy[
n] ;
753 const double dpsidz_n = dpsidz[
n] ;
755 detJ_weight_coeff_k * ( dpsidx_m * dpsidx_n +
756 dpsidy_m * dpsidy_n +
757 dpsidz_m * dpsidz_n ) +
758 bases_val_m * ( advection_x * dpsidx_n +
759 advection_y * dpsidy_n +
760 advection_z * dpsidz_n +
761 source_deriv * bases_vals[
n] ) ;
770 double x[ FunctionCount ] ;
771 double y[ FunctionCount ] ;
772 double z[ FunctionCount ] ;
773 unsigned node_index[ ElemNodeCount ];
777 scalar_type elem_mat[ FunctionCount ][ FunctionCount ] ;
780 gatherSolution(ielem,
val, node_index, x, y, z, elem_res, elem_mat);
783 computeElementResidualJacobian(
val, x, y, z, elem_res , elem_mat );
786 scatterResidual( ielem, node_index, elem_res, elem_mat );
791 class CoordinateMap ,
typename ScalarType >
794 CrsMatrix< ScalarType , ExecutionSpace > ,
805 static const unsigned FunctionCount = base_type::FunctionCount;
806 static const unsigned IntegrationCount = base_type::IntegrationCount;
807 static const unsigned ElemNodeCount = base_type::ElemNodeCount;
819 base_type(arg_mesh, arg_solution, arg_elem_graph,
820 arg_jacobian, arg_residual) {}
826 const size_t nelem = this->elem_node_ids.dimension_0();
827 parallel_for( nelem , *
this );
833 unsigned node_index[],
834 double x[],
double y[],
double z[],
837 for (
unsigned i = 0 ; i < ElemNodeCount ; ++i ) {
838 const unsigned ni = this->elem_node_ids( ielem , i );
842 x[i] = this->node_coords( ni , 0 );
843 y[i] = this->node_coords( ni , 1 );
844 z[i] = this->node_coords( ni , 2 );
846 val[i].val() = this->solution( ni );
847 val[i].diff( i, FunctionCount );
853 const unsigned node_index[],
856 for(
unsigned i = 0 ; i < FunctionCount ; i++ ) {
857 const unsigned row = node_index[i] ;
858 if ( row < this->residual.dimension_0() ) {
859 atomic_add( & this->residual( row ) , res[i].
val() );
861 for(
unsigned j = 0 ; j < FunctionCount ; j++ ) {
862 const unsigned entry = this->elem_graph( ielem , i , j );
863 if ( entry != ~0u ) {
864 atomic_add( & this->jacobian.coeff( entry ) ,
865 res[i].fastAccessDx(j) );
879 double coeff_k = 3.456;
880 double coeff_src = 1.234;
881 double advection[] = { 1.1, 1.2, 1.3 };
882 double dpsidx[ FunctionCount ] ;
883 double dpsidy[ FunctionCount ] ;
884 double dpsidz[ FunctionCount ] ;
885 for (
unsigned i = 0 ; i < IntegrationCount ; ++i ) {
887 const double integ_weight = this->elem_data.weights[i];
888 const double* bases_vals = this->elem_data.values[i];
890 this->transform_gradients( this->elem_data.gradients[i] ,
892 dpsidx , dpsidy , dpsidz );
893 const double detJ_weight = detJ * integ_weight;
894 const double detJ_weight_coeff_k = detJ_weight * coeff_k;
900 for (
unsigned m = 0 ; m < FunctionCount ; m++ ) {
901 value_at_pt += dof_values[m] * bases_vals[m] ;
902 gradx_at_pt += dof_values[m] * dpsidx[m] ;
903 grady_at_pt += dof_values[m] * dpsidy[m] ;
904 gradz_at_pt += dof_values[m] * dpsidz[m] ;
908 coeff_src * value_at_pt * value_at_pt ;
911 advection[0]*gradx_at_pt +
912 advection[1]*grady_at_pt +
913 advection[2]*gradz_at_pt;
915 for (
unsigned m = 0; m < FunctionCount; ++m) {
916 const double bases_val_m = bases_vals[m] * detJ_weight ;
917 const double dpsidx_m = dpsidx[m] ;
918 const double dpsidy_m = dpsidy[m] ;
919 const double dpsidz_m = dpsidz[m] ;
922 detJ_weight_coeff_k * ( dpsidx_m * gradx_at_pt +
923 dpsidy_m * grady_at_pt +
924 dpsidz_m * gradz_at_pt ) +
925 bases_val_m * ( advection_term + source_term ) ;
933 double x[ FunctionCount ] ;
934 double y[ FunctionCount ] ;
935 double z[ FunctionCount ] ;
936 unsigned node_index[ ElemNodeCount ];
942 gatherSolution( ielem,
val, node_index, x, y, z, elem_res );
945 computeElementResidual(
val, x, y, z, elem_res );
948 scatterResidual( ielem, node_index, elem_res );
953 class CoordinateMap ,
typename ScalarType >
956 CrsMatrix< ScalarType , ExecutionSpace > ,
958 public ElementComputation< Kokkos::Example::BoxElemFixture< ExecutionSpace , Order , CoordinateMap > ,
959 CrsMatrix< ScalarType , ExecutionSpace > ,
970 static const unsigned FunctionCount = base_type::FunctionCount;
971 static const unsigned IntegrationCount = base_type::IntegrationCount;
972 static const unsigned ElemNodeCount = base_type::ElemNodeCount;
984 base_type(arg_mesh, arg_solution, arg_elem_graph,
985 arg_jacobian, arg_residual) {}
991 const size_t nelem = this->elem_node_ids.dimension_0();
992 parallel_for( nelem , *
this );
998 unsigned node_index[],
999 double x[],
double y[],
double z[],
1002 for (
unsigned i = 0 ; i < ElemNodeCount ; ++i ) {
1003 const unsigned ni = this->elem_node_ids( ielem , i );
1005 node_index[i] = ni ;
1007 x[i] = this->node_coords( ni , 0 );
1008 y[i] = this->node_coords( ni , 1 );
1009 z[i] = this->node_coords( ni , 2 );
1011 val[i] = this->solution( ni );
1022 double coeff_k = 3.456;
1023 double coeff_src = 1.234;
1024 double advection[] = { 1.1, 1.2, 1.3 };
1025 double dpsidx[ FunctionCount ] ;
1026 double dpsidy[ FunctionCount ] ;
1027 double dpsidz[ FunctionCount ] ;
1028 for (
unsigned i = 0 ; i < IntegrationCount ; ++i ) {
1030 const double integ_weight = this->elem_data.weights[i];
1031 const double* bases_vals = this->elem_data.values[i];
1033 this->transform_gradients( this->elem_data.gradients[i] ,
1035 dpsidx , dpsidy , dpsidz );
1036 const double detJ_weight = detJ * integ_weight;
1037 const double detJ_weight_coeff_k = detJ_weight * coeff_k;
1043 for (
unsigned m = 0 ; m < FunctionCount ; m++ ) {
1044 value_at_pt.val() += dof_values[m] * bases_vals[m] ;
1045 value_at_pt.fastAccessDx(m) = bases_vals[m] ;
1047 gradx_at_pt.val() += dof_values[m] * dpsidx[m] ;
1048 gradx_at_pt.fastAccessDx(m) = dpsidx[m] ;
1050 grady_at_pt.val() += dof_values[m] * dpsidy[m] ;
1051 grady_at_pt.fastAccessDx(m) = dpsidy[m] ;
1053 gradz_at_pt.val() += dof_values[m] * dpsidz[m] ;
1054 gradz_at_pt.fastAccessDx(m) = dpsidz[m] ;
1058 coeff_src * value_at_pt * value_at_pt ;
1061 advection[0]*gradx_at_pt +
1062 advection[1]*grady_at_pt +
1063 advection[2]*gradz_at_pt;
1065 for (
unsigned m = 0; m < FunctionCount; ++m) {
1066 const double bases_val_m = bases_vals[m] * detJ_weight ;
1067 const double dpsidx_m = dpsidx[m] ;
1068 const double dpsidy_m = dpsidy[m] ;
1069 const double dpsidz_m = dpsidz[m] ;
1072 detJ_weight_coeff_k * ( dpsidx_m * gradx_at_pt +
1073 dpsidy_m * grady_at_pt +
1074 dpsidz_m * gradz_at_pt ) +
1075 bases_val_m * ( advection_term + source_term ) ;
1083 double x[ FunctionCount ] ;
1084 double y[ FunctionCount ] ;
1085 double z[ FunctionCount ] ;
1086 unsigned node_index[ ElemNodeCount ];
1092 gatherSolution( ielem,
val, node_index, x, y, z, elem_res );
1095 computeElementResidual(
val, x, y, z, elem_res );
1098 this->scatterResidual( ielem, node_index, elem_res );
1103 class CoordinateMap ,
typename ScalarType >
1106 CrsMatrix< ScalarType , ExecutionSpace > ,
1108 public ElementComputation< Kokkos::Example::BoxElemFixture< ExecutionSpace , Order , CoordinateMap > ,
1109 CrsMatrix< ScalarType , ExecutionSpace > ,
1120 static const unsigned FunctionCount = base_type::FunctionCount;
1121 static const unsigned IntegrationCount = base_type::IntegrationCount;
1122 static const unsigned ElemNodeCount = base_type::ElemNodeCount;
1134 base_type(arg_mesh, arg_solution, arg_elem_graph,
1135 arg_jacobian, arg_residual) {}
1141 const size_t nelem = this->elem_node_ids.dimension_0();
1142 parallel_for( nelem , *
this );
1154 double coeff_k = 3.456;
1155 double coeff_src = 1.234;
1156 double advection[] = { 1.1, 1.2, 1.3 };
1157 double dpsidx[ FunctionCount ] ;
1158 double dpsidy[ FunctionCount ] ;
1159 double dpsidz[ FunctionCount ] ;
1165 for (
unsigned i = 0 ; i < IntegrationCount ; ++i ) {
1167 const double integ_weight = this->elem_data.weights[i];
1168 const double* bases_vals = this->elem_data.values[i];
1170 this->transform_gradients( this->elem_data.gradients[i] ,
1172 dpsidx , dpsidy , dpsidz );
1173 const double detJ_weight = detJ * integ_weight;
1174 const double detJ_weight_coeff_k = detJ_weight * coeff_k;
1176 value_at_pt.val() = 0.0 ;
1177 gradx_at_pt.val() = 0.0 ;
1178 grady_at_pt.val() = 0.0 ;
1179 gradz_at_pt.val() = 0.0 ;
1180 for (
unsigned m = 0 ; m < FunctionCount ; m++ ) {
1181 value_at_pt.val() += dof_values[m] * bases_vals[m] ;
1182 gradx_at_pt.val() += dof_values[m] * dpsidx[m] ;
1183 grady_at_pt.val() += dof_values[m] * dpsidy[m] ;
1184 gradz_at_pt.val() += dof_values[m] * dpsidz[m] ;
1188 coeff_src * value_at_pt * value_at_pt ;
1191 advection[0]*gradx_at_pt +
1192 advection[1]*grady_at_pt +
1193 advection[2]*gradz_at_pt;
1195 for (
unsigned m = 0; m < FunctionCount; ++m) {
1196 const double bases_val_m = bases_vals[m] * detJ_weight ;
1198 detJ_weight_coeff_k * ( dpsidx[m] * gradx_at_pt +
1199 dpsidy[m] * grady_at_pt +
1200 dpsidz[m] * gradz_at_pt ) +
1201 bases_val_m * ( advection_term + source_term ) ;
1203 elem_res[m] += res.val();
1206 for(
unsigned n = 0;
n < FunctionCount;
n++) {
1207 mat[
n] += res.fastAccessDx(0) * bases_vals[
n] +
1208 res.fastAccessDx(1) * dpsidx[
n] +
1209 res.fastAccessDx(2) * dpsidy[
n] +
1210 res.fastAccessDx(3) * dpsidz[
n];
1219 double x[ FunctionCount ] ;
1220 double y[ FunctionCount ] ;
1221 double z[ FunctionCount ] ;
1222 unsigned node_index[ ElemNodeCount ];
1226 scalar_type elem_mat[ FunctionCount ][ FunctionCount ] ;
1229 this->gatherSolution( ielem,
val, node_index, x, y, z, elem_res, elem_mat );
1232 computeElementResidualJacobian(
val, x, y, z, elem_res, elem_mat );
1235 this->scatterResidual( ielem, node_index, elem_res, elem_mat );
KOKKOS_INLINE_FUNCTION void operator()(const unsigned iwork) const
KOKKOS_INLINE_FUNCTION void operator()(const unsigned ielem) const
KOKKOS_INLINE_FUNCTION void scatterResidual(const unsigned ielem, const unsigned node_index[], const scalar_type res[], const scalar_type mat[][FunctionCount]) const
ElementComputation(const typename base_type::mesh_type &arg_mesh, const typename base_type::vector_type &arg_solution, const typename base_type::elem_graph_type &arg_elem_graph, const typename base_type::sparse_matrix_type &arg_jacobian, const typename base_type::vector_type &arg_residual)
KOKKOS_INLINE_FUNCTION void join(volatile unsigned &update, volatile const unsigned &input) const
static const unsigned TensorDim
base_type::scalar_type scalar_type
base_type::execution_space execution_space
const elem_vectors_type elem_residuals
KOKKOS_INLINE_FUNCTION void computeElementResidual(const scalar_type dof_values[], const double x[], const double y[], const double z[], fad_scalar_type elem_res[]) const
static const unsigned SpatialDim
Kokkos::View< scalar_type *[FunctionCount], execution_space > elem_vectors_type
static const unsigned element_node_count
CrsMatrix< ScalarType, ExecutionSpace > sparse_matrix_type
KOKKOS_INLINE_FUNCTION void gatherSolution(const unsigned ielem, scalar_type val[], unsigned node_index[], double x[], double y[], double z[], fad_scalar_type res[]) const
const elem_matrices_type elem_jacobians
static const unsigned function_count
KOKKOS_INLINE_FUNCTION void operator()(const unsigned ielem) const
ElementComputation< Kokkos::Example::BoxElemFixture< ExecutionSpace, Order, CoordinateMap >, CrsMatrix< ScalarType, ExecutionSpace >, Analytic > base_type
Kokkos::View< const double *[SpaceDim], Device > node_coord_type
Kokkos::View< const unsigned *[ElemNode], Device > elem_node_type
ElementComputationBase< ExecutionSpace, Order, CoordinateMap, ScalarType > base_type
ExecutionSpace execution_space
Kokkos::View< scalar_type *, Kokkos::LayoutLeft, execution_space > vector_type
base_type::execution_space execution_space
const vector_type residual
Sacado::Fad::SFad< scalar_type, 4 > fad_scalar_type
static const unsigned IntegrationCount
const elem_node_type elem_node_ids
KOKKOS_INLINE_FUNCTION void init(const TagFillNodeSet &, unsigned &update) const
const unsigned node_count
KOKKOS_INLINE_FUNCTION void join(const TagFillNodeSet &, volatile unsigned &update, volatile const unsigned &input) const
KOKKOS_INLINE_FUNCTION void fill_graph_entries(const unsigned iset) const
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
NodeNodeGraph< elem_node_type, sparse_graph_type, ElemNodeCount >::ElemGraphType elem_graph_type
const vector_type solution
base_type::execution_space execution_space
#define KOKKOS_INLINE_FUNCTION
ElemNodeIdView::execution_space execution_space
KOKKOS_INLINE_FUNCTION void computeElementResidualJacobian(const scalar_type dof_values[], const double x[], const double y[], const double z[], scalar_type elem_res[], scalar_type elem_mat[][FunctionCount]) const
KOKKOS_INLINE_FUNCTION void operator()(const unsigned ielem) const
ElementComputationBase< ExecutionSpace, Order, CoordinateMap, ScalarType > base_type
double sort_graph_entries
double fill_graph_entries
KOKKOS_INLINE_FUNCTION void operator()(const TagFillNodeSet &, unsigned ielem, unsigned &count) const
ElementComputationBase(const ElementComputationBase &rhs)
KOKKOS_INLINE_FUNCTION void operator()(const unsigned irow, unsigned &update, const bool final) const
Kokkos::Example::HexElement_Data< mesh_type::ElemNode > element_data_type
CrsGraphType::row_map_type::non_const_type RowMapType
KOKKOS_INLINE_FUNCTION void scatterResidual(const unsigned ielem, const unsigned node_index[], fad_scalar_type res[]) const
const element_data_type elem_data
KOKKOS_INLINE_FUNCTION void sort_graph_entries(const unsigned irow) const
ElementComputationBase(const mesh_type &arg_mesh, const vector_type &arg_solution, const elem_graph_type &arg_elem_graph, const sparse_matrix_type &arg_jacobian, const vector_type &arg_residual)
ElementComputation(const ElementComputation &rhs)
CrsMatrix(const StaticCrsGraphType &arg_graph)
mesh_type::elem_node_type elem_node_type
Do not initialize the derivative array.
Kokkos::StaticCrsGraph< unsigned, Space, void, unsigned > StaticCrsGraphType
KOKKOS_INLINE_FUNCTION void gatherSolution(const unsigned ielem, scalar_type val[], unsigned node_index[], double x[], double y[], double z[], scalar_type res[], scalar_type mat[][FunctionCount]) const
static const unsigned ElemNodeCount
const ElemNodeIdView elem_node_id
static const unsigned integration_count
mesh_type::node_coord_type node_coord_type
const elem_graph_type elem_graph
ElementComputation(const ElementComputation &rhs)
KOKKOS_INLINE_FUNCTION void operator()(const unsigned ielem) const
Kokkos::Example::BoxElemFixture< ExecutionSpace, Order, CoordinateMap > mesh_type
NodeNodeGraph(const ElemNodeIdView &arg_elem_node_id, const unsigned arg_node_count, Times &results)
const sparse_matrix_type jacobian
const node_coord_type node_coords
base_type::execution_space execution_space
Sacado::Fad::SFad< scalar_type, FunctionCount > fad_scalar_type
KOKKOS_INLINE_FUNCTION void fill_elem_graph_map(const unsigned ielem) const
View< ValueType *, Space > coeff_type
sparse_matrix_type::StaticCrsGraphType sparse_graph_type
Kokkos::View< unsigned, execution_space > UnsignedValue
ElementComputation(const typename base_type::mesh_type &arg_mesh, const typename base_type::vector_type &arg_solution, const typename base_type::elem_graph_type &arg_elem_graph, const typename base_type::sparse_matrix_type &arg_jacobian, const typename base_type::vector_type &arg_residual)
ElementComputation< Kokkos::Example::BoxElemFixture< ExecutionSpace, Order, CoordinateMap >, CrsMatrix< ScalarType, ExecutionSpace >, FadElement > base_type
KOKKOS_INLINE_FUNCTION void computeElementResidualJacobian(const scalar_type dof_values[], const double x[], const double y[], const double z[], scalar_type elem_res[], scalar_type elem_mat[][FunctionCount]) const
ElementComputation(const typename base_type::mesh_type &arg_mesh, const typename base_type::vector_type &arg_solution, const typename base_type::elem_graph_type &arg_elem_graph, const typename base_type::sparse_matrix_type &arg_jacobian, const typename base_type::vector_type &arg_residual)
ElementComputation(const ElementComputation &rhs)
KOKKOS_INLINE_FUNCTION void computeElementResidual(const fad_scalar_type dof_values[], const double x[], const double y[], const double z[], fad_scalar_type elem_res[]) const
ElementComputation(const ElementComputation &rhs)
ElementComputation(const typename base_type::mesh_type &arg_mesh, const typename base_type::vector_type &arg_solution, const typename base_type::elem_graph_type &arg_elem_graph, const typename base_type::sparse_matrix_type &arg_jacobian, const typename base_type::vector_type &arg_residual)
base_type::scalar_type scalar_type
base_type::scalar_type scalar_type
double fill_element_graph
Kokkos::View< scalar_type *[FunctionCount][FunctionCount], execution_space > elem_matrices_type
Sacado::Fad::SFad< scalar_type, FunctionCount > fad_scalar_type
Kokkos::View< unsigned *[ElemNode][ElemNode], execution_space > ElemGraphType
static const unsigned FunctionCount
Generate a distributed unstructured finite element mesh from a partitioned NX*NY*NZ box of elements...
KOKKOS_INLINE_FUNCTION void gatherSolution(const unsigned ielem, fad_scalar_type val[], unsigned node_index[], double x[], double y[], double z[], fad_scalar_type res[]) const
KOKKOS_INLINE_FUNCTION void init(unsigned &update) const
base_type::scalar_type scalar_type
Kokkos::UnorderedMap< key_type, void, execution_space > SetType
pair< unsigned, unsigned > key_type
static const unsigned spatial_dimension