Compadre  1.2.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Classes | Protected Attributes | List of all members
Compadre::GMLS Class Reference

Generalized Moving Least Squares (GMLS) More...

Detailed Description

Generalized Moving Least Squares (GMLS)

This class sets up a batch of GMLS problems from a given set of neighbor lists, target sites, and source sites. GMLS requires a target functional, reconstruction space, and sampling functional to be specified. For a given choice of reconstruction space and sampling functional, multiple targets can be generated with very little additional computation, which is why this class allows for multiple target functionals to be specified.

Definition at line 25 of file Compadre_GMLS.hpp.

#include <Compadre_GMLS.hpp>

Classes

struct  ApplyCurvatureTargets
 Tag for functor to evaluate curvature targets and apply to coefficients of curvature reconstruction. More...
 
struct  ApplyManifoldTargets
 Tag for functor to evaluate targets, apply target evaluation to polynomial coefficients to store in _alphas. More...
 
struct  ApplyStandardTargets
 Tag for functor to evaluate targets, apply target evaluation to polynomial coefficients to store in _alphas. More...
 
struct  AssembleCurvaturePsqrtW
 Tag for functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity for curvature. More...
 
struct  AssembleManifoldPsqrtW
 Tag for functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity. More...
 
struct  AssembleStandardPsqrtW
 Tag for functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity. More...
 
struct  ComputeCoarseTangentPlane
 Tag for functor to create a coarse tangent approximation from a given neighborhood of points. More...
 
struct  ComputePrestencilWeights
 Tag for functor to calculate prestencil weights to apply to data to transform into a format expected by a GMLS stencil. More...
 
struct  FixTangentDirectionOrdering
 Tag for functor to determine if tangent directions need reordered, and to reorder them if needed. More...
 
struct  GetAccurateTangentDirections
 Tag for functor to evaluate curvature targets and construct accurate tangent direction approximation for manifolds. More...
 

Public Member Functions

Instantiation / Destruction
 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) More...
 
 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 functional. More...
 
 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 used. More...
 
 ~GMLS ()
 Destructor. More...
 
Functors

Member functions that perform operations on the entire batch

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. More...
 
KOKKOS_INLINE_FUNCTION void operator() (const ApplyStandardTargets &, const member_type &teamMember) const
 Functor to evaluate targets, apply target evaluation to polynomial coefficients to store in _alphas. More...
 
KOKKOS_INLINE_FUNCTION void operator() (const ComputeCoarseTangentPlane &, const member_type &teamMember) const
 Functor to create a coarse tangent approximation from a given neighborhood of points. More...
 
KOKKOS_INLINE_FUNCTION void operator() (const AssembleCurvaturePsqrtW &, const member_type &teamMember) const
 Functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity for curvature. More...
 
KOKKOS_INLINE_FUNCTION void operator() (const GetAccurateTangentDirections &, const member_type &teamMember) const
 Functor to evaluate curvature targets and construct accurate tangent direction approximation for manifolds. More...
 
KOKKOS_INLINE_FUNCTION void operator() (const FixTangentDirectionOrdering &, const member_type &teamMember) const
 Functor to determine if tangent directions need reordered, and to reorder them if needed We require that the normal is consistent with a right hand rule on the tangent vectors. More...
 
KOKKOS_INLINE_FUNCTION void operator() (const ApplyCurvatureTargets &, const member_type &teamMember) const
 Functor to evaluate curvature targets and apply to coefficients of curvature reconstruction. More...
 
KOKKOS_INLINE_FUNCTION void operator() (const AssembleManifoldPsqrtW &, const member_type &teamMember) const
 Functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity. More...
 
KOKKOS_INLINE_FUNCTION void operator() (const ApplyManifoldTargets &, const member_type &teamMember) const
 Functor to evaluate targets, apply target evaluation to polynomial coefficients to store in _alphas. More...
 
KOKKOS_INLINE_FUNCTION void operator() (const ComputePrestencilWeights &, const member_type &teamMember) const
 Functor to calculate prestencil weights to apply to data to transform into a format expected by a GMLS stencil. More...
 
Accessors

Retrieve member variables through public member functions

host_managed_local_index_type getPolynomialCoefficientsDomainRangeSize () const
 Returns (size of the basis used in instance's polynomial reconstruction) x (data input dimension) More...
 
int getPolynomialCoefficientsSize () const
 Returns size of the basis used in instance's polynomial reconstruction. More...
 
host_managed_local_index_type getPolynomialCoefficientsMemorySize () const
 Returns 2D array size in memory on which coefficients are stored. More...
 
int getDimensions () const
 Dimension of the GMLS problem, set only at class instantiation. More...
 
int getGlobalDimensions () const
 Dimension of the GMLS problem's point data (spatial description of points in ambient space), set only at class instantiation. More...
 
int getLocalDimensions () const
 Local dimension of the GMLS problem (less than global dimension if on a manifold), set only at class instantiation. More...
 
DenseSolverType getDenseSolverType ()
 Get dense solver type. More...
 
ProblemType getProblemType ()
 Get problem type. More...
 
ConstraintType getConstraintType ()
 Get constraint type. More...
 
WeightingFunctionType getWeightingType () const
 Type for weighting kernel for GMLS problem. More...
 
WeightingFunctionType getManifoldWeightingType () const
 Type for weighting kernel for curvature. More...
 
int getWeightingPower () const
 Power for weighting kernel for GMLS problem. More...
 
int getManifoldWeightingPower () const
 Power for weighting kernel for curvature. More...
 
int getNumberOfQuadraturePoints () const
 Number of quadrature points. More...
 
int getOrderOfQuadraturePoints () const
 Order of quadrature points. More...
 
int getDimensionOfQuadraturePoints () const
 Dimensions of quadrature points. More...
 
std::string getQuadratureType () const
 Type of quadrature points. More...
 
decltype(_neighbor_lists)* getNeighborLists ()
 Get neighbor list accessor. More...
 
decltype(_TgetTangentDirections () const
 Get a view (device) of all tangent direction bundles. More...
 
decltype(_ref_NgetReferenceNormalDirections () const
 Get a view (device) of all reference outward normal directions. More...
 
double getTangentBundle (const int target_index, const int direction, const int component) const
 Get component of tangent or normal directions for manifold problems. More...
 
double getReferenceNormalDirection (const int target_index, const int component) const
 Get component of tangent or normal directions for manifold problems. More...
 
int getTargetOperationLocalIndex (TargetOperation lro) const
 Get the local index (internal) to GMLS for a particular TargetOperation Every TargetOperation has a global index which can be readily found in Compadre::TargetOperation but this function returns the index used inside of the GMLS class. More...
 
decltype(_prestencil_weightsgetPrestencilWeights () const
 Get a view (device) of all rank 2 preprocessing tensors This is a rank 5 tensor that is able to provide data transformation into a form that GMLS is able to operate on. More...
 
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 multiplied by the number of neighbors for each row and then needs a neighbor number added to this returned value to be meaningful) More...
 
decltype(_alphasgetAlphas () const
 Get a view (device) of all alphas. More...
 
decltype(_RHSgetFullPolynomialCoefficientsBasis () const
 Get a view (device) of all polynomial coefficients basis. More...
 
SamplingFunctional getPolynomialSamplingFunctional () const
 Get the polynomial sampling functional specified at instantiation. More...
 
SamplingFunctional getDataSamplingFunctional () const
 Get the data sampling functional specified at instantiation (often the same as the polynomial sampling functional) More...
 
ReconstructionSpace getReconstructionSpace () const
 Get the reconstruction space specified at instantiation. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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 access into alphas from a rank 1 view into a rank 3 view. More...
 
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 access into alphas from a rank 1 view into a rank 3 view. More...
 
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. More...
 
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 functionals. More...
 
int getOutputDimensionOfOperation (TargetOperation lro, bool ambient=false) const
 Dimensions ^ output rank for target operation. More...
 
int getInputDimensionOfOperation (TargetOperation lro) const
 Dimensions ^ input rank for target operation (always in local chart if on a manifold, never ambient space) More...
 
int getOutputDimensionOfSampling (SamplingFunctional sro) const
 Dimensions ^ output rank for sampling operation (always in local chart if on a manifold, never ambient space) More...
 
int getInputDimensionOfSampling (SamplingFunctional sro) const
 Dimensions ^ output rank for sampling operation (always in ambient space, never local chart on a manifold) More...
 
int calculateBasisMultiplier (const ReconstructionSpace rs) const
 Calculate basis_multiplier. More...
 
int calculateSamplingMultiplier (const ReconstructionSpace rs, const SamplingFunctional sro) const
 Calculate sampling_multiplier. More...
 
Modifiers

Changed member variables through public member functions

void resetCoefficientData ()
 
template<typename view_type_1 , typename view_type_2 , typename view_type_3 , typename view_type_4 >
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) More...
 
template<typename view_type_1 , typename view_type_2 , typename view_type_3 , typename view_type_4 >
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) More...
 
template<typename view_type_1 , typename view_type_2 >
void setAdditionalEvaluationSitesData (view_type_1 additional_evaluation_indices, view_type_2 additional_evaluation_coordinates)
 (OPTIONAL) Sets additional evaluation sites for each target site More...
 
template<typename view_type >
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). More...
 
template<typename view_type >
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). More...
 
template<typename view_type >
std::enable_if
< view_type::rank==2, void >
::type 
setNeighborLists (view_type neighbor_lists)
 Sets neighbor list information. More...
 
template<typename view_type >
void setSourceSites (view_type source_coordinates)
 Sets source coordinate information. More...
 
template<typename view_type >
void setSourceSites (decltype(_source_coordinates) source_coordinates)
 Sets source coordinate information. More...
 
template<typename view_type >
void setTargetSites (view_type target_coordinates)
 Sets target coordinate information. Rows of this 2D-array should correspond to rows of the neighbor lists. More...
 
template<typename view_type >
void setTargetSites (decltype(_target_coordinates) target_coordinates)
 Sets target coordinate information. Rows of this 2D-array should correspond to rows of the neighbor lists. More...
 
template<typename view_type >
void setWindowSizes (view_type epsilons)
 Sets window sizes, also called the support of the kernel. More...
 
template<typename view_type >
void setWindowSizes (decltype(_epsilons) epsilons)
 Sets window sizes, also called the support of the kernel (device) More...
 
template<typename view_type >
void setTangentBundle (view_type tangent_directions)
 (OPTIONAL) Sets orthonormal tangent directions for reconstruction on a manifold. More...
 
template<typename view_type >
void setReferenceOutwardNormalDirection (view_type outward_normal_directions, bool use_to_orient_surface=true)
 (OPTIONAL) Sets outward normal direction. More...
 
template<typename view_type >
void setSourceExtraData (view_type extra_data)
 (OPTIONAL) Sets extra data to be used by sampling functionals in certain instances. More...
 
template<typename view_type >
void setSourceExtraData (decltype(_source_extra_data) extra_data)
 (OPTIONAL) Sets extra data to be used by sampling functionals in certain instances. More...
 
template<typename view_type >
void setTargetExtraData (view_type extra_data)
 (OPTIONAL) Sets extra data to be used by target operations in certain instances. More...
 
template<typename view_type >
void setTargetExtraData (decltype(_target_extra_data) extra_data)
 (OPTIONAL) Sets extra data to be used by target operations in certain instances. More...
 
template<typename view_type >
void setAuxiliaryEvaluationCoordinates (view_type evaluation_coordinates)
 (OPTIONAL) Sets additional points for evaluation of target operation on polynomial reconstruction. More...
 
template<typename view_type >
void setAuxiliaryEvaluationCoordinates (decltype(_additional_evaluation_coordinates) evaluation_coordinates)
 (OPTIONAL) Sets additional points for evaluation of target operation on polynomial reconstruction. More...
 
template<typename view_type >
void setAuxiliaryEvaluationIndicesLists (view_type indices_lists)
 (OPTIONAL) Sets the additional target evaluation coordinate indices list information. More...
 
template<typename view_type >
void setAuxiliaryEvaluationIndicesLists (decltype(_additional_evaluation_indices) indices_lists)
 (OPTIONAL) Sets the additional target evaluation coordinate indices list information. More...
 
void setWeightingType (const std::string &wt)
 Type for weighting kernel for GMLS problem. More...
 
void setWeightingType (const WeightingFunctionType wt)
 Type for weighting kernel for GMLS problem. More...
 
void setCurvatureWeightingType (const std::string &wt)
 Type for weighting kernel for curvature. More...
 
void setCurvatureWeightingType (const WeightingFunctionType wt)
 Type for weighting kernel for curvature. More...
 
void setPolynomialOrder (const int poly_order)
 Sets basis order to be used when reoncstructing any function. More...
 
void setCurvaturePolynomialOrder (const int manifold_poly_order)
 Sets basis order to be used when reoncstructing curvature. More...
 
void setWeightingPower (int wp)
 Power for weighting kernel for GMLS problem. More...
 
void setCurvatureWeightingPower (int wp)
 Power for weighting kernel for curvature. More...
 
void setOrderOfQuadraturePoints (int order)
 Number quadrature points to use. More...
 
void setDimensionOfQuadraturePoints (int dim)
 Dimensions of quadrature points to use. More...
 
void setQuadratureType (std::string quadrature_type)
 Type of quadrature points. More...
 
void addTargets (TargetOperation lro)
 Adds a target to the vector of target functional to be applied to the reconstruction. More...
 
void addTargets (std::vector< TargetOperation > lro)
 Adds a vector of target functionals to the vector of target functionals already to be applied to the reconstruction. More...
 
void clearTargets ()
 Empties the vector of target functionals to apply to the reconstruction. More...
 
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 batch of GMLS problems to be solved for. Provides alpha values ! that can later be contracted against data or degrees of freedom to form a ! global linear system. ! More...
 
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 polynomial coefficients. But now that is inside of generatePolynomialCoefficients because ! it must be to handle number_of_batches>1. Effectively, this just calls generatePolynomialCoefficients. ! More...
 

Static Public Member Functions

Public Utility
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..3 and polynomial order m The divfree options will return the divergence-free basis if true. More...
 
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. More...
 
static KOKKOS_INLINE_FUNCTION
double 
Wab (const double r, const double h, const WeightingFunctionType &weighting_type, const int power)
 Evaluates the weighting kernel. More...
 
static KOKKOS_INLINE_FUNCTION
double 
EuclideanVectorLength (const XYZ &delta_vector, const int dimension)
 Returns Euclidean norm of a vector. More...
 
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. More...
 
static int getSamplingOutputIndex (const SamplingFunctional sf, const int output_component_axis_1, const int output_component_axis_2)
 Helper function for finding alpha coefficients. More...
 
static int getOutputRankOfSampling (SamplingFunctional sro)
 Output rank for sampling operation. More...
 
static int getInputRankOfSampling (SamplingFunctional sro)
 Input rank for sampling operation. More...
 

Protected Member Functions

Private Modifiers

Private function because information lives on the device

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 P. More...
 
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. More...
 
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. More...
 
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) More...
 
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. More...
 
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. More...
 
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. More...
 
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 curvature. More...
 
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. More...
 
