9 &&
"keep_coefficients is set to true, but number of batches exceeds 1.");
41 &&
"No target operations added to GMLS class before calling generatePolynomialCoefficients().");
58 &&
"Cannot solve GMLS problems with the NEUMANN_GRAD_SCALAR constraint using QR Factorization.");
74 }
catch(std::exception &e) {
75 printf(
"Insufficient memory to store alphas: \n\n%s", e.what());
85 std::pow(2,sro.use_target_site_weights),
90 max_num_neighbors : 1,
95 }
catch(std::exception &e) {
96 printf(
"Insufficient memory to store prestencil weights: \n\n%s", e.what());
127 int team_scratch_size_a = 0;
130 int team_scratch_size_b = 0;
131 int thread_scratch_size_a = 0;
132 int thread_scratch_size_b = 0;
150 const int max_manifold_NP = (manifold_NP >
_NP) ? manifold_NP : _NP;
159 team_scratch_size_b += scratch_vector_type::shmem_size( (
_dimensions-1)*max_num_neighbors );
161 thread_scratch_size_a += scratch_vector_type::shmem_size(this_num_cols);
162 thread_scratch_size_a += scratch_vector_type::shmem_size(
165 team_scratch_size_b += scratch_vector_type::shmem_size(max_num_neighbors);
166 team_scratch_size_b += scratch_vector_type::shmem_size(max_num_neighbors);
173 Kokkos::deep_copy(
_T, 0.0);
185 thread_scratch_size_a += scratch_vector_type::shmem_size(this_num_cols);
186 thread_scratch_size_a += scratch_vector_type::shmem_size(
199 int RHS_dim_0, RHS_dim_1;
202 int P_dim_0, P_dim_1;
213 _w = Kokkos::View<double*>(
"w", max_batch_size*
TO_GLOBAL(max_num_rows));
216 }
catch (std::exception &e) {
217 printf(
"Failed to allocate space for RHS, P, and w. Consider increasing number_of_batches: \n\n%s", e.what());
228 &&
"Normal vectors are required for solving GMLS problems with the NEUMANN_GRAD_SCALAR constraint.");
232 for (
int batch_num=0; batch_num<number_of_batches; ++batch_num) {
235 Kokkos::deep_copy(
_RHS, 0.0);
236 Kokkos::deep_copy(
_P, 0.0);
237 Kokkos::deep_copy(
_w, 0.0);
238 Kokkos::deep_copy(
_Z, 0.0);
261 Kokkos::parallel_for(
"ComputeCoarseTangentPlane", tp, functor_compute_coarse_tangent_plane);
269 Kokkos::parallel_for(
"FixTangentDirectionOrdering", tp, functor_fix_tangent_direction_ordering);
274 Kokkos::parallel_for(
"AssembleCurvaturePsqrtW", tp, functor_assemble_curvature_psqrtw);
278 Kokkos::Profiling::pushRegion(
"Curvature LU Factorization");
282 GMLS_LinearAlgebra::batchQRPivotingSolve<layout_right,layout_left,layout_right>(
_pm,
_RHS.data(), RHS_dim_0, RHS_dim_1,
_P.data(), P_dim_1, P_dim_0, manifold_NP, manifold_NP, max_num_neighbors, this_batch_size, implicit_RHS);
283 Kokkos::Profiling::popRegion();
286 Kokkos::Profiling::pushRegion(
"Curvature QR+Pivoting Factorization");
287 GMLS_LinearAlgebra::batchQRPivotingSolve<layout_right,layout_right,layout_right>(
_pm,
_P.data(), P_dim_0, P_dim_1,
_RHS.data(), RHS_dim_0, RHS_dim_1, max_num_neighbors, manifold_NP, max_num_neighbors, this_batch_size, implicit_RHS);
288 Kokkos::Profiling::popRegion();
293 Kokkos::parallel_for(
"GetAccurateTangentDirections", tp, functor_get_accurate_tangent_directions);
296 Kokkos::deep_copy(
_P, 0.0);
298 if (batch_num==number_of_batches-1) {
300 _host_T = Kokkos::create_mirror_view(
_T);
308 Kokkos::parallel_for(
"AssembleCurvaturePsqrtW", tp, functor_assemble_curvature_psqrtw);
312 Kokkos::Profiling::pushRegion(
"Curvature LU Factorization");
313 GMLS_LinearAlgebra::batchQRPivotingSolve<layout_right,layout_left,layout_right>(
_pm,
_RHS.data(), RHS_dim_0, RHS_dim_1,
_P.data(), P_dim_1, P_dim_0, manifold_NP, manifold_NP, max_num_neighbors, this_batch_size, implicit_RHS);
314 Kokkos::Profiling::popRegion();
317 Kokkos::Profiling::pushRegion(
"Curvature QR+Pivoting Factorization");
318 GMLS_LinearAlgebra::batchQRPivotingSolve<layout_right,layout_right,layout_right>(
_pm,
_P.data(), P_dim_0, P_dim_1,
_RHS.data(), RHS_dim_0, RHS_dim_1, max_num_neighbors, manifold_NP, max_num_neighbors, this_batch_size, implicit_RHS);
319 Kokkos::Profiling::popRegion();
324 Kokkos::parallel_for(
"ApplyCurvatureTargets", tp, functor_apply_curvature_targets);
332 Kokkos::parallel_for(
"ComputePrestencilWeights", tp, functor_compute_prestencil_weights);
335 Kokkos::deep_copy(
_P, 0.0);
339 Kokkos::parallel_for(
"AssembleManifoldPsqrtW", tp, functor_assemble_manifold_psqrtw);
343 Kokkos::Profiling::pushRegion(
"Manifold LU Factorization");
344 GMLS_LinearAlgebra::batchQRPivotingSolve<layout_right,layout_left,layout_right>(
_pm,
_RHS.data(), RHS_dim_0, RHS_dim_1,
_P.data(), P_dim_1, P_dim_0, this_num_cols, this_num_cols, max_num_rows, this_batch_size, implicit_RHS);
345 Kokkos::Profiling::popRegion();
347 Kokkos::Profiling::pushRegion(
"Manifold QR+Pivoting Factorization");
348 GMLS_LinearAlgebra::batchQRPivotingSolve<layout_right,layout_right,layout_right>(
_pm,
_P.data(), P_dim_0, P_dim_1,
_RHS.data(), RHS_dim_0, RHS_dim_1, max_num_rows, this_num_cols, max_num_rows, this_batch_size, implicit_RHS);
349 Kokkos::Profiling::popRegion();
362 Kokkos::parallel_for(
"AssembleStandardPsqrtW", tp, functor_assemble_standard_psqrtw);
367 Kokkos::Profiling::pushRegion(
"LU Factorization");
368 GMLS_LinearAlgebra::batchQRPivotingSolve<layout_right,layout_left,layout_right>(
_pm,
_RHS.data(), RHS_dim_0, RHS_dim_1,
_P.data(), P_dim_1, P_dim_0, this_num_cols + added_coeff_size, this_num_cols + added_coeff_size, max_num_rows +
_d_ss.
_added_alpha_size, this_batch_size, implicit_RHS);
369 Kokkos::Profiling::popRegion();
371 Kokkos::Profiling::pushRegion(
"QR+Pivoting Factorization");
373 GMLS_LinearAlgebra::batchQRPivotingSolve<layout_right,layout_right,layout_right>(
_pm,
_RHS.data(), RHS_dim_0, RHS_dim_1,
_P.data(), P_dim_1, P_dim_0, this_num_cols + added_coeff_size, this_num_cols + added_coeff_size, max_num_rows +
_d_ss.
_added_alpha_size, this_batch_size, implicit_RHS);
375 GMLS_LinearAlgebra::batchQRPivotingSolve<layout_right,layout_right,layout_right>(
_pm,
_P.data(), P_dim_0, P_dim_1,
_RHS.data(), RHS_dim_0, RHS_dim_1, max_num_rows, this_num_cols, max_num_rows, this_batch_size, implicit_RHS);
377 Kokkos::Profiling::popRegion();
381 Kokkos::parallel_for(
"ComputePrestencilWeights", tp, functor_compute_prestencil_weights);
398 Kokkos::parallel_for(
"EvaluateManifoldTargets", tp, functor_evaluate_manifold_targets);
408 Kokkos::parallel_for(
"EvaluateStandardTargets", tp, functor_evaluate_standard_targets);
415 const auto work_item_property = Kokkos::Experimental::WorkItemProperty::HintLightWeight;
416 const auto tp2 = Kokkos::Experimental::require(tp, work_item_property);
417 auto functor_apply_targets =
ApplyTargets(gmls_solution_data);
419 Kokkos::parallel_for(
"ApplyTargets", tp2, functor_apply_targets);
428 _w = Kokkos::View<double*>(
"w",0);
429 _Z = Kokkos::View<double*>(
"Z",0);
430 if (number_of_batches > 1) {
431 _RHS = Kokkos::View<double*>(
"RHS",0);
432 _P = Kokkos::View<double*>(
"P",0);
436 _RHS = Kokkos::View<double*>(
"RHS", 0);
437 if (!keep_coefficients)
_P = Kokkos::View<double*>(
"P", 0);
440 _P = Kokkos::View<double*>(
"P", 0);
441 if (!keep_coefficients)
_RHS = Kokkos::View<double*>(
"RHS", 0);
443 _RHS = Kokkos::View<double*>(
"RHS", 0);
444 if (!keep_coefficients)
_P = Kokkos::View<double*>(
"P", 0);
473 data._sampling_multiplier = gmls._sampling_multiplier;
474 data._initial_index_for_batch = gmls._initial_index_for_batch;
475 data._d_ss = gmls._d_ss;
478 const int max_num_neighbors = gmls._pc._nla.getMaxNumNeighbors();
479 const int max_num_rows = gmls._sampling_multiplier*max_num_neighbors;
487 data.this_num_cols = gmls._basis_multiplier*gmls._NP;
489 int RHS_dim_0, RHS_dim_1;
490 getRHSDims(gmls._dense_solver_type, gmls._constraint_type, gmls._reconstruction_space,
491 gmls._dimensions, max_num_rows, data.this_num_cols, RHS_dim_0, RHS_dim_1);
492 int P_dim_0, P_dim_1;
493 getPDims(gmls._dense_solver_type, gmls._constraint_type, gmls._reconstruction_space,
494 gmls._dimensions, max_num_rows, data.this_num_cols, P_dim_0, P_dim_1);
497 data.Coeffs_data = gmls._RHS.data();
498 data.Coeffs_dim_0 = RHS_dim_0;
499 data.Coeffs_dim_1 = RHS_dim_1;
501 data.Coeffs_data = gmls._P.data();
502 data.Coeffs_dim_0 = P_dim_1;
503 data.Coeffs_dim_1 = P_dim_0;
506 data.P_target_row_dim_0 = gmls._d_ss._total_alpha_values*gmls._d_ss._max_evaluation_sites_per_target;
507 data.P_target_row_dim_1 = data.this_num_cols;
508 data.P_target_row_data = gmls._Z.data();
510 data.operations_size = gmls._operations.size();
512 data.number_of_neighbors_list = gmls._pc._nla._number_of_neighbors_list;
513 data.additional_number_of_neighbors_list = gmls._additional_pc._nla._number_of_neighbors_list;
521 data._source_extra_data = gmls._source_extra_data;
522 data._target_extra_data = gmls._target_extra_data;
524 data._epsilons = gmls._epsilons ;
525 data._prestencil_weights = gmls._prestencil_weights ;
526 data._additional_pc = gmls._additional_pc;
527 data._poly_order = gmls._poly_order ;
528 data._curvature_poly_order = gmls._curvature_poly_order;
530 data._global_dimensions = gmls._global_dimensions;
531 data._local_dimensions = gmls._local_dimensions;
532 data._dimensions = gmls._dimensions;
533 data._reconstruction_space = gmls._reconstruction_space;
534 data._reconstruction_space_rank = gmls._reconstruction_space_rank;
535 data._dense_solver_type = gmls._dense_solver_type;
536 data._problem_type = gmls._problem_type;
537 data._constraint_type = gmls._constraint_type;
538 data._polynomial_sampling_functional = gmls._polynomial_sampling_functional;
539 data._data_sampling_functional = gmls._data_sampling_functional;
540 data._curvature_support_operations = gmls._curvature_support_operations;
541 data._operations = gmls._operations;
542 data._weighting_type = gmls._weighting_type;
543 data._curvature_weighting_type = gmls._curvature_weighting_type;
544 data._weighting_p = gmls._weighting_p;
545 data._weighting_n = gmls._weighting_n;
546 data._curvature_weighting_p = gmls._curvature_weighting_p;
547 data._curvature_weighting_n = gmls._curvature_weighting_n;
548 data._basis_multiplier = gmls._basis_multiplier;
549 data._sampling_multiplier = gmls._sampling_multiplier;
550 data._data_sampling_multiplier = gmls._data_sampling_multiplier;
551 data._initial_index_for_batch = gmls._initial_index_for_batch;
553 data._order_of_quadrature_points = gmls._order_of_quadrature_points;
554 data._dimension_of_quadrature_points = gmls._dimension_of_quadrature_points;
556 data._d_ss = gmls._d_ss;
558 data.max_num_neighbors = gmls._pc._nla.getMaxNumNeighbors();
559 data.max_num_rows = gmls._sampling_multiplier*data.max_num_neighbors;
560 int basis_powers_space_multiplier = (gmls._reconstruction_space ==
BernsteinPolynomial) ? 2 : 1;
562 data.manifold_NP =
GMLS::getNP(gmls._curvature_poly_order, gmls._dimensions-1,
564 data.max_manifold_NP = (data.manifold_NP > gmls._NP) ? data.manifold_NP : gmls._NP;
565 data.this_num_cols = gmls._basis_multiplier*data.max_manifold_NP;
566 data.max_poly_order = (gmls._poly_order > gmls._curvature_poly_order) ?
567 gmls._poly_order : gmls._curvature_poly_order;
569 data.ref_N_data = gmls._ref_N.data();
570 data.ref_N_dim = gmls._dimensions;
572 data.thread_workspace_dim = (data.max_poly_order+1)*gmls._global_dimensions*basis_powers_space_multiplier;
573 data.manifold_gradient_dim = (gmls._dimensions-1)*data.max_num_neighbors;
575 data.manifold_curvature_coefficients_data = gmls._manifold_curvature_coefficients.data();
578 data.manifold_NP = 0;
579 data.this_num_cols = gmls._basis_multiplier*gmls._NP;
580 data.thread_workspace_dim = (gmls._poly_order+1)*gmls._global_dimensions*basis_powers_space_multiplier;
581 data.manifold_gradient_dim = 0;
585 data.T_data = gmls._T.data();
587 data.P_target_row_dim_0 = gmls._d_ss._total_alpha_values*gmls._d_ss._max_evaluation_sites_per_target;
588 data.P_target_row_dim_1 = data.this_num_cols;
589 data.P_target_row_data = gmls._Z.data();
591 data.RHS_data = gmls._RHS.data();
592 getRHSDims(gmls._dense_solver_type, gmls._constraint_type, gmls._reconstruction_space,
593 gmls._dimensions, data.max_num_rows, data.this_num_cols, data.RHS_dim_0, data.RHS_dim_1);
595 data.P_data = gmls._P.data();
596 getPDims(gmls._dense_solver_type, gmls._constraint_type, gmls._reconstruction_space,
597 gmls._dimensions, data.max_num_rows, data.this_num_cols, data.P_dim_0, data.P_dim_1);
599 data.w_data = gmls._w.data();
602 data.Coeffs_data = gmls._RHS.data();
603 data.Coeffs_dim_0 = data.RHS_dim_0;
604 data.Coeffs_dim_1 = data.RHS_dim_1;
606 data.Coeffs_data = gmls._P.data();
607 data.Coeffs_dim_0 = data.P_dim_1;
608 data.Coeffs_dim_1 = data.P_dim_0;
Kokkos::TeamPolicy< device_execution_space > TeamPolicyThreadsAndVectors(const global_index_type batch_size, const int threads_per_team=-1, const int vector_lanes_per_thread=-1) const
Creates a team policy for a parallel_for parallel_for will break out over loops over teams with each ...
KOKKOS_INLINE_FUNCTION int getOutputDimensionOfSampling(SamplingFunctional sro, const int local_dimensions)
Dimensions ^ output rank for sampling operation (always in local chart if on a manifold, never ambient space)
std::size_t global_index_type
KOKKOS_INLINE_FUNCTION int getAdditionalAlphaSizeFromConstraint(DenseSolverType dense_solver_type, ConstraintType constraint_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::View< double * > _manifold_curvature_coefficients
curvature polynomial coefficients for all problems
Neumann Gradient Scalar Type.
SolutionSet< host_memory_space > _h_ss
Solution Set (contains all alpha values from solution and alpha layout methods)
Functor to evaluate curvature targets and apply to coefficients of curvature reconstruction.
KOKKOS_INLINE_FUNCTION int getNumberOfTargets() const
Get number of total targets having neighborhoods (host/device).
bool _reference_outward_normal_direction_provided
whether or not the reference outward normal directions were provided by the user. ...
KOKKOS_INLINE_FUNCTION int getAdditionalCoeffSizeFromConstraintAndSpace(DenseSolverType dense_solver_type, ConstraintType constraint_type, ReconstructionSpace reconstruction_space, const int dimension)
int _NP
dimension of basis for polynomial reconstruction
bool _orthonormal_tangent_space_provided
whether or not the orthonormal tangent directions were provided by the user.
const GMLSBasisData extractBasisData() const
Get GMLS basis data.
Functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity.
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
KOKKOS_INLINE_FUNCTION int calculateSamplingMultiplier(const ReconstructionSpace rs, const SamplingFunctional sro, const int local_dimensions)
Calculate sampling_multiplier.
int _sampling_multiplier
actual dimension of the sampling functional e.g.
int _curvature_poly_order
order of basis for curvature reconstruction
int _default_vector_lanes
Functor to evaluate curvature targets and construct accurate tangent direction approximation for mani...
constexpr SamplingFunctional VaryingManifoldVectorPointSample
For integrating polynomial dotted with normal over an edge.
SamplingFunctional _polynomial_sampling_functional
polynomial sampling functional used to construct P matrix, set at GMLS class instantiation NOTE: ca...
int _data_sampling_multiplier
effective dimension of the data sampling functional e.g.
Quadrature _qm
manages and calculates quadrature
Kokkos::View< double * >::HostMirror _host_T
tangent vectors information (host)
Functor to evaluate targets operations on the basis.
constexpr SamplingFunctional StaggeredEdgeAnalyticGradientIntegralSample
Analytical integral of a gradient source vector is just a difference of the scalar source at neighbor...
Kokkos::View< TargetOperation *, memory_space > _lro
vector of user requested target operations
bool _contains_valid_alphas
whether internal alpha values are valid (set externally on a solve)
bool verifyPointConnections(bool assert_valid=false)
Verify whether _pc is valid.
void generatePolynomialCoefficients(const int number_of_batches=1, const bool keep_coefficients=false, const bool clear_cache=true)
Generates polynomial coefficients by setting up and solving least squares problems Sets up the batch ...
DenseSolverType _dense_solver_type
solver type for GMLS problem - can be QR, SVD or LU
neighbor_lists_type * getNeighborLists() const
Get neighbor list accessor.
Kokkos::View< double *, layout_right, memory_space > _alphas
generated alpha coefficients (device)
Each target applies a different transform for each neighbor.
Solve GMLS problem on a manifold (will use QR or SVD to solve the resultant GMLS problem dependent on...
device_mirror_target_view_type _target_coordinates
Kokkos::View< double *****, layout_right > _prestencil_weights
generated weights for nontraditional samples required to transform data into expected sampling functi...
Functor to evaluate targets on a manifold.
global_index_type getTotalNeighborsOverAllListsHost() const
Get the sum of the number of neighbors of all targets' neighborhoods (host)
#define TO_GLOBAL(variable)
int _initial_index_for_batch
initial index for current batch
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...
Functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity.
LU factorization performed on P^T*W*P matrix.
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)
Scalar polynomial basis centered at the target site and scaled by sum of basis powers e...
QR+Pivoting factorization performed on P*sqrt(w) matrix.
point_connections_type _additional_pc
(OPTIONAL) connections between additional points and neighbors
KOKKOS_INLINE_FUNCTION int getMaxNumNeighbors() const
Get the maximum number of neighbors of all targets' neighborhoods (host/device)
int _dimension_of_quadrature_points
dimension of quadrature rule
int _added_alpha_size
additional alpha coefficients due to constraints
const GMLSSolutionData extractSolutionData() const
Get GMLS solution data.
SamplingFunctional _data_sampling_functional
generally the same as _polynomial_sampling_functional, but can differ if specified at ...
template void batchQRPivotingSolve< layout_right, layout_right, layout_right >(ParallelManager, double *, int, int, double *, int, int, int, int, int, const int, const bool)
ParallelManager _pm
determines scratch level spaces and is used to call kernels
Kokkos::View< double * > _Z
stores evaluations of targets applied to basis
int _default_threads
largest team size
ProblemType _problem_type
problem type for GMLS problem, can also be set to STANDARD for normal or MANIFOLD for manifold proble...
int _basis_multiplier
dimension of the reconstructed function e.g.
int _global_dimensions
spatial dimension of the points, set at class instantiation only
SolutionSet< device_memory_space > _d_ss
Kokkos::View< double * > _T
Rank 3 tensor for high order approximation of tangent vectors for all problems.
void setTeamScratchSize(const int level, const int value)
Functor to calculate prestencil weights to apply to data to transform into a format expected by a GML...
bool verifyAdditionalPointConnections(bool assert_valid=false)
Verify whether _additional_pc is valid.
bool _store_PTWP_inv_PTW
whether polynomial coefficients were requested to be stored (in a state not yet applied to data) ...
Functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity for curvature...
void setThreadScratchSize(const int level, const int value)
Functor to determine if tangent directions need reordered, and to reorder them if needed We require t...
ReconstructionSpace _reconstruction_space
reconstruction space for GMLS problems, set at GMLS class instantiation
Kokkos::View< double * > _P
P*sqrt(w) matrix for all problems.
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)
constexpr SamplingFunctional PointSample
Available sampling functionals.
int _dimensions
dimension of the problem, set at class instantiation only
std::string _quadrature_type
quadrature rule type
Kokkos::View< TargetOperation * >::HostMirror _host_operations
vector containing target functionals to be applied for reconstruction problem (host) ...
struct SubviewND< T, T2, enable_if_t<(T::rank< 2)> >{T _data_in;T2 _data_original_view;bool _scalar_as_vector_if_needed;SubviewND(T data_in, T2 data_original_view, bool scalar_as_vector_if_needed){_data_in=data_in;_data_original_view=data_original_view;_scalar_as_vector_if_needed=scalar_as_vector_if_needed;}auto get1DView(const int column_num) -> decltype(Kokkos::subview(_data_in, Kokkos::ALL))
Creates 1D subviews of data from a 1D view, generally constructed with CreateNDSliceOnDeviceView.
Each target applies a different data transform, but the same to each neighbor.
point_connections_type _pc
connections between points and neighbors
ConstraintType _constraint_type
constraint type for GMLS problem
int _local_dimensions
dimension of the problem, set at class instantiation only. For manifolds, generally _global_dimension...
Functor to apply target evaluation to polynomial coefficients to store in _alphas.
int _max_evaluation_sites_per_target
maximum number of evaluation sites for each target (includes target site)
int _order_of_quadrature_points
order of exact polynomial integration for quadrature rule
void generateAlphas(const int number_of_batches=1, const bool keep_coefficients=false, const bool clear_cache=true)
Meant to calculate target operations and apply the evaluations to the previously constructed polynomi...
#define compadre_assert_release(condition)
compadre_assert_release is used for assertions that should always be checked, but generally are not e...
Functor to create a coarse tangent approximation from a given neighborhood of points.
NeighborLists< Kokkos::View< int * > > _neighbor_lists
Accessor to get neighbor list data, offset data, and number of neighbors per target.
int _total_alpha_values
used for sizing P_target_row and the _alphas view
int _poly_order
order of basis for polynomial reconstruction
Bernstein polynomial basis.
template void batchQRPivotingSolve< layout_right, layout_left, layout_right >(ParallelManager, double *, int, int, double *, int, int, int, int, int, const int, const bool)
KOKKOS_INLINE_FUNCTION int calculateBasisMultiplier(const ReconstructionSpace rs, const int local_dimensions)
Calculate basis_multiplier.