PolyFEM
Loading...
Searching...
No Matches
polyfem::basis::RBFWithQuadratic Class Reference

#include <RBFWithQuadratic.hpp>

Public Member Functions

 RBFWithQuadratic (const assembler::LinearAssembler &assembler, const Eigen::MatrixXd &centers, const Eigen::MatrixXd &collocation_points, const Eigen::MatrixXd &local_basis_integral, const quadrature::Quadrature &quadr, Eigen::MatrixXd &rhs, bool with_constraints=true)
 Initialize RBF functions over a polytope element.
 
void basis (const int local_index, const Eigen::MatrixXd &uv, Eigen::MatrixXd &val) const
 Evaluates one RBF function over a list of coordinates.
 
void grad (const int local_index, const Eigen::MatrixXd &uv, Eigen::MatrixXd &val) const
 Evaluates the gradient of one RBF function over a list of coordinates.
 
void bases_values (const Eigen::MatrixXd &samples, Eigen::MatrixXd &val) const
 Batch evaluates the RBF + polynomials on a set of sample points.
 
void bases_grads (const int axis, const Eigen::MatrixXd &samples, Eigen::MatrixXd &val) const
 Batch evaluates the gradient of the RBF + polynomials on a set of sample points.
 

Static Public Member Functions

static int index_mapping (const int alpha, const int beta, const int d, const int ass_dim)
 
static void setup_monomials_vals_2d (const int star_index, const Eigen::MatrixXd &pts, assembler::ElementAssemblyValues &vals)
 
static void setup_monomials_strong_2d (const int dim, const assembler::LinearAssembler &assembler, const Eigen::MatrixXd &pts, const QuadratureVector &da, std::array< Eigen::MatrixXd, 5 > &strong)
 

Private Member Functions

bool is_volume () const
 
void compute_kernels_matrix (const Eigen::MatrixXd &samples, Eigen::MatrixXd &A) const
 
void compute_constraints_matrix_2d_old (const int num_bases, const quadrature::Quadrature &quadr, const Eigen::MatrixXd &local_basis_integral, Eigen::MatrixXd &L, Eigen::MatrixXd &t) const
 
void compute_constraints_matrix_2d (const assembler::LinearAssembler &assembler, const int num_bases, const quadrature::Quadrature &quadr, const Eigen::MatrixXd &local_basis_integral, Eigen::MatrixXd &L, Eigen::MatrixXd &t) const
 
void compute_constraints_matrix_3d (const assembler::LinearAssembler &assembler, const int num_bases, const quadrature::Quadrature &quadr, const Eigen::MatrixXd &local_basis_integral, Eigen::MatrixXd &L, Eigen::MatrixXd &t) const
 
void compute_weights (const assembler::LinearAssembler &assembler, const Eigen::MatrixXd &collocation_points, const Eigen::MatrixXd &local_basis_integral, const quadrature::Quadrature &quadr, Eigen::MatrixXd &rhs, bool with_constraints)
 

Private Attributes

Eigen::MatrixXd centers_
 
Eigen::MatrixXd weights_
 

Detailed Description

Definition at line 13 of file RBFWithQuadratic.hpp.

Constructor & Destructor Documentation

◆ RBFWithQuadratic()

RBFWithQuadratic::RBFWithQuadratic ( const assembler::LinearAssembler assembler,
const Eigen::MatrixXd &  centers,
const Eigen::MatrixXd &  collocation_points,
const Eigen::MatrixXd &  local_basis_integral,
const quadrature::Quadrature quadr,
Eigen::MatrixXd &  rhs,
bool  with_constraints = true 
)

Initialize RBF functions over a polytope element.

Parameters
[in]centers#C x dim positions of the kernels used to define functions over the polytope. The centers are placed at a small offset distance from the boundary of the element, due to the singularity at the centers
[in]collocation_points#S x dim positions of the collocation points, used to approximate the RBF functions over the boundary of the element
[in]local_basis_integral#B x dim+dim*(dim+1)/2 of the constant right-hand side for the integral constraint for each basis over the polytope
[in]quadrQuadrature points and weights inside the polytope
[in]rhs#S x #B of boundary conditions. Each column defines how the i-th basis of the mesh should evaluate on the collocation points sampled on the boundary of the polytope
[in]with_constraintsImpose integral constraints to guarantee linear reproduction for the Poisson equation