Private Accessors

Private function because information lives on the device

KOKKOS_INLINE_FUNCTION int getNNeighbors (const int target_index) const
 Returns number of neighbors for a particular target. More...
 
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 that neighbor. More...
 
KOKKOS_INLINE_FUNCTION int getMaxNNeighbors () const
 Returns the maximum neighbor lists size over all target sites. More...
 
KOKKOS_INLINE_FUNCTION int getMaxEvaluationSitesPerTarget () const
 Returns the maximum number of evaluation sites over all target sites (target sites are included in total) More...
 
KOKKOS_INLINE_FUNCTION int getNEvaluationSitesPerTarget (const int target_index) const
 (OPTIONAL) Returns number of additional evaluation sites for a particular target More...
 
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 contains the coordinates for that evaluation More...
 
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. More...
 
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. More...
 
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. More...
 
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. More...
 
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 orthonormal basis V. More...
 
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 orthonormal basis V^T. More...
 
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. More...
 
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. More...
 

Static Protected Member Functions

Private Utility
static DenseSolverType parseSolverType (const std::string &dense_solver_type)
 Parses a string to determine solver type. More...
 
static ProblemType parseProblemType (const std::string &problem_type)
 Parses a string to determine problem type. More...
 
static ConstraintType parseConstraintType (const std::string &constraint_type)
 Parses a string to determine constraint type. More...
 

Protected Attributes

pool_type _random_number_pool
 
Kokkos::View< double * > _w
 contains weights for all problems More...
 
Kokkos::View< double * > _P
 P*sqrt(w) matrix for all problems. More...
 
Kokkos::View< double * > _RHS
 sqrt(w)*Identity matrix for all problems, later holds polynomial coefficients for all problems More...
 
Kokkos::View< double * > _T
 Rank 3 tensor for high order approximation of tangent vectors for all problems. More...
 
Kokkos::View< double * > _ref_N
 Rank 2 tensor for high order approximation of tangent vectors for all problems. More...
 
Kokkos::View< double * >
::HostMirror 
_host_T
 tangent vectors information (host) More...
 
Kokkos::View< double * >
::HostMirror 
_host_ref_N
 reference outward normal vectors information (host) More...
 
Kokkos::View< double * > _manifold_metric_tensor_inverse
 metric tensor inverse for all problems More...
 
Kokkos::View< double * > _manifold_curvature_coefficients
 curvature polynomial coefficients for all problems More...
 
Kokkos::View< double * > _manifold_curvature_gradient
 _dimension-1 gradient values for curvature for all problems More...
 
Kokkos::View< double
**, layout_right
_source_extra_data
 Extra data available to basis functions (optional) More...
 
Kokkos::View< double
**, layout_right
_target_extra_data
 Extra data available to target operations (optional) More...
 
NeighborLists< Kokkos::View
< int * > > 
_neighbor_lists
 Accessor to get neighbor list data, offset data, and number of neighbors per target. More...
 
Kokkos::View< int
*, host_memory_space
_host_number_of_neighbors_list
 convenient copy on host of number of neighbors More...
 
Kokkos::View< double
**, layout_right
_source_coordinates
 all coordinates for the source for which _neighbor_lists refers (device) More...
 
Kokkos::View< double
**, layout_right
_target_coordinates
 coordinates for target sites for reconstruction (device) More...
 
Kokkos::View< double * > _epsilons
 h supports determined through neighbor search (device) More...
 
Kokkos::View< double * >
::HostMirror 
_host_epsilons
 h supports determined through neighbor search (host) More...
 
Kokkos::View< double
*, layout_right
_alphas
 generated alpha coefficients (device) More...
 
Kokkos::View< const double
*, layout_right >::HostMirror 
_host_alphas
 generated alpha coefficients (host) More...
 
Kokkos::View< double
*****, layout_right
_prestencil_weights
 generated weights for nontraditional samples required to transform data into expected sampling functional form (device). More...
 
Kokkos::View< const double
*****, layout_right >
::HostMirror 
_host_prestencil_weights
 generated weights for nontraditional samples required to transform data into expected sampling functional form (host) More...
 
Kokkos::View< double
**, layout_right
_additional_evaluation_coordinates
 (OPTIONAL) user provided additional coordinates for target operation evaluation (device) More...
 
Kokkos::View< int
**, layout_right
_additional_evaluation_indices
 (OPTIONAL) contains indices of entries in the _additional_evaluation_coordinates view (device) More...
 
