1 #ifndef _COMPADRE_GMLS_HPP_ 2 #define _COMPADRE_GMLS_HPP_ 4 #include "Compadre_Config.h" 35 Kokkos::View<double*>
_w;
38 Kokkos::View<double*>
_P;
41 Kokkos::View<double*>
_RHS;
46 Kokkos::View<double*>
_T;
53 Kokkos::View<double*>::HostMirror
_host_T;
92 Kokkos::View<double*, layout_right>
_alphas;
95 Kokkos::View<const double*, layout_right>::HostMirror
_host_alphas;
226 std::vector<TargetOperation>
_lro;
303 KOKKOS_INLINE_FUNCTION
304 void calcPij(
const member_type& teamMember,
double* delta,
double* thread_workspace,
const int target_index,
int neighbor_index,
const double alpha,
const int dimension,
const int poly_order,
bool specific_order_only =
false,
const scratch_matrix_right_type* V = NULL,
const ReconstructionSpace reconstruction_space =
ReconstructionSpace::ScalarTaylorPolynomial,
const SamplingFunctional sampling_strategy =
PointSample,
const int additional_evaluation_local_index = 0)
const;
321 KOKKOS_INLINE_FUNCTION
322 void calcGradientPij(
const member_type& teamMember,
double* delta,
double* thread_workspace,
const int target_index,
int neighbor_index,
const double alpha,
const int partial_direction,
const int dimension,
const int poly_order,
bool specific_order_only,
const scratch_matrix_right_type* V,
const ReconstructionSpace reconstruction_space,
const SamplingFunctional sampling_strategy,
const int additional_evaluation_local_index = 0)
const;
340 KOKKOS_INLINE_FUNCTION
341 void calcHessianPij(
const member_type& teamMember,
double* delta,
double* thread_workspace,
const int target_index,
int neighbor_index,
const double alpha,
const int partial_direction_1,
const int partial_direction_2,
const int dimension,
const int poly_order,
bool specific_order_only,
const scratch_matrix_right_type* V,
const ReconstructionSpace reconstruction_space,
const SamplingFunctional sampling_strategy,
const int additional_evaluation_local_index = 0)
const;
356 KOKKOS_INLINE_FUNCTION
357 void createWeightsAndP(
const member_type& teamMember,
scratch_vector_type delta,
scratch_vector_type thread_workspace,
scratch_matrix_right_type P,
scratch_vector_type w,
const int dimension,
int polynomial_order,
bool weight_p =
false,
scratch_matrix_right_type* V = NULL,
const ReconstructionSpace reconstruction_space =
ReconstructionSpace::ScalarTaylorPolynomial,
const SamplingFunctional sampling_strategy =
PointSample)
const;
372 KOKKOS_INLINE_FUNCTION
381 KOKKOS_INLINE_FUNCTION
394 KOKKOS_INLINE_FUNCTION
410 KOKKOS_INLINE_FUNCTION
414 KOKKOS_INLINE_FUNCTION
427 KOKKOS_INLINE_FUNCTION
434 KOKKOS_INLINE_FUNCTION
440 KOKKOS_INLINE_FUNCTION
446 KOKKOS_INLINE_FUNCTION
453 KOKKOS_INLINE_FUNCTION
461 KOKKOS_INLINE_FUNCTION
464 &&
"additional_list_num must be greater than or equal to 1, unlike neighbor lists which begin indexing at 0.");
470 KOKKOS_INLINE_FUNCTION
486 KOKKOS_INLINE_FUNCTION
502 KOKKOS_INLINE_FUNCTION
515 KOKKOS_INLINE_FUNCTION
517 XYZ coordinate_delta;
523 return coordinate_delta;
527 KOKKOS_INLINE_FUNCTION
531 val += global_coord.
x * (*V)(dim, 0);
532 val += global_coord.
y * (*V)(dim, 1);
533 if (
_dimensions>2) val += global_coord.
z * (*V)(dim, 2);
538 KOKKOS_INLINE_FUNCTION
543 val = local_coord.
x * (*V)(0, dim);
544 }
else if (dim == 0) {
545 val = local_coord.
x * ((*V)(0, dim) + (*V)(1, dim));
546 }
else if (dim == 1) {
547 val = local_coord.
y * ((*V)(0, dim) + (*V)(1, dim));
553 int getTargetOffsetIndexHost(
const int lro_num,
const int input_component,
const int output_component,
const int additional_evaluation_local_index = 0)
const {
557 + output_component );
561 KOKKOS_INLINE_FUNCTION
562 int getTargetOffsetIndexDevice(
const int lro_num,
const int input_component,
const int output_component,
const int additional_evaluation_local_index = 0)
const {
566 + output_component );
578 std::string solver_type_to_lower = dense_solver_type;
579 transform(solver_type_to_lower.begin(), solver_type_to_lower.end(), solver_type_to_lower.begin(), ::tolower);
580 if (solver_type_to_lower ==
"lu") {
589 std::string problem_type_to_lower = problem_type;
590 transform(problem_type_to_lower.begin(), problem_type_to_lower.end(), problem_type_to_lower.begin(), ::tolower);
591 if (problem_type_to_lower ==
"standard") {
593 }
else if (problem_type_to_lower ==
"manifold") {
602 std::string constraint_type_to_lower = constraint_type;
603 transform(constraint_type_to_lower.begin(), constraint_type_to_lower.end(), constraint_type_to_lower.begin(), ::tolower);
604 if (constraint_type_to_lower ==
"none") {
606 }
else if (constraint_type_to_lower ==
"neumann_grad_scalar") {
626 const int poly_order,
627 const int dimensions = 3,
628 const std::string dense_solver_type = std::string(
"QR"),
629 const std::string problem_type = std::string(
"STANDARD"),
630 const std::string constraint_type = std::string(
"NO_CONSTRAINT"),
631 const int manifold_curvature_poly_order = 2) :
654 (
"operations needed for manifold gradient reconstruction", 1);
655 auto curvature_support_operations_mirror =
657 curvature_support_operations_mirror(0) =
663 _lro = std::vector<TargetOperation>();
704 const int dimensions = 3,
705 const std::string dense_solver_type = std::string(
"QR"),
706 const std::string problem_type = std::string(
"STANDARD"),
707 const std::string constraint_type = std::string(
"NO_CONSTRAINT"),
708 const int manifold_curvature_poly_order = 2)
715 const int poly_order,
716 const int dimensions = 3,
717 const std::string dense_solver_type = std::string(
"QR"),
718 const std::string problem_type = std::string(
"STANDARD"),
719 const std::string constraint_type = std::string(
"NO_CONSTRAINT"),
720 const int manifold_curvature_poly_order = 2)
721 :
GMLS(reconstruction_space, dual_sampling_strategy, dual_sampling_strategy, poly_order, dimensions, dense_solver_type, problem_type, constraint_type, manifold_curvature_poly_order) {}
768 KOKKOS_INLINE_FUNCTION
772 KOKKOS_INLINE_FUNCTION
776 KOKKOS_INLINE_FUNCTION
780 KOKKOS_INLINE_FUNCTION
784 KOKKOS_INLINE_FUNCTION
789 KOKKOS_INLINE_FUNCTION
793 KOKKOS_INLINE_FUNCTION
797 KOKKOS_INLINE_FUNCTION
801 KOKKOS_INLINE_FUNCTION
805 KOKKOS_INLINE_FUNCTION
818 KOKKOS_INLINE_FUNCTION
828 KOKKOS_INLINE_FUNCTION
831 const int np =
getNP(m, dimension, r_space);
835 nn = np * (1.7 + m*0.1);
838 nn = np * (1.4 + m*0.03);
852 KOKKOS_INLINE_FUNCTION
855 return std::pow(1.0-std::abs(r/h), power) * double(1.0-std::abs(r/h)>0.0);
858 double h_over_3 = h/3.0;
859 return 1./( h_over_3 * 2.5066282746310002416124 ) * std::exp(-.5*r*r/(h_over_3*h_over_3));
862 return ((1-x)+x*(1-x)*(1-2*x)) *
double(r<=h);
870 KOKKOS_INLINE_FUNCTION
872 double inside_val = delta_vector.
x*delta_vector.
x;
875 inside_val += delta_vector.
z*delta_vector.
z;
878 inside_val += delta_vector.
y*delta_vector.
y;
883 return std::sqrt(inside_val);
888 static int getTargetOutputIndex(
const int operation_num,
const int output_component_axis_1,
const int output_component_axis_2,
const int dimensions) {
890 return axis_1_size*output_component_axis_1 + output_component_axis_2;
896 return axis_1_size*output_component_axis_1 + output_component_axis_2;
936 &&
"Entire batch not computed at once, so getFullPolynomialCoefficientsBasis() can not be called.");
938 &&
"generateAlphas() called with keep_coefficients set to false.");
1000 double getTangentBundle(
const int target_index,
const int direction,
const int component)
const {
1003 scratch_matrix_right_type::HostMirror
1005 return T(direction, component);
1011 "getRefenceNormalDirection called, but reference outwrad normal directions were never provided.");
1012 scratch_vector_type::HostMirror
1014 return ref_N(component);
1046 const int output_component_axis_2,
const int input_component_axis_1,
1047 const int input_component_axis_2,
const int additional_evaluation_local_index = 0)
const {
1050 compadre_assert_debug((lro_number >= 0) &&
"getAlphaColumnOffset called for a TargetOperation that was not registered.");
1068 &&
"Entire batch not computed at once, so getFullPolynomialCoefficientsBasis() can not be called.");
1070 &&
"generateAlphas() called with keep_coefficients set to false.");
1090 return getAlpha(lro, target_index, 0, 0, neighbor_index, 0, 0, additional_evaluation_site);
1096 return getAlpha(lro, target_index, output_component, 0, neighbor_index, 0, 0, additional_evaluation_site);
1100 double getAlpha0TensorTo2Tensor(
TargetOperation lro,
const int target_index,
const int output_component_axis_1,
const int output_component_axis_2,
const int neighbor_index,
const int additional_evaluation_site = 0)
const {
1101 return getAlpha(lro, target_index, output_component_axis_1, output_component_axis_2, neighbor_index, 0, 0, additional_evaluation_site);
1107 return getAlpha(lro, target_index, 0, 0, neighbor_index, input_component, 0, additional_evaluation_site);
1113 return getAlpha(lro, target_index, output_component, 0, neighbor_index, input_component, 0, additional_evaluation_site);
1117 double getAlpha1TensorTo2Tensor(
TargetOperation lro,
const int target_index,
const int output_component_axis_1,
const int output_component_axis_2,
const int neighbor_index,
const int input_component,
const int additional_evaluation_site = 0)
const {
1119 return getAlpha(lro, target_index, output_component_axis_1, output_component_axis_2, neighbor_index, input_component, 0, additional_evaluation_site);
1123 double getAlpha2TensorTo0Tensor(
TargetOperation lro,
const int target_index,
const int neighbor_index,
const int input_component_axis_1,
const int input_component_axis_2,
const int additional_evaluation_site = 0)
const {
1124 return getAlpha(lro, target_index, 0, 0, neighbor_index, input_component_axis_1, input_component_axis_2, additional_evaluation_site);
1128 double getAlpha2TensorTo1Tensor(
TargetOperation lro,
const int target_index,
const int output_component,
const int neighbor_index,
const int input_component_axis_1,
const int input_component_axis_2,
const int additional_evaluation_site = 0)
const {
1129 return getAlpha(lro, target_index, output_component, 0, neighbor_index, input_component_axis_1, input_component_axis_2, additional_evaluation_site);
1133 double getAlpha2TensorTo2Tensor(
TargetOperation lro,
const int target_index,
const int output_component_axis_1,
const int output_component_axis_2,
const int neighbor_index,
const int input_component_axis_1,
const int input_component_axis_2,
const int additional_evaluation_site = 0)
const {
1134 return getAlpha(lro, target_index, output_component_axis_1, output_component_axis_2, neighbor_index, input_component_axis_1, input_component_axis_2, additional_evaluation_site);
1139 KOKKOS_INLINE_FUNCTION
1147 return (total_neighbors_before_target+
TO_GLOBAL(total_added_alphas_before_target))
1149 +
TO_GLOBAL(alpha_column_offset*alphas_per_tile_per_target);
1162 return (total_neighbors_before_target+
TO_GLOBAL(total_added_alphas_before_target))
1164 +
TO_GLOBAL(alpha_column_offset*alphas_per_tile_per_target);
1169 double getAlpha(
TargetOperation lro,
const int target_index,
const int output_component_axis_1,
const int output_component_axis_2,
const int neighbor_index,
const int input_component_axis_1,
const int input_component_axis_2,
const int additional_evaluation_site = 0)
const {
1191 output_component_axis_2, input_component_axis_1, input_component_axis_2, additional_evaluation_site);
1203 if (for_target)
return 0;
else return 1;
1209 const int target_index_in_weights =
1213 const int neighbor_index_in_weights =
1218 output_component, input_component);
1231 return this->_host_lro_input_tile_size[
_lro_lookup[(int)lro]];
1264 sm = std::min(bm,sm);
1279 if (
_RHS.extent(0) > 0)
1280 _RHS = Kokkos::View<double*>(
"RHS",0);
1284 template<
typename view_type_1,
typename view_type_2,
typename view_type_3,
typename view_type_4>
1286 view_type_1 neighbor_lists,
1287 view_type_2 source_coordinates,
1288 view_type_3 target_coordinates,
1289 view_type_4 epsilons) {
1290 this->setNeighborLists<view_type_1>(neighbor_lists);
1291 this->setSourceSites<view_type_2>(source_coordinates);
1292 this->setTargetSites<view_type_3>(target_coordinates);
1293 this->setWindowSizes<view_type_4>(epsilons);
1297 template<
typename view_type_1,
typename view_type_2,
typename view_type_3,
typename view_type_4>
1299 view_type_1 cr_neighbor_lists,
1300 view_type_1 number_of_neighbors_list,
1301 view_type_2 source_coordinates,
1302 view_type_3 target_coordinates,
1303 view_type_4 epsilons) {
1304 this->setNeighborLists<view_type_1>(cr_neighbor_lists, number_of_neighbors_list);
1305 this->setSourceSites<view_type_2>(source_coordinates);
1306 this->setTargetSites<view_type_3>(target_coordinates);
1307 this->setWindowSizes<view_type_4>(epsilons);
1311 template<
typename view_type_1,
typename view_type_2>
1313 view_type_1 additional_evaluation_indices,
1314 view_type_2 additional_evaluation_coordinates) {
1315 this->setAuxiliaryEvaluationIndicesLists<view_type_1>(additional_evaluation_indices);
1316 this->setAuxiliaryEvaluationCoordinates<view_type_2>(additional_evaluation_coordinates);
1320 template <
typename view_type>
1321 typename std::enable_if<view_type::rank==1&&std::is_same<decltype(_neighbor_lists)::internal_view_type,view_type>::value==1,
void>::type
1327 Kokkos::parallel_for(
"copy neighbor list sizes", Kokkos::RangePolicy<host_execution_space>(0,
_host_number_of_neighbors_list.extent(0)), KOKKOS_LAMBDA(
const int i) {
1336 template <
typename view_type>
1337 typename std::enable_if<view_type::rank==1&&std::is_same<decltype(_neighbor_lists)::internal_view_type,view_type>::value==0,
void>::type
1341 gmls_view_type d_neighbor_lists(
"compressed row neighbor lists data", neighbor_lists.extent(0));
1342 gmls_view_type d_number_of_neighbors_list(
"number of neighbors list", number_of_neighbors_list.extent(0));
1343 Kokkos::deep_copy(d_neighbor_lists, neighbor_lists);
1344 Kokkos::deep_copy(d_number_of_neighbors_list, number_of_neighbors_list);
1349 Kokkos::parallel_for(
"copy neighbor list sizes", Kokkos::RangePolicy<host_execution_space>(0,
_host_number_of_neighbors_list.extent(0)), KOKKOS_LAMBDA(
const int i) {
1359 template <
typename view_type>
1360 typename std::enable_if<view_type::rank==2, void>::type
setNeighborLists(view_type neighbor_lists) {
1362 _neighbor_lists = Convert2DToCompressedRowNeighborLists<decltype(neighbor_lists), Kokkos::View<int*> >(neighbor_lists);
1365 Kokkos::parallel_for(
"copy neighbor list sizes", Kokkos::RangePolicy<host_execution_space>(0,
_host_number_of_neighbors_list.extent(0)), KOKKOS_LAMBDA(
const int i) {
1375 template<
typename view_type>
1380 source_coordinates.extent(0), source_coordinates.extent(1));
1382 typedef typename view_type::memory_space input_array_memory_space;
1383 if (std::is_same<input_array_memory_space, device_memory_space>::value) {
1394 Kokkos::deep_copy(host_source_coordinates, source_coordinates);
1403 template<
typename view_type>
1411 template<
typename view_type>
1415 target_coordinates.extent(0), target_coordinates.extent(1));
1417 typedef typename view_type::memory_space input_array_memory_space;
1418 if (std::is_same<input_array_memory_space, device_memory_space>::value) {
1429 Kokkos::deep_copy(host_target_coordinates, target_coordinates);
1440 template<
typename view_type>
1451 template<
typename view_type>
1465 template<
typename view_type>
1477 template<
typename view_type>
1487 "Memory space does not match between _T and tangent_directions");
1492 Kokkos::parallel_for(
"copy tangent vectors", Kokkos::RangePolicy<device_execution_space>(0,
_target_coordinates.extent(0)), KOKKOS_LAMBDA(
const int i) {
1494 for (
int j=0; j<this_dimensions; ++j) {
1495 for (
int k=0; k<this_dimensions; ++k) {
1496 T(j,k) = tangent_directions(i, j, k);
1503 _host_T = Kokkos::create_mirror_view(
_T);
1511 template<
typename view_type>
1518 auto this_ref_N = this->
_ref_N;
1522 Kokkos::parallel_for(
"copy normal vectors", Kokkos::RangePolicy<device_execution_space>(0,
_target_coordinates.extent(0)), KOKKOS_LAMBDA(
const int i) {
1523 for (
int j=0; j<this_dimensions; ++j) {
1524 this_ref_N(i*this_dimensions + j) = outward_normal_directions(i, j);
1539 template<
typename view_type>
1546 Kokkos::deep_copy(host_extra_data, extra_data);
1554 template<
typename view_type>
1563 template<
typename view_type>
1570 Kokkos::deep_copy(host_extra_data, extra_data);
1578 template<
typename view_type>
1589 template <
typename view_type>
1593 evaluation_coordinates.extent(0), evaluation_coordinates.extent(1));
1595 typedef typename view_type::memory_space input_array_memory_space;
1596 if (std::is_same<input_array_memory_space, device_memory_space>::value) {
1607 Kokkos::deep_copy(host_additional_evaluation_coordinates, evaluation_coordinates);
1618 template <
typename view_type>
1627 template <
typename view_type>
1631 indices_lists.extent(0), indices_lists.extent(1));
1635 typedef typename view_type::memory_space input_array_memory_space;
1636 if (std::is_same<input_array_memory_space, device_memory_space>::value) {
1655 Kokkos::parallel_reduce(
"additional evaluation indices",
1657 KOKKOS_LAMBDA(
const int i,
int& t_max_evaluation_sites_per_target) {
1658 number_of_additional_evaluation_indices(i) = additional_evaluation_indices(i,0);
1659 t_max_evaluation_sites_per_target = (t_max_evaluation_sites_per_target > number_of_additional_evaluation_indices(i)+1)
1660 ? t_max_evaluation_sites_per_target : number_of_additional_evaluation_indices(i)+1;
1669 template <
typename view_type>
1681 Kokkos::parallel_reduce(
"additional evaluation indices",
1683 KOKKOS_LAMBDA(
const int i,
int& t_max_evaluation_sites_per_target) {
1684 number_of_additional_evaluation_indices(i) = additional_evaluation_indices(i,0);
1685 t_max_evaluation_sites_per_target = (t_max_evaluation_sites_per_target > number_of_additional_evaluation_indices(i)+1)
1686 ? t_max_evaluation_sites_per_target : number_of_additional_evaluation_indices(i)+1;
1695 std::string wt_to_lower = wt;
1696 transform(wt_to_lower.begin(), wt_to_lower.end(), wt_to_lower.begin(), ::tolower);
1697 if (wt_to_lower ==
"power") {
1699 }
else if (wt_to_lower ==
"gaussian") {
1701 }
else if (wt_to_lower ==
"cubicspline") {
1718 std::string wt_to_lower = wt;
1719 transform(wt_to_lower.begin(), wt_to_lower.end(), wt_to_lower.begin(), ::tolower);
1720 if (wt_to_lower ==
"power") {
1722 }
else if (wt_to_lower ==
"gaussian") {
1724 }
else if (wt_to_lower ==
"cubicspline") {
1784 std::vector<TargetOperation> temporary_lro_vector(1, lro);
1795 for (
size_t i=0; i<lro.size(); ++i) {
1797 bool operation_found =
false;
1799 for (
size_t j=0; j<
_lro.size(); ++j) {
1802 if (
_lro[j]==lro[i]) {
1804 operation_found =
true;
1814 if (!operation_found) {
1816 _lro.push_back(lro[i]);
1832 int total_offset = 0;
1833 int output_offset = 0;
1834 int input_offset = 0;
1835 for (
size_t i=0; i<
_lro.size(); ++i) {
1847 total_offset += input_tile_size * output_tile_size;
1848 output_offset += output_tile_size;
1849 input_offset += input_tile_size;
1898 void generateAlphas(
const int number_of_batches = 1,
const bool keep_coefficients =
false);
Kokkos::View< double * > _manifold_metric_tensor_inverse
metric tensor inverse for all problems
KOKKOS_INLINE_FUNCTION double convertLocalToGlobalCoordinate(const XYZ local_coord, const int dim, const scratch_matrix_right_type *V) const
Returns a component of the global coordinate after transformation from local to global under the orth...
Divergence-free vector polynomial basis.
Kokkos::View< int * > _lro_input_tile_size
dimensions ^ rank of tensor of output for each sampling functional (device)
Kokkos::View< double *, layout_right > _alphas
generated alpha coefficients (device)
KOKKOS_INLINE_FUNCTION int getNumberOfNeighborsDevice(int target_index) const
Get number of neighbors for a given target (device)
decltype(_RHS) getFullPolynomialCoefficientsBasis() const
Get a view (device) of all polynomial coefficients basis.
Standard GMLS problem type.
std::size_t global_index_type
Kokkos::View< const double *****, layout_right >::HostMirror _host_prestencil_weights
generated weights for nontraditional samples required to transform data into expected sampling functi...
Kokkos::View< TargetOperation * > _operations
vector containing target functionals to be applied for reconstruction problem (device) ...
bool _entire_batch_computed_at_once
whether entire calculation was computed at once the alternative is that it was broken up over many sm...
KOKKOS_INLINE_FUNCTION double getNeighborCoordinate(const int target_index, const int neighbor_list_num, const int dim, const scratch_matrix_right_type *V=NULL) const
Returns one component of the neighbor coordinate for a particular target.
int _max_evaluation_sites_per_target
maximum number of evaluation sites for each target (includes target site)
ProblemType getProblemType()
Get problem type.
Kokkos::View< double * > _manifold_curvature_coefficients
curvature polynomial coefficients for all problems
int _curvature_weighting_power
power to be used for weighting kernel for curvature
KOKKOS_INLINE_FUNCTION void createWeightsAndP(const member_type &teamMember, scratch_vector_type delta, scratch_vector_type thread_workspace, scratch_matrix_right_type P, scratch_vector_type w, const int dimension, int polynomial_order, bool weight_p=false, scratch_matrix_right_type *V=NULL, const ReconstructionSpace reconstruction_space=ReconstructionSpace::ScalarTaylorPolynomial, const SamplingFunctional sampling_strategy=PointSample) const
Fills the _P matrix with either P or P*sqrt(w)
int output_rank
Rank of sampling functional output for each SamplingFunctional.
void generatePolynomialCoefficients(const int number_of_batches=1, const bool keep_coefficients=false)
Generates polynomial coefficients by setting up and solving least squares problems ! Sets up the batc...
static KOKKOS_INLINE_FUNCTION double EuclideanVectorLength(const XYZ &delta_vector, const int dimension)
Returns Euclidean norm of a vector.
Kokkos::View< int **, layout_right > _additional_evaluation_indices
(OPTIONAL) contains indices of entries in the _additional_evaluation_coordinates view (device) ...
Tag for functor to evaluate curvature targets and apply to coefficients of curvature reconstruction...
Kokkos::View< double * > _ref_N
Rank 2 tensor for high order approximation of tangent vectors for all problems.
KOKKOS_INLINE_FUNCTION void calcGradientPij(const member_type &teamMember, double *delta, double *thread_workspace, const int target_index, int neighbor_index, const double alpha, const int partial_direction, const int dimension, const int poly_order, bool specific_order_only, const scratch_matrix_right_type *V, const ReconstructionSpace reconstruction_space, const SamplingFunctional sampling_strategy, const int additional_evaluation_local_index=0) const
Evaluates the gradient of a polynomial basis under the Dirac Delta (pointwise) sampling function...
KOKKOS_INLINE_FUNCTION int getSize(const int degree, const int dimension)
Returns size of basis.
#define compadre_kernel_assert_release(condition)
compadre_kernel_assert_release is similar to compadre_assert_release, but is a call on the device...
constexpr SamplingFunctional ManifoldVectorPointSample
Point evaluations of the entire vector source function (but on a manifold, so it includes a transform...
GMLS(const int poly_order, const int dimensions=3, const std::string dense_solver_type=std::string("QR"), const std::string problem_type=std::string("STANDARD"), const std::string constraint_type=std::string("NO_CONSTRAINT"), const int manifold_curvature_poly_order=2)
Constructor for the case when the data sampling functional does not match the polynomial sampling fun...
bool _reference_outward_normal_direction_provided
whether or not the reference outward normal directions were provided by the user. ...
KOKKOS_INLINE_FUNCTION int getMaxNumNeighbors() const
Get the maximum number of neighbors of all targets' neighborhoods (host/device)
team_policy::member_type member_type
int getOrderOfQuadraturePoints() const
Order of quadrature points.
Neumann Gradient Scalar Type.
double getAlpha2TensorTo0Tensor(TargetOperation lro, const int target_index, const int neighbor_index, const int input_component_axis_1, const int input_component_axis_2, const int additional_evaluation_site=0) const
Helper function for getting alphas for scalar reconstruction from matrix data.
int _NP
dimension of basis for polynomial reconstruction
int _total_alpha_values
used for sizing P_target_row and the _alphas view
Kokkos::View< double * > _epsilons
h supports determined through neighbor search (device)
bool _orthonormal_tangent_space_provided
whether or not the orthonormal tangent directions were provided by the user.
void setProblemData(view_type_1 cr_neighbor_lists, view_type_1 number_of_neighbors_list, view_type_2 source_coordinates, view_type_3 target_coordinates, view_type_4 epsilons)
Sets basic problem data (neighbor lists data, number of neighbors list, source coordinates, and target coordinates)
Kokkos::View< int *, host_memory_space > _host_number_of_neighbors_list
convenient copy on host of number of neighbors
void setQuadratureType(std::string quadrature_type)
Type of quadrature points.
static int getTargetOutputIndex(const int operation_num, const int output_component_axis_1, const int output_component_axis_2, const int dimensions)
Helper function for finding alpha coefficients.
bool _use_reference_outward_normal_direction_provided_to_orient_surface
whether or not to use reference outward normal directions to orient the surface in a manifold problem...
Kokkos::View< double * > _RHS
sqrt(w)*Identity matrix for all problems, later holds polynomial coefficients for all problems ...
Kokkos::View< double * > _w
contains weights for all problems
void setWindowSizes(view_type epsilons)
Sets window sizes, also called the support of the kernel.
int _sampling_multiplier
actual dimension of the sampling functional e.g.
int _reconstruction_space_rank
actual rank of reconstruction basis
int _curvature_poly_order
order of basis for curvature reconstruction
int getWeightingPower() const
Power for weighting kernel for GMLS problem.
Kokkos::View< double **, layout_right > _target_coordinates
coordinates for target sites for reconstruction (device)
Each target applies a different data transform, but the same to each neighbor.
KOKKOS_INLINE_FUNCTION int getAdditionalEvaluationIndex(const int target_index, const int additional_list_num) const
(OPTIONAL) Mapping from [0,number of additional evaluation sites for a target] to the row that contai...
std::enable_if< view_type::rank==2, void >::type setNeighborLists(view_type neighbor_lists)
Sets neighbor list information.
int _data_sampling_multiplier
effective dimension of the data sampling functional e.g.
std::vector< int > _lro_lookup
vector containing a mapping from a target functionals enum value to the its place in the list of targ...
int getInputDimensionOfSampling(SamplingFunctional sro) const
Dimensions ^ output rank for sampling operation (always in ambient space, never local chart on a mani...
int getNumberOfQuadraturePoints() const
Number of quadrature points.
Each target applies a different transform for each neighbor.
NeighborLists< Kokkos::View< int * > > _neighbor_lists
Accessor to get neighbor list data, offset data, and number of neighbors per target.
double getPreStencilWeight(SamplingFunctional sro, const int target_index, const int neighbor_index, bool for_target, const int output_component=0, const int input_component=0) const
Returns a stencil to transform data from its existing state into the input expected for some sampling...
Quadrature _qm
manages and calculates quadrature
Kokkos::View< double * >::HostMirror _host_T
tangent vectors information (host)
decltype(_ref_N) getReferenceNormalDirections() const
Get a view (device) of all reference outward normal directions.
Kokkos::View< int * >::HostMirror _host_lro_input_tensor_rank
tensor rank of sampling functional (host)
double getAlpha0TensorTo0Tensor(TargetOperation lro, const int target_index, const int neighbor_index, const int additional_evaluation_site=0) const
Helper function for getting alphas for scalar reconstruction from scalar data.
KOKKOS_INLINE_FUNCTION double getTargetCoordinate(const int target_index, const int dim, const scratch_matrix_right_type *V=NULL) const
Returns one component of the target coordinate for a particular target.
void clearTargets()
Empties the vector of target functionals to apply to the reconstruction.
std::vector< TargetOperation > _lro
vector of user requested target operations
void setTargetExtraData(view_type extra_data)
(OPTIONAL) Sets extra data to be used by target operations in certain instances.
int getTargetOperationLocalIndex(TargetOperation lro) const
Get the local index (internal) to GMLS for a particular TargetOperation Every TargetOperation has a g...
void setReferenceOutwardNormalDirection(view_type outward_normal_directions, bool use_to_orient_surface=true)
(OPTIONAL) Sets outward normal direction.
void setAuxiliaryEvaluationIndicesLists(decltype(_additional_evaluation_indices) indices_lists)
(OPTIONAL) Sets the additional target evaluation coordinate indices list information.
void setDimensionOfQuadraturePoints(int dim)
Dimensions of quadrature points to use.
Scalar polynomial basis centered at the target site and scaled by sum of basis powers e...
double getAlpha1TensorTo0Tensor(TargetOperation lro, const int target_index, const int neighbor_index, const int input_component, const int additional_evaluation_site=0) const
Helper function for getting alphas for scalar reconstruction from vector data.
Solve GMLS problem on a manifold (will use QR or SVD to solve the resultant GMLS problem dependent on...
KOKKOS_INLINE_FUNCTION void calcPij(const member_type &teamMember, double *delta, double *thread_workspace, const int target_index, int neighbor_index, const double alpha, const int dimension, const int poly_order, bool specific_order_only=false, const scratch_matrix_right_type *V=NULL, const ReconstructionSpace reconstruction_space=ReconstructionSpace::ScalarTaylorPolynomial, const SamplingFunctional sampling_strategy=PointSample, const int additional_evaluation_local_index=0) const
Evaluates the polynomial basis under a particular sampling function. Generally used to fill a row of ...
int _max_num_neighbors
maximum number of neighbors over all target sites
Kokkos::View< int * > _lro_input_tensor_rank
tensor rank of sampling functional (device)
Should be the total count of all available target functionals.
KOKKOS_INLINE_FUNCTION double convertGlobalToLocalCoordinate(const XYZ global_coord, const int dim, const scratch_matrix_right_type *V) const
Returns a component of the local coordinate after transformation from global to local under the ortho...
SamplingFunctional getPolynomialSamplingFunctional() const
Get the polynomial sampling functional specified at instantiation.
int getLocalDimensions() const
Local dimension of the GMLS problem (less than global dimension if on a manifold), set only at class instantiation.
Kokkos::View< int * >::HostMirror _host_lro_output_tensor_rank
tensor rank of target functional (host)
WeightingFunctionType _curvature_weighting_type
weighting kernel type for curvature problem
WeightingFunctionType _weighting_type
weighting kernel type for GMLS
Kokkos::View< int * > _lro_output_tile_size
dimensions ^ rank of tensor of output for each target functional (device)
static KOKKOS_INLINE_FUNCTION int getNN(const int m, const int dimension=3, const ReconstructionSpace r_space=ReconstructionSpace::ScalarTaylorPolynomial)
Returns number of neighbors needed for unisolvency for a given basis order and dimension.
const SamplingFunctional _data_sampling_functional
generally the same as _polynomial_sampling_functional, but can differ if specified at GMLS class inst...
constexpr int TargetOutputTensorRank[]
Rank of target functional output for each TargetOperation Rank of target functional input for each Ta...
Tag for functor to evaluate targets, apply target evaluation to polynomial coefficients to store in _...
double getAlpha0TensorTo2Tensor(TargetOperation lro, const int target_index, const int output_component_axis_1, const int output_component_axis_2, const int neighbor_index, const int additional_evaluation_site=0) const
Helper function for getting alphas for matrix reconstruction from scalar data.
static ConstraintType parseConstraintType(const std::string &constraint_type)
Parses a string to determine constraint type.
KOKKOS_INLINE_FUNCTION int getNNeighbors(const int target_index) const
Returns number of neighbors for a particular target.
int calculateBasisMultiplier(const ReconstructionSpace rs) const
Calculate basis_multiplier.
DenseSolverType _dense_solver_type
solver type for GMLS problem - can be QR, SVD or LU
static ProblemType parseProblemType(const std::string &problem_type)
Parses a string to determine problem type.
KOKKOS_INLINE_FUNCTION int getNEvaluationSitesPerTarget(const int target_index) const
(OPTIONAL) Returns number of additional evaluation sites for a particular target
QR+Pivoting factorization performed on P*sqrt(w) matrix.
void addTargets(TargetOperation lro)
Adds a target to the vector of target functional to be applied to the reconstruction.
Tag for functor to create a coarse tangent approximation from a given neighborhood of points...
decltype(_alphas) getAlphas() const
Get a view (device) of all alphas.
pool_type _random_number_pool
Kokkos::View< int * > _lro_total_offsets
index for where this operation begins the for _alpha coefficients (device)
Kokkos::View< double * > _manifold_curvature_gradient
_dimension-1 gradient values for curvature for all problems
void setSourceExtraData(view_type extra_data)
(OPTIONAL) Sets extra data to be used by sampling functionals in certain instances.
ProblemType
Problem type, that optionally can handle manifolds.
ConstraintType getConstraintType()
Get constraint type.
Kokkos::View< int *, host_execution_space > host_managed_local_index_type
void setAuxiliaryEvaluationCoordinates(decltype(_additional_evaluation_coordinates) evaluation_coordinates)
(OPTIONAL) Sets additional points for evaluation of target operation on polynomial reconstruction...
void setAuxiliaryEvaluationIndicesLists(view_type indices_lists)
(OPTIONAL) Sets the additional target evaluation coordinate indices list information.
Kokkos::View< double *****, layout_right > _prestencil_weights
generated weights for nontraditional samples required to transform data into expected sampling functi...
void addTargets(std::vector< TargetOperation > lro)
Adds a vector of target functionals to the vector of target functionals already to be applied to the ...
#define TO_GLOBAL(variable)
KOKKOS_INLINE_FUNCTION int getSize(const int degree, const int dimension)
Returns size of basis.
void setCurvatureWeightingType(const std::string &wt)
Type for weighting kernel for curvature.
std::enable_if< view_type::rank==1 &&std::is_same< decltype(_neighbor_lists)::internal_view_type, view_type >::value==1, void >::type setNeighborLists(view_type neighbor_lists, view_type number_of_neighbors_list)
Sets neighbor list information from compressed row neighborhood lists data (if same view_type)...
int _initial_index_for_batch
initial index for current batch
int input_rank
Rank of sampling functional input for each SamplingFunctional.
KOKKOS_INLINE_FUNCTION int getTargetOffsetIndexDevice(const int lro_num, const int input_component, const int output_component, const int additional_evaluation_local_index=0) const
Handles offset from operation input/output + extra evaluation sites.
static KOKKOS_INLINE_FUNCTION int getNP(const int m, const int dimension=3, const ReconstructionSpace r_space=ReconstructionSpace::ScalarTaylorPolynomial)
Returns size of the basis for a given polynomial order and dimension General to dimension 1...
SamplingFunctional getDataSamplingFunctional() const
Get the data sampling functional specified at instantiation (often the same as the polynomial samplin...
void setSourceSites(decltype(_source_coordinates) source_coordinates)
Sets source coordinate information.
Kokkos::View< double **, layout_right > _source_extra_data
Extra data available to basis functions (optional)
KOKKOS_INLINE_FUNCTION int getMaxEvaluationSitesPerTarget() const
Returns the maximum number of evaluation sites over all target sites (target sites are included in to...
void setTargetSites(decltype(_target_coordinates) target_coordinates)
Sets target coordinate information. Rows of this 2D-array should correspond to rows of the neighbor l...
Scalar basis reused as many times as there are components in the vector resulting in a much cheaper p...
Kokkos::View< int **, layout_right >::HostMirror _host_additional_evaluation_indices
(OPTIONAL) contains indices of entries in the _additional_evaluation_coordinates view (host) ...
KOKKOS_INLINE_FUNCTION int getMaxNNeighbors() const
Returns the maximum neighbor lists size over all target sites.
void setAdditionalEvaluationSitesData(view_type_1 additional_evaluation_indices, view_type_2 additional_evaluation_coordinates)
(OPTIONAL) Sets additional evaluation sites for each target site
std::string getQuadratureType() const
Type of quadrature points.
ReconstructionSpace
Space in which to reconstruct polynomial.
const SamplingFunctional _polynomial_sampling_functional
polynomial sampling functional used to construct P matrix, set at GMLS class instantiation ...
Tag for functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity.
std::enable_if< view_type::rank==1 &&std::is_same< decltype(_neighbor_lists)::internal_view_type, view_type >::value==0, void >::type setNeighborLists(view_type neighbor_lists, view_type number_of_neighbors_list)
Sets neighbor list information from compressed row neighborhood lists data (if different view_type)...
DenseSolverType getDenseSolverType()
Get dense solver type.
GMLS(ReconstructionSpace reconstruction_space, SamplingFunctional dual_sampling_strategy, const int poly_order, const int dimensions=3, const std::string dense_solver_type=std::string("QR"), const std::string problem_type=std::string("STANDARD"), const std::string constraint_type=std::string("NO_CONSTRAINT"), const int manifold_curvature_poly_order=2)
Constructor for the case when nonstandard sampling functionals or reconstruction spaces are to be use...
decltype(_T) getTangentDirections() const
Get a view (device) of all tangent direction bundles.
Tag for functor to calculate prestencil weights to apply to data to transform into a format expected ...
KOKKOS_INLINE_FUNCTION void getPDims(DenseSolverType dense_solver_type, ConstraintType constraint_type, ReconstructionSpace reconstruction_space, const int dimension, const int M, const int N, int &out_row, int &out_col)
Kokkos::View< double * >::HostMirror _host_epsilons
h supports determined through neighbor search (host)
host_managed_local_index_type getPolynomialCoefficientsMemorySize() const
Returns 2D array size in memory on which coefficients are stored.
double getAlpha0TensorTo1Tensor(TargetOperation lro, const int target_index, const int output_component, const int neighbor_index, const int additional_evaluation_site=0) const
Helper function for getting alphas for vector reconstruction from scalar data.
Kokkos::View< double **, layout_right, Kokkos::MemoryTraits< Kokkos::Unmanaged > > scratch_matrix_right_type
Kokkos::View< double *, Kokkos::MemoryTraits< Kokkos::Unmanaged > > scratch_vector_type
double getAlpha2TensorTo1Tensor(TargetOperation lro, const int target_index, const int output_component, const int neighbor_index, const int input_component_axis_1, const int input_component_axis_2, const int additional_evaluation_site=0) const
Helper function for getting alphas for vector reconstruction from matrix data.
KOKKOS_INLINE_FUNCTION int getNumberOfQuadraturePoints() const
int getAlphaColumnOffset(TargetOperation lro, const int output_component_axis_1, const int output_component_axis_2, const int input_component_axis_1, const int input_component_axis_2, const int additional_evaluation_local_index=0) const
Retrieves the offset for an operator based on input and output component, generic to row (but still m...
void setWeightingType(const WeightingFunctionType wt)
Type for weighting kernel for GMLS problem.
double getAlpha1TensorTo2Tensor(TargetOperation lro, const int target_index, const int output_component_axis_1, const int output_component_axis_2, const int neighbor_index, const int input_component, const int additional_evaluation_site=0) const
Helper function for getting alphas for matrix reconstruction from vector data.
WeightingFunctionType
Available weighting kernel function types.
int getOutputDimensionOfSampling(SamplingFunctional sro) const
Dimensions ^ output rank for sampling operation (always in local chart if on a manifold, never ambient space)
NeighborLists assists in accessing entries of compressed row neighborhood lists.
Kokkos::View< double **, layout_right > _source_coordinates
all coordinates for the source for which _neighbor_lists refers (device)
Kokkos::View< int * > _lro_output_tensor_rank
tensor rank of target functional (device)
void setCurvaturePolynomialOrder(const int manifold_poly_order)
Sets basis order to be used when reoncstructing curvature.
double getAlpha(TargetOperation lro, const int target_index, const int output_component_axis_1, const int output_component_axis_2, const int neighbor_index, const int input_component_axis_1, const int input_component_axis_2, const int additional_evaluation_site=0) const
Underlying function all interface helper functions call to retrieve alpha values. ...
void generateAlphas(const int number_of_batches=1, const bool keep_coefficients=false)
Meant to calculate target operations and apply the evaluations to the previously ! constructed polyno...
Kokkos::View< double **, layout_right > _target_extra_data
Extra data available to target operations (optional)
int _dimension_of_quadrature_points
dimension of quadrature rule
global_index_type getAlphaIndexHost(const int target_index, const int alpha_column_offset) const
Gives index into alphas given two axes, which when incremented by the neighbor number transforms acce...
double getReferenceNormalDirection(const int target_index, const int component) const
Get component of tangent or normal directions for manifold problems.
Kokkos::View< double * >::HostMirror _host_ref_N
reference outward normal vectors information (host)
KOKKOS_INLINE_FUNCTION double getTargetAuxiliaryCoordinate(const int target_index, const int additional_list_num, const int dim, const scratch_matrix_right_type *V=NULL) const
(OPTIONAL) Returns one component of the additional evaluation coordinates.
int getManifoldWeightingPower() const
Power for weighting kernel for curvature.
Point evaluation of the gradient of a scalar.
KOKKOS_INLINE_FUNCTION int getNeighborDevice(int target_index, int neighbor_num) const
Offers N(i,j) indexing where N(i,j) is the index of the jth neighbor of i (device) ...
GMLS(ReconstructionSpace reconstruction_space, const SamplingFunctional polynomial_sampling_strategy, const SamplingFunctional data_sampling_strategy, const int poly_order, const int dimensions=3, const std::string dense_solver_type=std::string("QR"), const std::string problem_type=std::string("STANDARD"), const std::string constraint_type=std::string("NO_CONSTRAINT"), const int manifold_curvature_poly_order=2)
Minimal constructor providing no data (neighbor lists, source sites, target sites) ...
int getOutputDimensionOfOperation(TargetOperation lro, bool ambient=false) const
Dimensions ^ output rank for target operation.
double getAlpha1TensorTo1Tensor(TargetOperation lro, const int target_index, const int output_component, const int neighbor_index, const int input_component, const int additional_evaluation_site=0) const
Helper function for getting alphas for vector reconstruction from vector data.
Kokkos::View< int * >::HostMirror _host_lro_input_tile_size
dimensions ^ rank of tensor of output for each sampling functional (host)
ParallelManager _pm
determines scratch level spaces and is used to call kernels
host_managed_local_index_type getPolynomialCoefficientsDomainRangeSize() const
Returns (size of the basis used in instance's polynomial reconstruction) x (data input dimension) ...
void setWindowSizes(decltype(_epsilons) epsilons)
Sets window sizes, also called the support of the kernel (device)
int getDimensions() const
Dimension of the GMLS problem, set only at class instantiation.
int getPolynomialCoefficientsSize() const
Returns size of the basis used in instance's polynomial reconstruction.
Kokkos::View< TargetOperation * > _curvature_support_operations
vector containing target functionals to be applied for curvature
static int getInputRankOfSampling(SamplingFunctional sro)
Input rank for sampling operation.
Generalized Moving Least Squares (GMLS)
KOKKOS_INLINE_FUNCTION XYZ getRelativeCoord(const int target_index, const int neighbor_list_num, const int dimension, const scratch_matrix_right_type *V=NULL) const
Returns the relative coordinate as a vector between the target site and the neighbor site...
KOKKOS_INLINE_FUNCTION global_index_type getRowOffsetDevice(int target_index) const
Get offset into compressed row neighbor lists (device)
void setTargetExtraData(decltype(_target_extra_data) extra_data)
(OPTIONAL) Sets extra data to be used by target operations in certain instances.
int _added_alpha_size
additional alpha coefficients due to constraints
ProblemType _problem_type
problem type for GMLS problem, can also be set to STANDARD for normal or MANIFOLD for manifold proble...
int getInputDimensionOfOperation(TargetOperation lro) const
Dimensions ^ input rank for target operation (always in local chart if on a manifold, never ambient space)
int _basis_multiplier
dimension of the reconstructed function e.g.
WeightingFunctionType getWeightingType() const
Type for weighting kernel for GMLS problem.
int _global_dimensions
spatial dimension of the points, set at class instantiation only
Kokkos::View< double **, layout_right > _additional_evaluation_coordinates
(OPTIONAL) user provided additional coordinates for target operation evaluation (device) ...
DenseSolverType
Dense solver type.
Kokkos::View< double * > _T
Rank 3 tensor for high order approximation of tangent vectors for all problems.
double getTangentBundle(const int target_index, const int direction, const int component) const
Get component of tangent or normal directions for manifold problems.
void resetCoefficientData()
KOKKOS_INLINE_FUNCTION void computeTargetFunctionals(const member_type &teamMember, scratch_vector_type t1, scratch_vector_type t2, scratch_matrix_right_type P_target_row) const
Evaluates a polynomial basis with a target functional applied to each member of the basis...
int getDimensionOfQuadraturePoints() const
Dimensions of quadrature points.
ReconstructionSpace getReconstructionSpace() const
Get the reconstruction space specified at instantiation.
void setSourceSites(view_type source_coordinates)
Sets source coordinate information.
LU factorization performed on P^T*W*P matrix.
bool _store_PTWP_inv_PTW
whether polynomial coefficients were requested to be stored (in a state not yet applied to data) ...
static int getOutputRankOfSampling(SamplingFunctional sro)
Output rank for sampling operation.
int transform_type
Describes the SamplingFunction relationship to targets, neighbors.
Kokkos::View< int * > _number_of_additional_evaluation_indices
(OPTIONAL) contains the # of additional coordinate indices for each target
void setPolynomialOrder(const int poly_order)
Sets basis order to be used when reoncstructing any function.
constexpr int ActualReconstructionSpaceRank[]
Number of actual components in the ReconstructionSpace.
void setTargetSites(view_type target_coordinates)
Sets target coordinate information. Rows of this 2D-array should correspond to rows of the neighbor l...
ReconstructionSpace _reconstruction_space
reconstruction space for GMLS problems, set at GMLS class instantiation
void setAuxiliaryEvaluationCoordinates(view_type evaluation_coordinates)
(OPTIONAL) Sets additional points for evaluation of target operation on polynomial reconstruction...
Kokkos::View< const double *, layout_right >::HostMirror _host_alphas
generated alpha coefficients (host)
int getTargetOffsetIndexHost(const int lro_num, const int input_component, const int output_component, const int additional_evaluation_local_index=0) const
Handles offset from operation input/output + extra evaluation sites.
Kokkos::Random_XorShift64_Pool pool_type
void setCurvatureWeightingPower(int wp)
Power for weighting kernel for curvature.
int _weighting_power
power to be used for weighting kernel
Kokkos::View< double * > _P
P*sqrt(w) matrix for all problems.
Tag for functor to evaluate targets, apply target evaluation to polynomial coefficients to store in _...
Tag for functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity.
void setWeightingType(const std::string &wt)
Type for weighting kernel for GMLS problem.
KOKKOS_INLINE_FUNCTION void getRHSDims(DenseSolverType dense_solver_type, ConstraintType constraint_type, ReconstructionSpace reconstruction_space, const int dimension, const int M, const int N, int &RHS_row, int &RHS_col)
void setProblemData(view_type_1 neighbor_lists, view_type_2 source_coordinates, view_type_3 target_coordinates, view_type_4 epsilons)
Sets basic problem data (neighbor lists, source coordinates, and target coordinates) ...
constexpr SamplingFunctional PointSample
Available sampling functionals.
int _dimensions
dimension of the problem, set at class instantiation only
KOKKOS_INLINE_FUNCTION global_index_type getAlphaIndexDevice(const int target_index, const int alpha_column_offset) const
Gives index into alphas given two axes, which when incremented by the neighbor number transforms acce...
decltype(_neighbor_lists) * getNeighborLists()
Get neighbor list accessor.
std::string _quadrature_type
quadrature rule type
Kokkos::View< TargetOperation * >::HostMirror _host_operations
vector containing target functionals to be applied for reconstruction problem (host) ...
static KOKKOS_INLINE_FUNCTION double Wab(const double r, const double h, const WeightingFunctionType &weighting_type, const int power)
Evaluates the weighting kernel.
bool _nontrivial_nullspace
whether or not operator to be inverted for GMLS problem has a nontrivial nullspace (requiring SVD) ...
double getAlpha2TensorTo2Tensor(TargetOperation lro, const int target_index, const int output_component_axis_1, const int output_component_axis_2, const int neighbor_index, const int input_component_axis_1, const int input_component_axis_2, const int additional_evaluation_site=0) const
Helper function for getting alphas for matrix reconstruction from matrix data.
KOKKOS_INLINE_FUNCTION void calcHessianPij(const member_type &teamMember, double *delta, double *thread_workspace, const int target_index, int neighbor_index, const double alpha, const int partial_direction_1, const int partial_direction_2, const int dimension, const int poly_order, bool specific_order_only, const scratch_matrix_right_type *V, const ReconstructionSpace reconstruction_space, const SamplingFunctional sampling_strategy, const int additional_evaluation_local_index=0) const
Evaluates the Hessian of a polynomial basis under the Dirac Delta (pointwise) sampling function...
KOKKOS_INLINE_FUNCTION void applyTargetsToCoefficients(const member_type &teamMember, scratch_vector_type t1, scratch_vector_type t2, scratch_matrix_right_type Q, scratch_vector_type w, scratch_matrix_right_type P_target_row, const int target_NP) const
Helper function for applying the evaluations from a target functional to the polynomial coefficients...
Tag for functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity for curva...
#define compadre_assert_debug(condition)
compadre_assert_debug is used for assertions that are checked in loops, as these significantly impact...
Kokkos::View< int * >::HostMirror _host_lro_total_offsets
index for where this operation begins the for _alpha coefficients (host)
KOKKOS_INLINE_FUNCTION void operator()(const AssembleStandardPsqrtW &, const member_type &teamMember) const
Functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity.
int getNumberOfNeighborsHost(int target_index) const
Get number of neighbors for a given target (host)
void setWeightingPower(int wp)
Power for weighting kernel for GMLS problem.
void setCurvatureWeightingType(const WeightingFunctionType wt)
Type for weighting kernel for curvature.
int getGlobalDimensions() const
Dimension of the GMLS problem's point data (spatial description of points in ambient space)...
int calculateSamplingMultiplier(const ReconstructionSpace rs, const SamplingFunctional sro) const
Calculate sampling_multiplier.
ConstraintType _constraint_type
constraint type for GMLS problem
KOKKOS_INLINE_FUNCTION int getNeighborIndex(const int target_index, const int neighbor_list_num) const
Mapping from [0,number of neighbors for a target] to the row that contains the source coordinates for...
KOKKOS_INLINE_FUNCTION void computeCurvatureFunctionals(const member_type &teamMember, scratch_vector_type t1, scratch_vector_type t2, scratch_matrix_right_type P_target_row, const scratch_matrix_right_type *V, const local_index_type local_neighbor_index=-1) const
Evaluates a polynomial basis for the curvature with a gradient target functional applied.
int _local_dimensions
dimension of the problem, set at class instantiation only. For manifolds, generally _global_dimension...
Kokkos::View< int * >::HostMirror _host_lro_output_tile_size
dimensions ^ rank of tensor of output for each target functional (host)
global_index_type getRowOffsetHost(int target_index) const
Get offset into compressed row neighbor lists (host)
int _order_of_quadrature_points
order of exact polynomial integration for quadrature rule
#define compadre_assert_release(condition)
compadre_assert_release is used for assertions that should always be checked, but generally are not e...
Tag for functor to determine if tangent directions need reordered, and to reorder them if needed...
Tag for functor to evaluate curvature targets and construct accurate tangent direction approximation ...
decltype(_prestencil_weights) getPrestencilWeights() const
Get a view (device) of all rank 2 preprocessing tensors This is a rank 5 tensor that is able to provi...
KOKKOS_INLINE_FUNCTION void computeTargetFunctionalsOnManifold(const member_type &teamMember, scratch_vector_type t1, scratch_vector_type t2, scratch_matrix_right_type P_target_row, scratch_matrix_right_type V, scratch_matrix_right_type G_inv, scratch_vector_type curvature_coefficients, scratch_vector_type curvature_gradients) const
Evaluates a polynomial basis with a target functional applied, using information from the manifold cu...
WeightingFunctionType getManifoldWeightingType() const
Type for weighting kernel for curvature.
constexpr SamplingFunctional VectorPointSample
Point evaluations of the entire vector source function.
void setOrderOfQuadraturePoints(int order)
Number quadrature points to use.
static DenseSolverType parseSolverType(const std::string &dense_solver_type)
Parses a string to determine solver type.
int _poly_order
order of basis for polynomial reconstruction
TargetOperation
Available target functionals.
void setSourceExtraData(decltype(_source_extra_data) extra_data)
(OPTIONAL) Sets extra data to be used by sampling functionals in certain instances.
KOKKOS_INLINE_FUNCTION int getNumberOfTargets() const
Get number of total targets having neighborhoods (host/device).
void setTangentBundle(view_type tangent_directions)
(OPTIONAL) Sets orthonormal tangent directions for reconstruction on a manifold.
#define compadre_kernel_assert_debug(condition)
ConstraintType
Constraint type.
static int getSamplingOutputIndex(const SamplingFunctional sf, const int output_component_axis_1, const int output_component_axis_2)
Helper function for finding alpha coefficients.
KOKKOS_INLINE_FUNCTION void createWeightsAndPForCurvature(const member_type &teamMember, scratch_vector_type delta, scratch_vector_type thread_workspace, scratch_matrix_right_type P, scratch_vector_type w, const int dimension, bool only_specific_order, scratch_matrix_right_type *V=NULL) const
Fills the _P matrix with P*sqrt(w) for use in solving for curvature.