Definition at line 185 of file RBFWithQuadratic.cpp.

References compute_weights().

Here is the call graph for this function:

Member Function Documentation

◆ bases_grads()

void RBFWithQuadratic::bases_grads ( const int  axis,
const Eigen::MatrixXd &  samples,
Eigen::MatrixXd &  val 
) const

Batch evaluates the gradient of the RBF + polynomials on a set of sample points.

Parameters
[in]axisThe axis (0, 1, 2) with respect to which to compute the gradient
[in]uv#uv x dim matrix of points to evaluate
[out]val#uv x n_loc_bases of bases gradient wrt axis over the sample points

Definition at line 236 of file RBFWithQuadratic.cpp.

References centers_, is_volume(), val, weights_, and x.

Referenced by grad().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ bases_values()

void RBFWithQuadratic::bases_values ( const Eigen::MatrixXd &  samples,
Eigen::MatrixXd &  val 
) const

Batch evaluates the RBF + polynomials on a set of sample points.

Parameters
[in]uv#uv x dim matrix of points to evaluate
[out]val#uv x n_loc_bases of bases values over the sample points

Definition at line 224 of file RBFWithQuadratic.cpp.

References compute_kernels_matrix(), val, and weights_.

Referenced by basis().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ basis()

void RBFWithQuadratic::basis ( const int  local_index,
const Eigen::MatrixXd &  uv,
Eigen::MatrixXd &  val 
) const

Evaluates one RBF function over a list of coordinates.

Parameters
[in]local_indexi-th RBF function to evaluate
[in]uv#uv x dim matrix of coordinates to evaluate (in object domain)
[out]val#uv x 1 matrix of computed values

Definition at line 201 of file RBFWithQuadratic.cpp.

References bases_values(), and val.

Referenced by compute_weights().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_constraints_matrix_2d()

void RBFWithQuadratic::compute_constraints_matrix_2d ( const assembler::LinearAssembler assembler,
const int  num_bases,
const quadrature::Quadrature quadr,
const Eigen::MatrixXd &  local_basis_integral,
Eigen::MatrixXd &  L,
Eigen::MatrixXd &  t 
) const
private

Definition at line 503 of file RBFWithQuadratic.cpp.

References polyfem::assembler::LinearAssembler::assemble(), polyfem::assembler::ElementAssemblyValues::basis_values, centers_, polyfem::assembler::ElementAssemblyValues::has_parameterization, polyfem::assembler::Assembler::is_tensor(), is_volume(), polyfem::quadrature::Quadrature::points, setup_monomials_strong_2d(), setup_monomials_vals_2d(), and polyfem::quadrature::Quadrature::weights.

Referenced by compute_weights().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_constraints_matrix_2d_old()

void RBFWithQuadratic::compute_constraints_matrix_2d_old ( const int  num_bases,
const quadrature::Quadrature quadr,
const Eigen::MatrixXd &  local_basis_integral,
Eigen::MatrixXd &  L,
Eigen::MatrixXd &  t 
) const
private

Definition at line 412 of file RBFWithQuadratic.cpp.

References centers_, is_volume(), polyfem::quadrature::Quadrature::points, polyfem::quadrature::Quadrature::weights, and weights_.

Here is the call graph for this function:

◆ compute_constraints_matrix_3d()

void RBFWithQuadratic::compute_constraints_matrix_3d ( const assembler::LinearAssembler assembler,
const int  num_bases,
const quadrature::Quadrature quadr,
const Eigen::MatrixXd &  local_basis_integral,
Eigen::MatrixXd &  L,
Eigen::MatrixXd &  t 
) const
private

Definition at line 614 of file RBFWithQuadratic.cpp.

References centers_, is_volume(), polyfem::quadrature::Quadrature::points, polyfem::quadrature::Quadrature::weights, and weights_.

Referenced by compute_weights().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_kernels_matrix()

void RBFWithQuadratic::compute_kernels_matrix ( const Eigen::MatrixXd &  samples,
Eigen::MatrixXd &  A 
) const
private

Definition at line 383 of file RBFWithQuadratic.cpp.

References centers_, is_volume(), and x.

Referenced by bases_values(), and compute_weights().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ compute_weights()

void RBFWithQuadratic::compute_weights ( const assembler::LinearAssembler assembler,
const Eigen::MatrixXd &  collocation_points,
const Eigen::MatrixXd &  local_basis_integral,
const quadrature::Quadrature quadr,
Eigen::MatrixXd &  rhs,
bool  with_constraints 
)
private

Definition at line 714 of file RBFWithQuadratic.cpp.

References basis(), centers_, compute_constraints_matrix_2d(), compute_constraints_matrix_3d(), compute_kernels_matrix(), grad(), is_volume(), polyfem::logger(), polyfem::quadrature::Quadrature::points, val, polyfem::quadrature::Quadrature::weights, weights_, and x.

Referenced by RBFWithQuadratic().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ grad()

void RBFWithQuadratic::grad ( const int  local_index,
const Eigen::MatrixXd &  uv,
Eigen::MatrixXd &  val 
) const

Evaluates the gradient of one RBF function over a list of coordinates.

Parameters
[in]local_indexi-th RBF function to evaluate
[in]uv#uv x dim matrix of coordinates to evaluate (in object domain)
[out]val#uv x dim matrix of computed gradients

Definition at line 210 of file RBFWithQuadratic.cpp.

References bases_grads(), centers_, and val.

Referenced by compute_weights().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ index_mapping()

static int polyfem::basis::RBFWithQuadratic::index_mapping ( const int  alpha,
const int  beta,
const int  d,
const int  ass_dim 
)
inlinestatic

Definition at line 16 of file RBFWithQuadratic.hpp.

Referenced by polyfem::basis::RBFWithQuadraticLagrange::compute_constraints_matrix_2d(), and polyfem::basis::PolygonalBasis2d::compute_integral_constraints().

Here is the caller graph for this function:

◆ is_volume()

bool polyfem::basis::RBFWithQuadratic::is_volume ( ) const
inlineprivate

Definition at line 74 of file RBFWithQuadratic.hpp.

References centers_.

Referenced by bases_grads(), compute_constraints_matrix_2d(), compute_constraints_matrix_2d_old(), compute_constraints_matrix_3d(), compute_kernels_matrix(), and compute_weights().

Here is the caller graph for this function:

◆ setup_monomials_strong_2d()

void RBFWithQuadratic::setup_monomials_strong_2d ( const int  dim,
const assembler::LinearAssembler assembler,
const Eigen::MatrixXd &  pts,
const QuadratureVector da,
std::array< Eigen::MatrixXd, 5 > &  strong 
)
static

Definition at line 83 of file RBFWithQuadratic.cpp.

References polyfem::assembler::Assembler::compute_rhs(), da, and DiffScalarBase::setVariableCount().

Referenced by compute_constraints_matrix_2d(), polyfem::basis::RBFWithQuadraticLagrange::compute_constraints_matrix_2d(), and polyfem::basis::PolygonalBasis2d::compute_integral_constraints().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ setup_monomials_vals_2d()

void RBFWithQuadratic::setup_monomials_vals_2d ( const int  star_index,
const Eigen::MatrixXd &  pts,
assembler::ElementAssemblyValues vals 
)
static

Definition at line 139 of file RBFWithQuadratic.cpp.

References polyfem::assembler::ElementAssemblyValues::basis_values, and vals.

Referenced by compute_constraints_matrix_2d(), polyfem::basis::RBFWithQuadraticLagrange::compute_constraints_matrix_2d(), and polyfem::basis::PolygonalBasis2d::compute_integral_constraints().

Here is the caller graph for this function:

Member Data Documentation

◆ centers_

Eigen::MatrixXd polyfem::basis::RBFWithQuadratic::centers_
private

◆ weights_

Eigen::MatrixXd polyfem::basis::RBFWithQuadratic::weights_
private

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