Kokkos::View< int
**, layout_right >::HostMirror 
_host_additional_evaluation_indices
 (OPTIONAL) contains indices of entries in the _additional_evaluation_coordinates view (host) More...
 
Kokkos::View< int * > _number_of_additional_evaluation_indices
 (OPTIONAL) contains the # of additional coordinate indices for each target More...
 
int _poly_order
 order of basis for polynomial reconstruction More...
 
int _curvature_poly_order
 order of basis for curvature reconstruction More...
 
int _NP
 dimension of basis for polynomial reconstruction More...
 
int _global_dimensions
 spatial dimension of the points, set at class instantiation only More...
 
int _local_dimensions
 dimension of the problem, set at class instantiation only. For manifolds, generally _global_dimensions-1 More...
 
int _dimensions
 dimension of the problem, set at class instantiation only More...
 
ReconstructionSpace _reconstruction_space
 reconstruction space for GMLS problems, set at GMLS class instantiation More...
 
int _reconstruction_space_rank
 actual rank of reconstruction basis More...
 
DenseSolverType _dense_solver_type
 solver type for GMLS problem - can be QR, SVD or LU More...
 
ProblemType _problem_type
 problem type for GMLS problem, can also be set to STANDARD for normal or MANIFOLD for manifold problems More...
 
ConstraintType _constraint_type
 constraint type for GMLS problem More...
 
const SamplingFunctional _polynomial_sampling_functional
 polynomial sampling functional used to construct P matrix, set at GMLS class instantiation More...
 
const SamplingFunctional _data_sampling_functional
 generally the same as _polynomial_sampling_functional, but can differ if specified at GMLS class instantiation More...
 
Kokkos::View< TargetOperation * > _curvature_support_operations
 vector containing target functionals to be applied for curvature More...
 
Kokkos::View< TargetOperation * > _operations
 vector containing target functionals to be applied for reconstruction problem (device) More...
 
Kokkos::View< TargetOperation * >
::HostMirror 
_host_operations
 vector containing target functionals to be applied for reconstruction problem (host) More...
 
WeightingFunctionType _weighting_type
 weighting kernel type for GMLS More...
 
WeightingFunctionType _curvature_weighting_type
 weighting kernel type for curvature problem More...
 
int _weighting_power
 power to be used for weighting kernel More...
 
int _curvature_weighting_power
 power to be used for weighting kernel for curvature More...
 
int _basis_multiplier
 dimension of the reconstructed function e.g. More...
 
int _sampling_multiplier
 actual dimension of the sampling functional e.g. More...
 
int _data_sampling_multiplier
 effective dimension of the data sampling functional e.g. More...
 
bool _nontrivial_nullspace
 whether or not operator to be inverted for GMLS problem has a nontrivial nullspace (requiring SVD) More...
 
bool _orthonormal_tangent_space_provided
 whether or not the orthonormal tangent directions were provided by the user. More...
 
bool _reference_outward_normal_direction_provided
 whether or not the reference outward normal directions were provided by the user. More...
 
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. More...
 
bool _entire_batch_computed_at_once
 whether entire calculation was computed at once the alternative is that it was broken up over many smaller groups, in which case this is false, and so the _RHS matrix can not be stored or requested More...
 
bool _store_PTWP_inv_PTW
 whether polynomial coefficients were requested to be stored (in a state not yet applied to data) More...
 
int _initial_index_for_batch
 initial index for current batch More...
 
int _max_num_neighbors
 maximum number of neighbors over all target sites More...
 
int _max_evaluation_sites_per_target
 maximum number of evaluation sites for each target (includes target site) More...
 
std::vector< TargetOperation_lro
 vector of user requested target operations More...
 
std::vector< int > _lro_lookup
 vector containing a mapping from a target functionals enum value to the its place in the list of target functionals to be applied More...
 
Kokkos::View< int * > _lro_total_offsets
 index for where this operation begins the for _alpha coefficients (device) More...
 
Kokkos::View< int * >::HostMirror _host_lro_total_offsets
 index for where this operation begins the for _alpha coefficients (host) More...
 
Kokkos::View< int * > _lro_output_tile_size
 dimensions ^ rank of tensor of output for each target functional (device) More...
 
Kokkos::View< int * >::HostMirror _host_lro_output_tile_size
 dimensions ^ rank of tensor of output for each target functional (host) More...
 
Kokkos::View< int * > _lro_input_tile_size
 dimensions ^ rank of tensor of output for each sampling functional (device) More...
 
Kokkos::View< int * >::HostMirror _host_lro_input_tile_size
 dimensions ^ rank of tensor of output for each sampling functional (host) More...
 
Kokkos::View< int * > _lro_output_tensor_rank
 tensor rank of target functional (device) More...
 
Kokkos::View< int * >::HostMirror _host_lro_output_tensor_rank
 tensor rank of target functional (host) More...
 
Kokkos::View< int * > _lro_input_tensor_rank
 tensor rank of sampling functional (device) More...
 
Kokkos::View< int * >::HostMirror _host_lro_input_tensor_rank
 tensor rank of sampling functional (host) More...
 
int _total_alpha_values
 used for sizing P_target_row and the _alphas view More...
 
int _added_alpha_size
 additional alpha coefficients due to constraints More...
 
ParallelManager _pm
 determines scratch level spaces and is used to call kernels More...
 
int _order_of_quadrature_points
 order of exact polynomial integration for quadrature rule More...
 
int _dimension_of_quadrature_points
 dimension of quadrature rule More...
 
std::string _quadrature_type
 quadrature rule type More...
 
Quadrature _qm
 manages and calculates quadrature More...
 

Constructor & Destructor Documentation

Compadre::GMLS::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 
)
inline

Minimal constructor providing no data (neighbor lists, source sites, target sites)

Definition at line 625 of file Compadre_GMLS.hpp.

Compadre::GMLS::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 
)
inline

Constructor for the case when the data sampling functional does not match the polynomial sampling functional.

Only case anticipated is staggered Laplacian.

Definition at line 705 of file Compadre_GMLS.hpp.

Compadre::GMLS::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 
)
inline

Constructor for the case when nonstandard sampling functionals or reconstruction spaces are to be used.

Reconstruction space and sampling strategy can only be set at instantiation.

Definition at line 715 of file Compadre_GMLS.hpp.

Compadre::GMLS::~GMLS ( )
inline

Destructor.

Definition at line 726 of file Compadre_GMLS.hpp.

Member Function Documentation

void Compadre::GMLS::addTargets ( TargetOperation  lro)
inline

Adds a target to the vector of target functional to be applied to the reconstruction.

Definition at line 1787 of file Compadre_GMLS.hpp.

void Compadre::GMLS::addTargets ( std::vector< TargetOperation lro)
inline

Adds a vector of target functionals to the vector of target functionals already to be applied to the reconstruction.

Definition at line 1794 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::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
protected

Helper function for applying the evaluations from a target functional to the polynomial coefficients.

Definition at line 8 of file Compadre_GMLS_ApplyTargetEvaluations.hpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::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
protected

Evaluates the gradient of a polynomial basis under the Dirac Delta (pointwise) sampling function.

