Dune Core Modules (2.9.0)
bsplinebasis.hh
Go to the documentation of this file.
28 // A maze of dependencies between the different parts of this. We need a few forward declarations
402 while (preBasis_.knotVectors_[i][currentKnotSpan_[i]+1] < preBasis_.knotVectors_[i][currentKnotSpan_[i]]+1e-8)
410 while (preBasis_.knotVectors_[i][currentKnotSpan_[i]+1] < preBasis_.knotVectors_[i][currentKnotSpan_[i]]+1e-8)
416 localBasis_.scaling_[i][i] = preBasis_.knotVectors_[i][currentKnotSpan_[i]+1] - preBasis_.knotVectors_[i][currentKnotSpan_[i]];
601 assert( std::accumulate(elements_.begin(), elements_.end(), 1, std::multiplies<unsigned>()) == gridView_.size(0) );
653 assert( std::accumulate(elements_.begin(), elements_.end(), 1, std::multiplies<unsigned>()) == gridView_.size(0) );
751 globalIJK[i] = std::max((int)currentKnotSpan[i] - (int)order[i], 0) + localIJK[i]; // needs to be a signed type!
862 if (order_[i]==0) // order-zero functions are piecewise constant, hence all derivatives are zero
868 R derivativeAddend1 = lowOrderOneDValues[i][j] / (knotVectors_[i][j+order_[i]]-knotVectors_[i][j]);
869 R derivativeAddend2 = lowOrderOneDValues[i][j+1] / (knotVectors_[i][j+order_[i]+1]-knotVectors_[i][j+1]);
929 evaluateAll(in[i], oneDValues[i], true, oneDDerivatives[i], false, oneDSecondDerivatives[i], knotVectors_[i], order_[i], currentKnotSpan[i]);
932 evaluateAll(in[i], oneDValues[i], true, oneDDerivatives[i], true, oneDSecondDerivatives[i], knotVectors_[i], order_[i], currentKnotSpan[i]);
990 if (directions[0] == j || directions[1] == j) //the spline has to be derived (once) in this direction
1009 static std::array<unsigned int,dim> getIJK(typename GridView::IndexSet::IndexType idx, std::array<unsigned int,dim> elements)
1033 std::size_t outSize = order+1; // The 'standard' value away from the boundaries of the knot vector
1212 R derivativeAddend1 = lowOrderTwoDValues[j] / (knotVector[j+order]-knotVector[j]) / (knotVector[j+order-1]-knotVector[j]);
1213 R derivativeAddend2 = lowOrderTwoDValues[j+1] / (knotVector[j+order]-knotVector[j]) / (knotVector[j+order]-knotVector[j+1]);
1214 R derivativeAddend3 = lowOrderTwoDValues[j+1] / (knotVector[j+order+1]-knotVector[j+1]) / (knotVector[j+order]-knotVector[j+1]);
1215 R derivativeAddend4 = lowOrderTwoDValues[j+2] / (knotVector[j+order+1]-knotVector[j+1]) / (knotVector[j+1+order]-knotVector[j+2]);
1226 outHess[j-offset] = order * (order-1) * ( derivativeAddend1 - derivativeAddend2 -derivativeAddend3 + derivativeAddend4 );
1311 return BSplinePreBasis<std::decay_t<decltype(gridView)>>(gridView, knotVector, order, makeOpen);
void resize(size_type r, size_type c, value_type v=value_type())
resize matrix to r × c
Definition: dynmatrix.hh:106
LocalBasis class in the sense of dune-localfunctions, presenting the restriction of a B-spline patch ...
Definition: bsplinebasis.hh:46
LocalBasisTraits< D, dim, FieldVector< D, dim >, R, 1, FieldVector< R, 1 >, FieldMatrix< R, 1, dim > > Traits
export type traits for function signature
Definition: bsplinebasis.hh:55
unsigned int order() const
Polynomial order of the shape functions.
Definition: bsplinebasis.hh:140
std::size_t size() const
Return the number of basis functions on the current knot span.
Definition: bsplinebasis.hh:147
void evaluate(const typename std::array< int, k > &directions, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions and derivatives of any order.
Definition: bsplinebasis.hh:97
void evaluateFunction(const FieldVector< D, dim > &in, std::vector< FieldVector< R, 1 > > &out) const
Evaluate all shape functions.
Definition: bsplinebasis.hh:70
void evaluateJacobian(const FieldVector< D, dim > &in, std::vector< FieldMatrix< D, 1, dim > > &out) const
Evaluate Jacobian of all shape functions.
Definition: bsplinebasis.hh:82
BSplineLocalBasis(const BSplinePreBasis< GV > &preBasis, const BSplineLocalFiniteElement< GV, R > &lFE)
Constructor with a given B-spline patch.
Definition: bsplinebasis.hh:61
Attaches a shape function to an entity.
Definition: bsplinebasis.hh:178
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: bsplinebasis.hh:321
std::size_t size() const
number of coefficients
Definition: bsplinebasis.hh:315
LocalFiniteElement in the sense of dune-localfunctions, for the B-spline basis on tensor-product grid...
Definition: bsplinebasis.hh:362
BSplineLocalFiniteElement(const BSplineLocalFiniteElement &other)
Copy constructor.
Definition: bsplinebasis.hh:383
const BSplineLocalBasis< GV, R > & localBasis() const
Hand out a LocalBasis object.
Definition: bsplinebasis.hh:427
LocalFiniteElementTraits< BSplineLocalBasis< GV, R >, BSplineLocalCoefficients< dim >, BSplineLocalInterpolation< dim, BSplineLocalBasis< GV, R > > > Traits
Export various types related to this LocalFiniteElement.
Definition: bsplinebasis.hh:372
BSplineLocalFiniteElement(const BSplinePreBasis< GV > &preBasis)
Constructor with a given B-spline basis.
Definition: bsplinebasis.hh:376
void bind(const std::array< unsigned, dim > &elementIdx)
Bind LocalFiniteElement to a specific knot span of the spline patch.
Definition: bsplinebasis.hh:394
const BSplineLocalCoefficients< dim > & localCoefficients() const
Hand out a LocalCoefficients object.
Definition: bsplinebasis.hh:433
GeometryType type() const
Return the reference element that the local finite element is defined on (here, a hypercube)
Definition: bsplinebasis.hh:455
unsigned size() const
Number of shape functions in this finite element.
Definition: bsplinebasis.hh:445
unsigned int size(int i) const
Number of degrees of freedom for one coordinate direction.
Definition: bsplinebasis.hh:463
const BSplineLocalInterpolation< dim, BSplineLocalBasis< GV, R > > & localInterpolation() const
Hand out a LocalInterpolation object.
Definition: bsplinebasis.hh:439
Local interpolation in the sense of dune-localfunctions, for the B-spline basis on tensor-product gri...
Definition: bsplinebasis.hh:340
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: bsplinebasis.hh:344
std::array< unsigned, dim > elements_
Number of grid elements in the different coordinate directions.
Definition: bsplinebasis.hh:1240
static std::array< unsigned int, dim > getIJK(typename GridView::IndexSet::IndexType idx, std::array< unsigned int, dim > elements)
Compute integer element coordinates from the element index.
Definition: bsplinebasis.hh:1009
const GridView & gridView() const
Obtain the grid view that the basis is defined on.
Definition: bsplinebasis.hh:680
static void evaluateFunctionFull(const typename GV::ctype &in, DynamicMatrix< R > &out, const std::vector< R > &knotVector, unsigned int order, unsigned int currentKnotSpan)
Evaluate all one-dimensional B-spline functions for a given coordinate direction.
Definition: bsplinebasis.hh:1085
void evaluateFunction(const FieldVector< typename GV::ctype, dim > &in, std::vector< FieldVector< R, 1 > > &out, const std::array< unsigned, dim > ¤tKnotSpan) const
Evaluate all B-spline basis functions at a given point.
Definition: bsplinebasis.hh:781
std::array< unsigned int, dim > order_
Order of the B-spline for each space dimension.
Definition: bsplinebasis.hh:1234
size_type dimension() const
Get the total dimension of the space spanned by this basis.
Definition: bsplinebasis.hh:719
static void evaluateAll(const typename GV::ctype &in, std::vector< R > &out, bool evaluateJacobian, std::vector< R > &outJac, bool evaluateHessian, std::vector< R > &outHess, const std::vector< R > &knotVector, unsigned int order, unsigned int currentKnotSpan)
Evaluate the second derivatives of all one-dimensional B-spline functions for a given coordinate dire...
Definition: bsplinebasis.hh:1128
static void evaluateFunction(const typename GV::ctype &in, std::vector< R > &out, const std::vector< R > &knotVector, unsigned int order, unsigned int currentKnotSpan)
Evaluate all one-dimensional B-spline functions for a given coordinate direction.
Definition: bsplinebasis.hh:1028
unsigned int size(size_t d) const
Number of shape functions in one direction.
Definition: bsplinebasis.hh:774
GV GridView
The grid view that the FE space is defined on.
Definition: bsplinebasis.hh:560
It indices(const Node &node, It it) const
Maps from subtree index set [0..size-1] to a globally unique multi index in global basis.
Definition: bsplinebasis.hh:735
void evaluate(const typename std::array< int, k > &directions, const FieldVector< typename GV::ctype, dim > &in, std::vector< FieldVector< R, 1 > > &out, const std::array< unsigned, dim > ¤tKnotSpan) const
Evaluate Derivatives of all B-spline basis functions.
Definition: bsplinebasis.hh:913
void update(const GridView &gv)
Update the stored grid view, to be called if the grid has changed.
Definition: bsplinebasis.hh:686
void initializeIndices()
Initialize the global indices.
Definition: bsplinebasis.hh:676
size_type size(const Dune::ReservedVector< ST, i > &prefix) const
Return number of possible values for next position in multi index.
Definition: bsplinebasis.hh:712
unsigned int size() const
Total number of B-spline basis functions.
Definition: bsplinebasis.hh:765
Node makeNode() const
Create tree node.
Definition: bsplinebasis.hh:694
BSplinePreBasis(const GridView &gridView, const std::vector< double > &knotVector, unsigned int order, bool makeOpen=true)
Construct a B-spline basis for a given grid view and set of knot vectors.
Definition: bsplinebasis.hh:590
void evaluateJacobian(const FieldVector< typename GV::ctype, dim > &in, std::vector< FieldMatrix< R, 1, dim > > &out, const std::array< unsigned, dim > ¤tKnotSpan) const
Evaluate Jacobian of all B-spline basis functions.
Definition: bsplinebasis.hh:812
std::array< std::vector< double >, dim > knotVectors_
The knot vectors, one for each space dimension.
Definition: bsplinebasis.hh:1237
size_type maxNodeSize() const
Get the maximal number of DOFs associated to node for any element.
Definition: bsplinebasis.hh:725
BSplinePreBasis(const GridView &gridView, const FieldVector< double, dim > &lowerLeft, const FieldVector< double, dim > &upperRight, const std::array< unsigned int, dim > &elements, unsigned int order, bool makeOpen=true)
Construct a B-spline basis for a given grid view with uniform knot vectors.
Definition: bsplinebasis.hh:642
Global basis for given pre-basis.
Definition: defaultglobalbasis.hh:46
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:125
constexpr size_type size() const noexcept
Returns number of elements in the vector.
Definition: reservedvector.hh:387
This file implements a quadratic diagonal matrix of fixed size.
This file implements a dense matrix with dynamic numbers of rows and columns.
auto bSpline(const std::vector< double > &knotVector, unsigned int order, bool makeOpen=true)
Create a pre-basis factory that can create a B-spline pre-basis.
Definition: bsplinebasis.hh:1306
constexpr GeometryType cube(unsigned int dim)
Returns a GeometryType representing a hypercube of dimension dim.
Definition: type.hh:472
constexpr T accumulate(Range &&range, T value, F &&f)
Accumulate values.
Definition: hybridutilities.hh:291
auto max(ADLTag< 0 >, const V &v1, const V &v2)
implements binary Simd::max()
Definition: defaults.hh:81
constexpr HybridTreePath< T..., std::size_t > push_back(const HybridTreePath< T... > &tp, std::size_t i)
Appends a run time index to a HybridTreePath.
Definition: treepath.hh:281
A unique label for each type of element that can occur in a grid.
|
Legal Statements / Impressum |
Hosted by TU Dresden |
generated with Hugo v0.80.0
(Apr 25, 22:37, 2024)