Parameters
delta[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large is the _basis_multipler*the dimension of the polynomial basis.
thread_workspace[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _poly_order*the spatial dimension of the polynomial basis.
target_index[in] - target number
neighbor_index[in] - index of neighbor for this target with respect to local numbering [0,...,number of neighbors for target]
alpha[in] - double to determine convex combination of target and neighbor site at which to evaluate polynomials. (1-alpha)*neighbor + alpha*target
partial_direction[in] - direction that partial is taken with respect to, e.g. 0 is x direction, 1 is y direction
dimension[in] - spatial dimension of basis to evaluate. e.g. dimension two basis of order one is 1, x, y, whereas for dimension 3 it is 1, x, y, z
poly_order[in] - polynomial basis degree
specific_order_only[in] - boolean for only evaluating one degree of polynomial when true
V[in] - orthonormal basis matrix size _dimensions * _dimensions whose first _dimensions-1 columns are an approximation of the tangent plane
reconstruction_space[in] - space of polynomial that a sampling functional is to evaluate
sampling_strategy[in] - sampling functional specification
additional_evaluation_local_index[in] - local index for evaluation sites

Definition at line 468 of file Compadre_GMLS_Basis.hpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::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
protected

Evaluates the Hessian of a polynomial basis under the Dirac Delta (pointwise) sampling function.

Parameters
delta[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large is the _basis_multipler*the dimension of the polynomial basis.
thread_workspace[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _poly_order*the spatial dimension of the polynomial basis.
target_index[in] - target number
neighbor_index[in] - index of neighbor for this target with respect to local numbering [0,...,number of neighbors for target]
alpha[in] - double to determine convex combination of target and neighbor site at which to evaluate polynomials. (1-alpha)*neighbor + alpha*target
partial_direction_1[in] - first direction that partial is taken with respect to, e.g. 0 is x direction, 1 is y direction
partial_direction_2[in] - second direction that partial is taken with respect to, e.g. 0 is x direction, 1 is y direction
dimension[in] - spatial dimension of basis to evaluate. e.g. dimension two basis of order one is 1, x, y, whereas for dimension 3 it is 1, x, y, z
poly_order[in] - polynomial basis degree
specific_order_only[in] - boolean for only evaluating one degree of polynomial when true
V[in] - orthonormal basis matrix size _dimensions * _dimensions whose first _dimensions-1 columns are an approximation of the tangent plane
reconstruction_space[in] - space of polynomial that a sampling functional is to evaluate
sampling_strategy[in] - sampling functional specification
additional_evaluation_local_index[in] - local index for evaluation sites

Definition at line 531 of file Compadre_GMLS_Basis.hpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::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
protected

Evaluates the polynomial basis under a particular sampling function. Generally used to fill a row of P.

Parameters
delta[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _basis_multipler*the dimension of the polynomial basis.
thread_workspace[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _poly_order*the spatial dimension of the polynomial basis.
target_index[in] - target number
neighbor_index[in] - index of neighbor for this target with respect to local numbering [0,...,number of neighbors for target]
alpha[in] - double to determine convex combination of target and neighbor site at which to evaluate polynomials. (1-alpha)*neighbor + alpha*target
dimension[in] - spatial dimension of basis to evaluate. e.g. dimension two basis of order one is 1, x, y, whereas for dimension 3 it is 1, x, y, z
poly_order[in] - polynomial basis degree
specific_order_only[in] - boolean for only evaluating one degree of polynomial when true
V[in] - orthonormal basis matrix size _dimensions * _dimensions whose first _dimensions-1 columns are an approximation of the tangent plane
reconstruction_space[in] - space of polynomial that a sampling functional is to evaluate
sampling_strategy[in] - sampling functional specification
additional_evaluation_local_index[in] - local index for evaluation sites

Definition at line 9 of file Compadre_GMLS_Basis.hpp.

int Compadre::GMLS::calculateBasisMultiplier ( const ReconstructionSpace  rs) const
inline

Calculate basis_multiplier.

Definition at line 1252 of file Compadre_GMLS.hpp.

int Compadre::GMLS::calculateSamplingMultiplier ( const ReconstructionSpace  rs,
const SamplingFunctional  sro 
) const
inline

Calculate sampling_multiplier.

Definition at line 1259 of file Compadre_GMLS.hpp.

void Compadre::GMLS::clearTargets ( )
inline

Empties the vector of target functionals to apply to the reconstruction.

Definition at line 1879 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::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
protected

Evaluates a polynomial basis for the curvature with a gradient target functional applied.

_operations is used by this function which is set through a modifier function

Parameters
teamMember[in] - Kokkos::TeamPolicy member type (created by parallel_for)
delta[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large is the _basis_multipler*the dimension of the polynomial basis.
thread_workspace[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _curvature_poly_order*the spatial dimension of the polynomial basis.
P_target_row[out] - 1D Kokkos View where the evaluation of the polynomial basis is stored
V[in] - orthonormal basis matrix size _dimensions * _dimensions whose first _dimensions-1 columns are an approximation of the tangent plane

Definition at line 951 of file Compadre_GMLS_Targets.hpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::computeTargetFunctionals ( const member_type teamMember,
scratch_vector_type  t1,
scratch_vector_type  t2,
scratch_matrix_right_type  P_target_row 
) const
protected

Evaluates a polynomial basis with a target functional applied to each member of the basis.

Parameters
teamMember[in] - Kokkos::TeamPolicy member type (created by parallel_for)
delta[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large is the _basis_multipler*the dimension of the polynomial basis.
thread_workspace[in/out] - scratch space that is allocated so that each team has its own copy. Must be at least as large is the _poly_order*_global_dimensions.
P_target_row[out] - 1D Kokkos View where the evaluation of the polynomial basis is stored

Definition at line 10 of file Compadre_GMLS_Targets.hpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::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
protected

Evaluates a polynomial basis with a target functional applied, using information from the manifold curvature.

_operations is used by this function which is set through a modifier function

Parameters
teamMember[in] - Kokkos::TeamPolicy member type (created by parallel_for)
delta[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large is the _basis_multipler*the dimension of the polynomial basis.
thread_workspace[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _curvature_poly_order*the spatial dimension of the polynomial basis.
P_target_row[out] - 1D Kokkos View where the evaluation of the polynomial basis is stored
V[in] - orthonormal basis matrix size _dimensions * _dimensions whose first _dimensions-1 columns are an approximation of the tangent plane
G_inv[in] - (_dimensions-1)*(_dimensions-1) Kokkos View containing inverse of metric tensor
curvature_coefficients[in] - polynomial coefficients for curvature
curvature_gradients[in] - approximation of gradient of curvature, Kokkos View of size (_dimensions-1)

Definition at line 1001 of file Compadre_GMLS_Targets.hpp.

KOKKOS_INLINE_FUNCTION double Compadre::GMLS::convertGlobalToLocalCoordinate ( const XYZ  global_coord,
const int  dim,
const scratch_matrix_right_type V 
) const
inlineprotected

Returns a component of the local coordinate after transformation from global to local under the orthonormal basis V.

Definition at line 528 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION double Compadre::GMLS::convertLocalToGlobalCoordinate ( const XYZ  local_coord,
const int  dim,
const scratch_matrix_right_type V 
) const
inlineprotected

Returns a component of the global coordinate after transformation from local to global under the orthonormal basis V^T.

Definition at line 539 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::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
protected

Fills the _P matrix with either P or P*sqrt(w)

Parameters
teamMember[in] - Kokkos::TeamPolicy member type (created by parallel_for)
delta[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large is the _basis_multipler*the dimension of the polynomial basis.
thread_workspace[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _poly_order*the spatial dimension of the polynomial basis.
P[out] - 2D Kokkos View which will contain evaluation of sampling functional on polynomial basis for each neighbor the target has (stored column major)
w[out] - 1D Kokkos View which will contain weighting kernel values for the target with each neighbor if weight_p = true
dimension[in] - spatial dimension of basis to evaluate. e.g. dimension two basis of order one is 1, x, y, whereas for dimension 3 it is 1, x, y, z
polynomial_order[in] - polynomial basis degree
weight_p[in] - boolean whether to fill w with kernel weights
V[in] - orthonormal basis matrix size _dimensions * _dimensions whose first _dimensions-1 columns are an approximation of the tangent plane
reconstruction_space[in] - space of polynomial that a sampling functional is to evaluate
sampling_strategy[in] - sampling functional specification

Definition at line 592 of file Compadre_GMLS_Basis.hpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::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
protected

Fills the _P matrix with P*sqrt(w) for use in solving for curvature.

Uses _curvature_poly_order as the polynomial order of the basis

Parameters
teamMember[in] - Kokkos::TeamPolicy member type (created by parallel_for)
delta[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large is the _basis_multipler*the dimension of the polynomial basis.
thread_workspace[in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _poly_order*the spatial dimension of the polynomial basis.
P[out] - 2D Kokkos View which will contain evaluation of sampling functional on polynomial basis for each neighbor the target has (stored column major)
w[out] - 1D Kokkos View which will contain weighting kernel values for the target with each neighbor if weight_p = true
dimension[in] - spatial dimension of basis to evaluate. e.g. dimension two basis of order one is 1, x, y, whereas for dimension 3 it is 1, x, y, z
only_specific_order[in] - boolean for only evaluating one degree of polynomial when true
V[in] - orthonormal basis matrix size _dimensions * _dimensions whose first _dimensions-1 columns are an approximation of the tangent plane

Definition at line 667 of file Compadre_GMLS_Basis.hpp.

static KOKKOS_INLINE_FUNCTION double Compadre::GMLS::EuclideanVectorLength ( const XYZ delta_vector,
const int  dimension 
)
inlinestatic

Returns Euclidean norm of a vector.

Definition at line 873 of file Compadre_GMLS.hpp.

void Compadre::GMLS::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 polynomial coefficients. But now that is inside of generatePolynomialCoefficients because ! it must be to handle number_of_batches>1. Effectively, this just calls generatePolynomialCoefficients. !

Parameters
number_of_batches[in] - how many batches to break up the total workload into (for storage) !
keep_coefficients[in] - whether to store (P^T W P)^-1 * P^T * W

Definition at line 446 of file Compadre_GMLS.cpp.

void Compadre::GMLS::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 batch of GMLS problems to be solved for. Provides alpha values ! that can later be contracted against data or degrees of freedom to form a ! global linear system. !

Parameters
number_of_batches[in] - how many batches to break up the total workload into (for storage) !
keep_coefficients[in] - whether to store (P^T W P)^-1 * P^T * W

Definition at line 11 of file Compadre_GMLS.cpp.

KOKKOS_INLINE_FUNCTION int Compadre::GMLS::getAdditionalEvaluationIndex ( const int  target_index,
const int  additional_list_num 
) const
inlineprotected

(OPTIONAL) Mapping from [0,number of additional evaluation sites for a target] to the row that contains the coordinates for that evaluation

Definition at line 462 of file Compadre_GMLS.hpp.

double Compadre::GMLS::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
inline

Underlying function all interface helper functions call to retrieve alpha values.

Definition at line 1173 of file Compadre_GMLS.hpp.

double Compadre::GMLS::getAlpha0TensorTo0Tensor ( TargetOperation  lro,
const int  target_index,
const int  neighbor_index,
const int  additional_evaluation_site = 0 
) const
inline

Helper function for getting alphas for scalar reconstruction from scalar data.

Definition at line 1092 of file Compadre_GMLS.hpp.

double Compadre::GMLS::getAlpha0TensorTo1Tensor ( TargetOperation  lro,
const int  target_index,
const int  output_component,
const int  neighbor_index,
const int  additional_evaluation_site = 0 
) const
inline

Helper function for getting alphas for vector reconstruction from scalar data.

Definition at line 1098 of file Compadre_GMLS.hpp.

double Compadre::GMLS::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
inline

Helper function for getting alphas for matrix reconstruction from scalar data.

Definition at line 1104 of file Compadre_GMLS.hpp.

double Compadre::GMLS::getAlpha1TensorTo0Tensor ( TargetOperation  lro,
const int  target_index,
const int  neighbor_index,
const int  input_component,
const int  additional_evaluation_site = 0 
) const
inline

Helper function for getting alphas for scalar reconstruction from vector data.

Definition at line 1109 of file Compadre_GMLS.hpp.

double Compadre::GMLS::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
inline

Helper function for getting alphas for vector reconstruction from vector data.

Definition at line 1115 of file Compadre_GMLS.hpp.

double Compadre::GMLS::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
inline

Helper function for getting alphas for matrix reconstruction from vector data.

Definition at line 1121 of file Compadre_GMLS.hpp.

double Compadre::GMLS::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
inline

Helper function for getting alphas for scalar reconstruction from matrix data.

Definition at line 1127 of file Compadre_GMLS.hpp.

double Compadre::GMLS::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
inline

Helper function for getting alphas for vector reconstruction from matrix data.

Definition at line 1132 of file Compadre_GMLS.hpp.

double Compadre::GMLS::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
inline

Helper function for getting alphas for matrix reconstruction from matrix data.

Definition at line 1137 of file Compadre_GMLS.hpp.

int Compadre::GMLS::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
inline

Retrieves the offset for an operator based on input and output component, generic to row (but still multiplied by the number of neighbors for each row and then needs a neighbor number added to this returned value to be meaningful)

Definition at line 1049 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION global_index_type Compadre::GMLS::getAlphaIndexDevice ( const int  target_index,
const int  alpha_column_offset 
) const
inline

Gives index into alphas given two axes, which when incremented by the neighbor number transforms access into alphas from a rank 1 view into a rank 3 view.

Definition at line 1144 of file Compadre_GMLS.hpp.

global_index_type Compadre::GMLS::getAlphaIndexHost ( const int  target_index,
const int  alpha_column_offset 
) const
inline

Gives index into alphas given two axes, which when incremented by the neighbor number transforms access into alphas from a rank 1 view into a rank 3 view.

Definition at line 1159 of file Compadre_GMLS.hpp.

decltype(_alphas) Compadre::GMLS::getAlphas ( ) const
inline

Get a view (device) of all alphas.

Definition at line 1067 of file Compadre_GMLS.hpp.

ConstraintType Compadre::GMLS::getConstraintType ( )
inline

Get constraint type.

Definition at line 968 of file Compadre_GMLS.hpp.

SamplingFunctional Compadre::GMLS::getDataSamplingFunctional ( ) const
inline

Get the data sampling functional specified at instantiation (often the same as the polynomial sampling functional)

Definition at line 1086 of file Compadre_GMLS.hpp.

DenseSolverType Compadre::GMLS::getDenseSolverType ( )
inline

Get dense solver type.

Definition at line 962 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getDimensionOfQuadraturePoints ( ) const
inline

Dimensions of quadrature points.

Definition at line 989 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getDimensions ( ) const
inline

Dimension of the GMLS problem, set only at class instantiation.

Definition at line 953 of file Compadre_GMLS.hpp.

decltype(_RHS) Compadre::GMLS::getFullPolynomialCoefficientsBasis ( ) const
inline

Get a view (device) of all polynomial coefficients basis.

Definition at line 1070 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getGlobalDimensions ( ) const
inline

Dimension of the GMLS problem's point data (spatial description of points in ambient space), set only at class instantiation.

Definition at line 956 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getInputDimensionOfOperation ( TargetOperation  lro) const
inline

Dimensions ^ input rank for target operation (always in local chart if on a manifold, never ambient space)

Definition at line 1234 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getInputDimensionOfSampling ( SamplingFunctional  sro) const
inline

Dimensions ^ output rank for sampling operation (always in ambient space, never local chart on a manifold)

Definition at line 1247 of file Compadre_GMLS.hpp.

static int Compadre::GMLS::getInputRankOfSampling ( SamplingFunctional  sro)
inlinestatic

Input rank for sampling operation.

Definition at line 907 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getLocalDimensions ( ) const
inline

Local dimension of the GMLS problem (less than global dimension if on a manifold), set only at class instantiation.

Definition at line 959 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getManifoldWeightingPower ( ) const
inline

Power for weighting kernel for curvature.

Definition at line 980 of file Compadre_GMLS.hpp.

WeightingFunctionType Compadre::GMLS::getManifoldWeightingType ( ) const
inline

Type for weighting kernel for curvature.

Definition at line 974 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION int Compadre::GMLS::getMaxEvaluationSitesPerTarget ( ) const
inlineprotected

Returns the maximum number of evaluation sites over all target sites (target sites are included in total)

Definition at line 447 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION int Compadre::GMLS::getMaxNNeighbors ( ) const
inlineprotected

Returns the maximum neighbor lists size over all target sites.

Definition at line 441 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION double Compadre::GMLS::getNeighborCoordinate ( const int  target_index,
const int  neighbor_list_num,
const int  dim,
const scratch_matrix_right_type V = NULL 
) const
inlineprotected

Returns one component of the neighbor coordinate for a particular target.

Whether global or local coordinates depends upon V being specified

Definition at line 503 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION int Compadre::GMLS::getNeighborIndex ( const int  target_index,
const int  neighbor_list_num 
) const
inlineprotected

Mapping from [0,number of neighbors for a target] to the row that contains the source coordinates for that neighbor.

Definition at line 435 of file Compadre_GMLS.hpp.

decltype(_neighbor_lists) * Compadre::GMLS::getNeighborLists ( )
inline

Get neighbor list accessor.

Definition at line 995 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION int Compadre::GMLS::getNEvaluationSitesPerTarget ( const int  target_index) const
inlineprotected

(OPTIONAL) Returns number of additional evaluation sites for a particular target

Definition at line 454 of file Compadre_GMLS.hpp.

static KOKKOS_INLINE_FUNCTION int Compadre::GMLS::getNN ( const int  m,
const int  dimension = 3,
const ReconstructionSpace  r_space = ReconstructionSpace::ScalarTaylorPolynomial 
)
inlinestatic

Returns number of neighbors needed for unisolvency for a given basis order and dimension.

Definition at line 831 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION int Compadre::GMLS::getNNeighbors ( const int  target_index) const
inlineprotected

Returns number of neighbors for a particular target.

Definition at line 428 of file Compadre_GMLS.hpp.

static KOKKOS_INLINE_FUNCTION int Compadre::GMLS::getNP ( const int  m,
const int  dimension = 3,
const ReconstructionSpace  r_space = ReconstructionSpace::ScalarTaylorPolynomial 
)
inlinestatic

Returns size of the basis for a given polynomial order and dimension General to dimension 1..3 and polynomial order m The divfree options will return the divergence-free basis if true.

Examples:
GMLS Tutorial, and Manifold GMLS Tutorial.

Definition at line 821 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getNumberOfQuadraturePoints ( ) const
inline

Number of quadrature points.

Definition at line 983 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getOrderOfQuadraturePoints ( ) const
inline

Order of quadrature points.

Definition at line 986 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getOutputDimensionOfOperation ( TargetOperation  lro,
bool  ambient = false 
) const
inline

Dimensions ^ output rank for target operation.

Definition at line 1226 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getOutputDimensionOfSampling ( SamplingFunctional  sro) const
inline

Dimensions ^ output rank for sampling operation (always in local chart if on a manifold, never ambient space)

Definition at line 1241 of file Compadre_GMLS.hpp.

static int Compadre::GMLS::getOutputRankOfSampling ( SamplingFunctional  sro)
inlinestatic

Output rank for sampling operation.

Definition at line 902 of file Compadre_GMLS.hpp.

host_managed_local_index_type Compadre::GMLS::getPolynomialCoefficientsDomainRangeSize ( ) const
inline

Returns (size of the basis used in instance's polynomial reconstruction) x (data input dimension)

Definition at line 921 of file Compadre_GMLS.hpp.

host_managed_local_index_type Compadre::GMLS::getPolynomialCoefficientsMemorySize ( ) const
inline

Returns 2D array size in memory on which coefficients are stored.

Definition at line 935 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getPolynomialCoefficientsSize ( ) const
inline

Returns size of the basis used in instance's polynomial reconstruction.

Definition at line 929 of file Compadre_GMLS.hpp.

SamplingFunctional Compadre::GMLS::getPolynomialSamplingFunctional ( ) const
inline

Get the polynomial sampling functional specified at instantiation.

Definition at line 1083 of file Compadre_GMLS.hpp.

double Compadre::GMLS::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
inline

Returns a stencil to transform data from its existing state into the input expected for some sampling functionals.

Definition at line 1203 of file Compadre_GMLS.hpp.

decltype(_prestencil_weights) Compadre::GMLS::getPrestencilWeights ( ) const
inline

Get a view (device) of all rank 2 preprocessing tensors This is a rank 5 tensor that is able to provide data transformation into a form that GMLS is able to operate on.

The ranks are as follows:

1 - Either size 2 if it operates on the target site and neighbor site (staggered schemes) or 1 if it operates only on the neighbor sites (almost every scheme)

2 - If the data transform varies with each target site (but could be the same for each neighbor of that target site), then this is the number of target sites

3 - If the data transform varies with each neighbor of each target site, then this is the number of neighbors for each respective target (max number of neighbors for all target sites is its uniform size)

4 - Data transform resulting in rank 1 data for the GMLS operator will have size _local_dimensions, otherwise 1

5 - Data transform taking in rank 1 data will have size _global_dimensions, otherwise 1

Definition at line 1042 of file Compadre_GMLS.hpp.

ProblemType Compadre::GMLS::getProblemType ( )
inline

Get problem type.

Definition at line 965 of file Compadre_GMLS.hpp.

std::string Compadre::GMLS::getQuadratureType ( ) const
inline

Type of quadrature points.

Definition at line 992 of file Compadre_GMLS.hpp.

ReconstructionSpace Compadre::GMLS::getReconstructionSpace ( ) const
inline

Get the reconstruction space specified at instantiation.

Definition at line 1089 of file Compadre_GMLS.hpp.

double Compadre::GMLS::getReferenceNormalDirection ( const int  target_index,
const int  component 
) const
inline

Get component of tangent or normal directions for manifold problems.

Definition at line 1013 of file Compadre_GMLS.hpp.

decltype(_ref_N) Compadre::GMLS::getReferenceNormalDirections ( ) const
inline

Get a view (device) of all reference outward normal directions.

Definition at line 1001 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION XYZ Compadre::GMLS::getRelativeCoord ( const int  target_index,
const int  neighbor_list_num,
const int  dimension,
const scratch_matrix_right_type V = NULL 
) const
inlineprotected

Returns the relative coordinate as a vector between the target site and the neighbor site.

Whether global or local coordinates depends upon V being specified

Definition at line 516 of file Compadre_GMLS.hpp.

static int Compadre::GMLS::getSamplingOutputIndex ( const SamplingFunctional  sf,
const int  output_component_axis_1,
const int  output_component_axis_2 
)
inlinestatic

Helper function for finding alpha coefficients.

Definition at line 896 of file Compadre_GMLS.hpp.

double Compadre::GMLS::getTangentBundle ( const int  target_index,
const int  direction,
const int  component 
) const
inline

Get component of tangent or normal directions for manifold problems.

Definition at line 1004 of file Compadre_GMLS.hpp.

decltype(_T) Compadre::GMLS::getTangentDirections ( ) const
inline

Get a view (device) of all tangent direction bundles.

Definition at line 998 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION double Compadre::GMLS::getTargetAuxiliaryCoordinate ( const int  target_index,
const int  additional_list_num,
const int  dim,
const scratch_matrix_right_type V = NULL 
) const
inlineprotected

(OPTIONAL) Returns one component of the additional evaluation coordinates.

Whether global or local coordinates depends upon V being specified

Definition at line 487 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION double Compadre::GMLS::getTargetCoordinate ( const int  target_index,
const int  dim,
const scratch_matrix_right_type V = NULL 
) const
inlineprotected

Returns one component of the target coordinate for a particular target.

Whether global or local coordinates depends upon V being specified

Definition at line 471 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION int Compadre::GMLS::getTargetOffsetIndexDevice ( const int  lro_num,
const int  input_component,
const int  output_component,
const int  additional_evaluation_local_index = 0 
) const
inlineprotected

Handles offset from operation input/output + extra evaluation sites.

Definition at line 562 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getTargetOffsetIndexHost ( const int  lro_num,
const int  input_component,
const int  output_component,
const int  additional_evaluation_local_index = 0 
) const
inlineprotected

Handles offset from operation input/output + extra evaluation sites.

Definition at line 553 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getTargetOperationLocalIndex ( TargetOperation  lro) const
inline

Get the local index (internal) to GMLS for a particular TargetOperation Every TargetOperation has a global index which can be readily found in Compadre::TargetOperation but this function returns the index used inside of the GMLS class.

Definition at line 1024 of file Compadre_GMLS.hpp.

static int Compadre::GMLS::getTargetOutputIndex ( const int  operation_num,
const int  output_component_axis_1,
const int  output_component_axis_2,
const int  dimensions 
)
inlinestatic

Helper function for finding alpha coefficients.

Definition at line 890 of file Compadre_GMLS.hpp.

int Compadre::GMLS::getWeightingPower ( ) const
inline

Power for weighting kernel for GMLS problem.

Definition at line 977 of file Compadre_GMLS.hpp.

WeightingFunctionType Compadre::GMLS::getWeightingType ( ) const
inline

Type for weighting kernel for GMLS problem.

Definition at line 971 of file Compadre_GMLS.hpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::operator() ( const AssembleStandardPsqrtW ,
const member_type teamMember 
) const

Functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity.

Definition at line 454 of file Compadre_GMLS.cpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::operator() ( const ApplyStandardTargets ,
const member_type teamMember 
) const

Functor to evaluate targets, apply target evaluation to polynomial coefficients to store in _alphas.

Definition at line 533 of file Compadre_GMLS.cpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::operator() ( const ComputeCoarseTangentPlane ,
const member_type teamMember 
) const

Functor to create a coarse tangent approximation from a given neighborhood of points.

Definition at line 586 of file Compadre_GMLS.cpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::operator() ( const AssembleCurvaturePsqrtW ,
const member_type teamMember 
) const

Functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity for curvature.

Definition at line 640 of file Compadre_GMLS.cpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::operator() ( const GetAccurateTangentDirections ,
const member_type teamMember 
) const

Functor to evaluate curvature targets and construct accurate tangent direction approximation for manifolds.

Definition at line 715 of file Compadre_GMLS.cpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::operator() ( const FixTangentDirectionOrdering ,
const member_type teamMember 
) const

Functor to determine if tangent directions need reordered, and to reorder them if needed We require that the normal is consistent with a right hand rule on the tangent vectors.

Definition at line 868 of file Compadre_GMLS.cpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::operator() ( const ApplyCurvatureTargets ,
const member_type teamMember 
) const

Functor to evaluate curvature targets and apply to coefficients of curvature reconstruction.

Definition at line 913 of file Compadre_GMLS.cpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::operator() ( const AssembleManifoldPsqrtW ,
const member_type teamMember 
) const

Functor to assemble the P*sqrt(weights) matrix and construct sqrt(weights)*Identity.

Definition at line 1057 of file Compadre_GMLS.cpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::operator() ( const ApplyManifoldTargets ,
const member_type teamMember 
) const

Functor to evaluate targets, apply target evaluation to polynomial coefficients to store in _alphas.

Definition at line 1128 of file Compadre_GMLS.cpp.

KOKKOS_INLINE_FUNCTION void Compadre::GMLS::operator() ( const ComputePrestencilWeights ,
const member_type teamMember 
) const

Functor to calculate prestencil weights to apply to data to transform into a format expected by a GMLS stencil.

Definition at line 1191 of file Compadre_GMLS.cpp.

static ConstraintType Compadre::GMLS::parseConstraintType ( const std::string &  constraint_type)
inlinestaticprotected

Parses a string to determine constraint type.

Definition at line 603 of file Compadre_GMLS.hpp.

static ProblemType Compadre::GMLS::parseProblemType ( const std::string &  problem_type)
inlinestaticprotected

Parses a string to determine problem type.

Definition at line 590 of file Compadre_GMLS.hpp.

static DenseSolverType Compadre::GMLS::parseSolverType ( const std::string &  dense_solver_type)
inlinestaticprotected

Parses a string to determine solver type.

Definition at line 577 of file Compadre_GMLS.hpp.

void Compadre::GMLS::resetCoefficientData ( )
inline

Definition at line 1282 of file Compadre_GMLS.hpp.

template<typename view_type_1 , typename view_type_2 >
void Compadre::GMLS::setAdditionalEvaluationSitesData ( view_type_1  additional_evaluation_indices,
view_type_2  additional_evaluation_coordinates 
)
inline

(OPTIONAL) Sets additional evaluation sites for each target site

Definition at line 1316 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setAuxiliaryEvaluationCoordinates ( view_type  evaluation_coordinates)
inline

(OPTIONAL) Sets additional points for evaluation of target operation on polynomial reconstruction.

If this is never called, then the target sites are the only locations where the target operations will be evaluated and applied to polynomial reconstructions.

Definition at line 1594 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setAuxiliaryEvaluationCoordinates ( decltype(_additional_evaluation_coordinates evaluation_coordinates)
inline

(OPTIONAL) Sets additional points for evaluation of target operation on polynomial reconstruction.

If this is never called, then the target sites are the only locations where the target operations will be evaluated and applied to polynomial reconstructions. (device)

Definition at line 1623 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setAuxiliaryEvaluationIndicesLists ( view_type  indices_lists)
inline

(OPTIONAL) Sets the additional target evaluation coordinate indices list information.

Should be # targets x maximum number of indices evaluation indices for any target + 1. first entry in every row should be the number of indices for the corresponding target.

Definition at line 1632 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setAuxiliaryEvaluationIndicesLists ( decltype(_additional_evaluation_indices indices_lists)
inline

(OPTIONAL) Sets the additional target evaluation coordinate indices list information.

Should be # targets x maximum number of indices evaluation indices for any target + 1. first entry in every row should be the number of indices for the corresponding target.

Definition at line 1674 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setCurvaturePolynomialOrder ( const int  manifold_poly_order)
inline

Sets basis order to be used when reoncstructing curvature.

Definition at line 1751 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setCurvatureWeightingPower ( int  wp)
inline

Power for weighting kernel for curvature.

Definition at line 1763 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setCurvatureWeightingType ( const std::string &  wt)
inline

Type for weighting kernel for curvature.

Definition at line 1721 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setCurvatureWeightingType ( const WeightingFunctionType  wt)
inline

Type for weighting kernel for curvature.

Definition at line 1738 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setDimensionOfQuadraturePoints ( int  dim)
inline

Dimensions of quadrature points to use.

Definition at line 1775 of file Compadre_GMLS.hpp.

template<typename view_type >
std::enable_if<view_type::rank==1&&std::is_same<decltype(_neighbor_lists)::internal_view_type,view_type>::value==1, void>::type Compadre::GMLS::setNeighborLists ( view_type  neighbor_lists,
view_type  number_of_neighbors_list 
)
inline

Sets neighbor list information from compressed row neighborhood lists data (if same view_type).

Definition at line 1326 of file Compadre_GMLS.hpp.

template<typename view_type >
std::enable_if<view_type::rank==1&&std::is_same<decltype(_neighbor_lists)::internal_view_type,view_type>::value==0, void>::type Compadre::GMLS::setNeighborLists ( view_type  neighbor_lists,
view_type  number_of_neighbors_list 
)
inline

Sets neighbor list information from compressed row neighborhood lists data (if different view_type).

Definition at line 1342 of file Compadre_GMLS.hpp.

template<typename view_type >
std::enable_if<view_type::rank==2, void>::type Compadre::GMLS::setNeighborLists ( view_type  neighbor_lists)
inline

Sets neighbor list information.

Should be # targets x maximum number of neighbors for any target + 1. first entry in ever row should be the number of neighbors for the corresponding target.

Definition at line 1364 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setOrderOfQuadraturePoints ( int  order)
inline

Number quadrature points to use.

Definition at line 1769 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setPolynomialOrder ( const int  poly_order)
inline

Sets basis order to be used when reoncstructing any function.

Definition at line 1744 of file Compadre_GMLS.hpp.

template<typename view_type_1 , typename view_type_2 , typename view_type_3 , typename view_type_4 >
void Compadre::GMLS::setProblemData ( view_type_1  neighbor_lists,
view_type_2  source_coordinates,
view_type_3  target_coordinates,
view_type_4  epsilons 
)
inline

Sets basic problem data (neighbor lists, source coordinates, and target coordinates)

Definition at line 1289 of file Compadre_GMLS.hpp.

template<typename view_type_1 , typename view_type_2 , typename view_type_3 , typename view_type_4 >
void Compadre::GMLS::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 
)
inline

Sets basic problem data (neighbor lists data, number of neighbors list, source coordinates, and target coordinates)

Definition at line 1302 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setQuadratureType ( std::string  quadrature_type)
inline

Type of quadrature points.

Definition at line 1781 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setReferenceOutwardNormalDirection ( view_type  outward_normal_directions,
bool  use_to_orient_surface = true 
)
inline

(OPTIONAL) Sets outward normal direction.

For manifolds this may be used for orienting surface. It is also accessible for sampling operators that require a normal direction.

Definition at line 1516 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setSourceExtraData ( view_type  extra_data)
inline

(OPTIONAL) Sets extra data to be used by sampling functionals in certain instances.

Definition at line 1544 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setSourceExtraData ( decltype(_source_extra_data extra_data)
inline

(OPTIONAL) Sets extra data to be used by sampling functionals in certain instances.

(device)

Definition at line 1559 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setSourceSites ( view_type  source_coordinates)
inline

Sets source coordinate information.

Rows of this 2D-array should correspond to neighbor IDs contained in the entries of the neighbor lists 2D array.

Definition at line 1380 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setSourceSites ( decltype(_source_coordinates source_coordinates)
inline

Sets source coordinate information.

Rows of this 2D-array should correspond to neighbor IDs contained in the entries of the neighbor lists 2D array.

Definition at line 1408 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setTangentBundle ( view_type  tangent_directions)
inline

(OPTIONAL) Sets orthonormal tangent directions for reconstruction on a manifold.

The first rank of this 2D array corresponds to the target indices, i.e., rows of the neighbor lists 2D array. The second rank is the ordinal of the tangent direction (spatial dimensions-1 are tangent, last one is normal), and the third rank is indices into the spatial dimension.

Definition at line 1482 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setTargetExtraData ( view_type  extra_data)
inline

(OPTIONAL) Sets extra data to be used by target operations in certain instances.

Definition at line 1568 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setTargetExtraData ( decltype(_target_extra_data extra_data)
inline

(OPTIONAL) Sets extra data to be used by target operations in certain instances.

(device)

Definition at line 1583 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setTargetSites ( view_type  target_coordinates)
inline

Sets target coordinate information. Rows of this 2D-array should correspond to rows of the neighbor lists.

Definition at line 1416 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setTargetSites ( decltype(_target_coordinates target_coordinates)
inline

Sets target coordinate information. Rows of this 2D-array should correspond to rows of the neighbor lists.

Definition at line 1445 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setWeightingPower ( int  wp)
inline

Power for weighting kernel for GMLS problem.

Definition at line 1757 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setWeightingType ( const std::string &  wt)
inline

Type for weighting kernel for GMLS problem.

Definition at line 1698 of file Compadre_GMLS.hpp.

void Compadre::GMLS::setWeightingType ( const WeightingFunctionType  wt)
inline

Type for weighting kernel for GMLS problem.

Definition at line 1715 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setWindowSizes ( view_type  epsilons)
inline

Sets window sizes, also called the support of the kernel.

Definition at line 1456 of file Compadre_GMLS.hpp.

template<typename view_type >
void Compadre::GMLS::setWindowSizes ( decltype(_epsilons epsilons)
inline

Sets window sizes, also called the support of the kernel (device)

Definition at line 1470 of file Compadre_GMLS.hpp.

static KOKKOS_INLINE_FUNCTION double Compadre::GMLS::Wab ( const double  r,
const double  h,
const WeightingFunctionType weighting_type,
const int  power 
)
inlinestatic

Evaluates the weighting kernel.

Parameters
r[in] - Euclidean distance of relative vector. Euclidean distance of (target - neighbor) in some basis.
h[in] - window size. Kernel is guaranteed to take on a value of zero if it exceeds h.
weighting_type[in] - weighting type to be evaluated as the kernel. e,g. power, Gaussian, etc..
power[in] - power parameter to be given to the kernel.

Definition at line 855 of file Compadre_GMLS.hpp.

Member Data Documentation

int Compadre::GMLS::_added_alpha_size
protected

additional alpha coefficients due to constraints

Definition at line 266 of file Compadre_GMLS.hpp.

Kokkos::View<double**, layout_right> Compadre::GMLS::_additional_evaluation_coordinates
protected

(OPTIONAL) user provided additional coordinates for target operation evaluation (device)

Definition at line 106 of file Compadre_GMLS.hpp.

Kokkos::View<int**, layout_right> Compadre::GMLS::_additional_evaluation_indices
protected

(OPTIONAL) contains indices of entries in the _additional_evaluation_coordinates view (device)

Definition at line 109 of file Compadre_GMLS.hpp.

Kokkos::View<double*, layout_right> Compadre::GMLS::_alphas
protected

generated alpha coefficients (device)

Definition at line 92 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_basis_multiplier
protected

dimension of the reconstructed function e.g.

reconstruction of vector on a 2D manifold in 3D would have _basis_multiplier of 2

Definition at line 183 of file Compadre_GMLS.hpp.

ConstraintType Compadre::GMLS::_constraint_type
protected

constraint type for GMLS problem

Definition at line 149 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_curvature_poly_order
protected

order of basis for curvature reconstruction

Definition at line 122 of file Compadre_GMLS.hpp.

Kokkos::View<TargetOperation*> Compadre::GMLS::_curvature_support_operations
protected

vector containing target functionals to be applied for curvature

Definition at line 159 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_curvature_weighting_power
protected

power to be used for weighting kernel for curvature

Definition at line 179 of file Compadre_GMLS.hpp.

WeightingFunctionType Compadre::GMLS::_curvature_weighting_type
protected

weighting kernel type for curvature problem

Definition at line 173 of file Compadre_GMLS.hpp.

const SamplingFunctional Compadre::GMLS::_data_sampling_functional
protected

generally the same as _polynomial_sampling_functional, but can differ if specified at GMLS class instantiation

Definition at line 156 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_data_sampling_multiplier
protected

effective dimension of the data sampling functional e.g.

in 3D, a scalar will be 1, a vector will be 3, and a vector of reused scalars will be 3

Definition at line 192 of file Compadre_GMLS.hpp.

DenseSolverType Compadre::GMLS::_dense_solver_type
protected

solver type for GMLS problem - can be QR, SVD or LU

Definition at line 143 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_dimension_of_quadrature_points
protected

dimension of quadrature rule

Definition at line 275 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_dimensions
protected

dimension of the problem, set at class instantiation only

Definition at line 134 of file Compadre_GMLS.hpp.

bool Compadre::GMLS::_entire_batch_computed_at_once
protected

whether entire calculation was computed at once the alternative is that it was broken up over many smaller groups, in which case this is false, and so the _RHS matrix can not be stored or requested

Definition at line 211 of file Compadre_GMLS.hpp.

Kokkos::View<double*> Compadre::GMLS::_epsilons
protected

h supports determined through neighbor search (device)

Definition at line 86 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_global_dimensions
protected

spatial dimension of the points, set at class instantiation only

Definition at line 128 of file Compadre_GMLS.hpp.

Kokkos::View<int**, layout_right>::HostMirror Compadre::GMLS::_host_additional_evaluation_indices
protected

(OPTIONAL) contains indices of entries in the _additional_evaluation_coordinates view (host)

Definition at line 112 of file Compadre_GMLS.hpp.

Kokkos::View<const double*, layout_right>::HostMirror Compadre::GMLS::_host_alphas
protected

generated alpha coefficients (host)

Definition at line 95 of file Compadre_GMLS.hpp.

Kokkos::View<double*>::HostMirror Compadre::GMLS::_host_epsilons
protected

h supports determined through neighbor search (host)

Definition at line 89 of file Compadre_GMLS.hpp.

Kokkos::View<int*>::HostMirror Compadre::GMLS::_host_lro_input_tensor_rank
protected

tensor rank of sampling functional (host)

Definition at line 260 of file Compadre_GMLS.hpp.

Kokkos::View<int*>::HostMirror Compadre::GMLS::_host_lro_input_tile_size
protected

dimensions ^ rank of tensor of output for each sampling functional (host)

Definition at line 248 of file Compadre_GMLS.hpp.

Kokkos::View<int*>::HostMirror Compadre::GMLS::_host_lro_output_tensor_rank
protected

tensor rank of target functional (host)

Definition at line 254 of file Compadre_GMLS.hpp.

Kokkos::View<int*>::HostMirror Compadre::GMLS::_host_lro_output_tile_size
protected

dimensions ^ rank of tensor of output for each target functional (host)

Definition at line 242 of file Compadre_GMLS.hpp.

Kokkos::View<int*>::HostMirror Compadre::GMLS::_host_lro_total_offsets
protected

index for where this operation begins the for _alpha coefficients (host)

Definition at line 236 of file Compadre_GMLS.hpp.

Kokkos::View<int*, host_memory_space> Compadre::GMLS::_host_number_of_neighbors_list
protected

convenient copy on host of number of neighbors

Definition at line 77 of file Compadre_GMLS.hpp.

Kokkos::View<TargetOperation*>::HostMirror Compadre::GMLS::_host_operations
protected

vector containing target functionals to be applied for reconstruction problem (host)

Definition at line 165 of file Compadre_GMLS.hpp.

Kokkos::View<const double*****, layout_right>::HostMirror Compadre::GMLS::_host_prestencil_weights
protected

generated weights for nontraditional samples required to transform data into expected sampling functional form (host)

Definition at line 103 of file Compadre_GMLS.hpp.

Kokkos::View<double*>::HostMirror Compadre::GMLS::_host_ref_N
protected

reference outward normal vectors information (host)

Definition at line 56 of file Compadre_GMLS.hpp.

Kokkos::View<double*>::HostMirror Compadre::GMLS::_host_T
protected

tangent vectors information (host)

Definition at line 53 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_initial_index_for_batch
protected

initial index for current batch

Definition at line 217 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_local_dimensions
protected

dimension of the problem, set at class instantiation only. For manifolds, generally _global_dimensions-1

Definition at line 131 of file Compadre_GMLS.hpp.

std::vector<TargetOperation> Compadre::GMLS::_lro
protected

vector of user requested target operations

Definition at line 226 of file Compadre_GMLS.hpp.

Kokkos::View<int*> Compadre::GMLS::_lro_input_tensor_rank
protected

tensor rank of sampling functional (device)

Definition at line 257 of file Compadre_GMLS.hpp.

Kokkos::View<int*> Compadre::GMLS::_lro_input_tile_size
protected

dimensions ^ rank of tensor of output for each sampling functional (device)

Definition at line 245 of file Compadre_GMLS.hpp.

std::vector<int> Compadre::GMLS::_lro_lookup
protected

vector containing a mapping from a target functionals enum value to the its place in the list of target functionals to be applied

Definition at line 230 of file Compadre_GMLS.hpp.

Kokkos::View<int*> Compadre::GMLS::_lro_output_tensor_rank
protected

tensor rank of target functional (device)

Definition at line 251 of file Compadre_GMLS.hpp.

Kokkos::View<int*> Compadre::GMLS::_lro_output_tile_size
protected

dimensions ^ rank of tensor of output for each target functional (device)

Definition at line 239 of file Compadre_GMLS.hpp.

Kokkos::View<int*> Compadre::GMLS::_lro_total_offsets
protected

index for where this operation begins the for _alpha coefficients (device)

Definition at line 233 of file Compadre_GMLS.hpp.

Kokkos::View<double*> Compadre::GMLS::_manifold_curvature_coefficients
protected

curvature polynomial coefficients for all problems

Definition at line 62 of file Compadre_GMLS.hpp.

Kokkos::View<double*> Compadre::GMLS::_manifold_curvature_gradient
protected

_dimension-1 gradient values for curvature for all problems

Definition at line 65 of file Compadre_GMLS.hpp.

Kokkos::View<double*> Compadre::GMLS::_manifold_metric_tensor_inverse
protected

metric tensor inverse for all problems

Definition at line 59 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_max_evaluation_sites_per_target
protected

maximum number of evaluation sites for each target (includes target site)

Definition at line 223 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_max_num_neighbors
protected

maximum number of neighbors over all target sites

Definition at line 220 of file Compadre_GMLS.hpp.

NeighborLists<Kokkos::View<int*> > Compadre::GMLS::_neighbor_lists
protected

Accessor to get neighbor list data, offset data, and number of neighbors per target.

Definition at line 74 of file Compadre_GMLS.hpp.

bool Compadre::GMLS::_nontrivial_nullspace
protected

whether or not operator to be inverted for GMLS problem has a nontrivial nullspace (requiring SVD)

Definition at line 195 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_NP
protected

dimension of basis for polynomial reconstruction

Definition at line 125 of file Compadre_GMLS.hpp.

Kokkos::View<int*> Compadre::GMLS::_number_of_additional_evaluation_indices
protected

(OPTIONAL) contains the # of additional coordinate indices for each target

Definition at line 115 of file Compadre_GMLS.hpp.

Kokkos::View<TargetOperation*> Compadre::GMLS::_operations
protected

vector containing target functionals to be applied for reconstruction problem (device)

Definition at line 162 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_order_of_quadrature_points
protected

order of exact polynomial integration for quadrature rule

Definition at line 272 of file Compadre_GMLS.hpp.

bool Compadre::GMLS::_orthonormal_tangent_space_provided
protected

whether or not the orthonormal tangent directions were provided by the user.

If they are not, then for the case of calculations on manifolds, a GMLS approximation of the tangent space will be made and stored for use.

Definition at line 200 of file Compadre_GMLS.hpp.

Kokkos::View<double*> Compadre::GMLS::_P
protected

P*sqrt(w) matrix for all problems.

Definition at line 38 of file Compadre_GMLS.hpp.

ParallelManager Compadre::GMLS::_pm
protected

determines scratch level spaces and is used to call kernels

Definition at line 269 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_poly_order
protected

order of basis for polynomial reconstruction

Definition at line 119 of file Compadre_GMLS.hpp.

const SamplingFunctional Compadre::GMLS::_polynomial_sampling_functional
protected

polynomial sampling functional used to construct P matrix, set at GMLS class instantiation

Definition at line 152 of file Compadre_GMLS.hpp.

Kokkos::View<double*****, layout_right> Compadre::GMLS::_prestencil_weights
protected

generated weights for nontraditional samples required to transform data into expected sampling functional form (device).

Definition at line 99 of file Compadre_GMLS.hpp.

ProblemType Compadre::GMLS::_problem_type
protected

problem type for GMLS problem, can also be set to STANDARD for normal or MANIFOLD for manifold problems

Definition at line 146 of file Compadre_GMLS.hpp.

Quadrature Compadre::GMLS::_qm
protected

manages and calculates quadrature

Definition at line 281 of file Compadre_GMLS.hpp.

std::string Compadre::GMLS::_quadrature_type
protected

quadrature rule type

Definition at line 278 of file Compadre_GMLS.hpp.

pool_type Compadre::GMLS::_random_number_pool
protected

Definition at line 29 of file Compadre_GMLS.hpp.

ReconstructionSpace Compadre::GMLS::_reconstruction_space
protected

reconstruction space for GMLS problems, set at GMLS class instantiation

Definition at line 137 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_reconstruction_space_rank
protected

actual rank of reconstruction basis

Definition at line 140 of file Compadre_GMLS.hpp.

Kokkos::View<double*> Compadre::GMLS::_ref_N
protected

Rank 2 tensor for high order approximation of tangent vectors for all problems.

First rank is for the target index, the second is for the spatial dimension (_dimensions)

Definition at line 50 of file Compadre_GMLS.hpp.

bool Compadre::GMLS::_reference_outward_normal_direction_provided
protected

whether or not the reference outward normal directions were provided by the user.

Definition at line 203 of file Compadre_GMLS.hpp.

Kokkos::View<double*> Compadre::GMLS::_RHS
protected

sqrt(w)*Identity matrix for all problems, later holds polynomial coefficients for all problems

Definition at line 41 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_sampling_multiplier
protected

actual dimension of the sampling functional e.g.

reconstruction of vector on a 2D manifold in 3D would have _basis_multiplier of 2 e.g. in 3D, a scalar will be 1, a vector will be 3, and a vector of reused scalars will be 1

Definition at line 188 of file Compadre_GMLS.hpp.

Kokkos::View<double**, layout_right> Compadre::GMLS::_source_coordinates
protected

all coordinates for the source for which _neighbor_lists refers (device)

Definition at line 80 of file Compadre_GMLS.hpp.

Kokkos::View<double**, layout_right> Compadre::GMLS::_source_extra_data
protected

Extra data available to basis functions (optional)

Definition at line 68 of file Compadre_GMLS.hpp.

bool Compadre::GMLS::_store_PTWP_inv_PTW
protected

whether polynomial coefficients were requested to be stored (in a state not yet applied to data)

Definition at line 214 of file Compadre_GMLS.hpp.

Kokkos::View<double*> Compadre::GMLS::_T
protected

Rank 3 tensor for high order approximation of tangent vectors for all problems.

First rank is for the target index, the second is for the local direction to the manifolds 0..(_dimensions-1) are tangent, _dimensions is the normal, and the third is for the spatial dimension (_dimensions)

Definition at line 46 of file Compadre_GMLS.hpp.

Kokkos::View<double**, layout_right> Compadre::GMLS::_target_coordinates
protected

coordinates for target sites for reconstruction (device)

Definition at line 83 of file Compadre_GMLS.hpp.

Kokkos::View<double**, layout_right> Compadre::GMLS::_target_extra_data
protected

Extra data available to target operations (optional)

Definition at line 71 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_total_alpha_values
protected

used for sizing P_target_row and the _alphas view

Definition at line 263 of file Compadre_GMLS.hpp.

bool Compadre::GMLS::_use_reference_outward_normal_direction_provided_to_orient_surface
protected

whether or not to use reference outward normal directions to orient the surface in a manifold problem.

Definition at line 206 of file Compadre_GMLS.hpp.

Kokkos::View<double*> Compadre::GMLS::_w
protected

contains weights for all problems

Definition at line 35 of file Compadre_GMLS.hpp.

int Compadre::GMLS::_weighting_power
protected

power to be used for weighting kernel

Definition at line 176 of file Compadre_GMLS.hpp.

WeightingFunctionType Compadre::GMLS::_weighting_type
protected

weighting kernel type for GMLS

Definition at line 170 of file Compadre_GMLS.hpp.


The documentation for this class was generated from the following files: