G+Smo
24.08.0
Geometry + Simulation Modules
|
Interface for the set of functions defined on a domain (the total number of functions in the set equals to \(S\) )
All G+Smo objects that evaluate function[s] derive from this object. Examples are gsMatrix, gsFunction and gsGeometry.
The available evaluation procedures are:
Name of procedure | Evaluate what |
---|---|
eval (points) | value |
deriv (points) | first derivative(s) |
deriv2 (points) | second derivative(s) |
There are two classes of functions for evaluation that differ in the arguments and in the returned values. The argument can either be
1) a matrix;
2) a gsTransform object representing the parametrization of the domain.
In the first case the the matrix specifies a set of points of a size $N$ (one point per column) and it is the caller responsibility to distinguish between functions defined on the physical and on the parametric domain. In the second case the evaluator will fetch the appropriate points depending on whether it is defined on the parametric or in the physical space. In complex situations: in particular when the physical domain of a function is the parametric domain of another the automatic mechanism is not guaranteed to work.
The result is a matrix in which the i-th column contains the requested values or derivatives at the i-th point, i.e. the point whose coordinate are the i-th column of the input matrix or of the matrix supplied by the gsTransform object. On each column the data is grouped in blocks corresponding to different functions, so that that if the requested evaluation contains s values \( v_1, \ldots, v_s \) for each pair \( (f_i, p_j) \) , \( i = 1, \ldots, S \), \( j = 1, \ldots, N \) , of (function, point) the output matrix looks like
\[ \left[ \begin{array}{ccccc} v_1(f_1,p_1) & v_1(f_1,p_2) & \ldots & v_1(f_1,p_N)\\ v_2(f_1,p_1) & v_2(f_1,p_2) & \ldots & v_2(f_1,p_N)\\ \vdots & \vdots & & \vdots\\ v_s(f_1,p_1) & v_s(f_1,p_2) & \ldots & v_s(f_1,p_N)\\ v_1(f_2,p_1) & v_1(f_2,p_2) & \ldots & v_1(f_2,p_N)\\ \vdots & \vdots & & \vdots\\ v_s(f_S,p_1) & v_s(f_S,p_2) & \ldots & v_s(f_S,p_N) \end{array} \right] \]
Any implementation of the gsFunctionSet interface must:
1) overload the virtual function int size() that returns the number of functions in the set;
2) overload the needed evaluation functions with gsMatrix argument: eval_into, deriv_into, deriv2_into; the one that are not implemented will fail at runtime printing, which function need implementing to the console
and possibly:
3) write optimized versions of div, curl, laplacian, and ect. evaluation. By default they are computed by evaluating the derivatives and applying the definition.
T | type for real numbers |
Public Types | |
typedef memory::shared_ptr < gsFunctionSet > | Ptr |
Shared pointer for gsFunctionSet. | |
typedef memory::unique_ptr < gsFunctionSet > | uPtr |
Unique pointer for gsFunctionSet. | |
Public Member Functions | |
gsMatrix< index_t > | active (const gsMatrix< T > &u) const |
Returns the indices of active (nonzero) functions at points u, as a list of indices. More... | |
virtual void | active_into (const gsMatrix< T > &u, gsMatrix< index_t > &result) const |
Indices of active (non-zero) function(s) for each point. More... | |
const gsBasis< T > & | basis (const index_t k) const |
Helper which casts and returns the k-th piece of this function set as a gsBasis. | |
uPtr | clone () |
Clone methode. Produceds a deep copy inside a uPtr. | |
virtual void | compute (const gsMatrix< T > &in, gsFuncData< T > &out) const |
Computes function data. More... | |
gsMatrix< T > | deriv (const gsMatrix< T > &u) const |
Evaluate the derivatives,. More... | |
gsMatrix< T > | deriv2 (const gsMatrix< T > &u) const |
Evaluates the second derivatives of active (i.e., non-zero) functions at points u. More... | |
virtual void | deriv2_into (const gsMatrix< T > &u, gsMatrix< T > &result) const |
Second derivatives. More... | |
virtual void | deriv_into (const gsMatrix< T > &u, gsMatrix< T > &result) const |
First derivatives. More... | |
virtual short_t | domainDim () const =0 |
Dimension of the (source) domain. More... | |
gsMatrix< T > | eval (const gsMatrix< T > &u) const |
Evaluate the function,. More... | |
virtual void | eval_into (const gsMatrix< T > &u, gsMatrix< T > &result) const |
Evaluates the function(s). More... | |
std::vector< gsMatrix< T > > | evalAllDers (const gsMatrix< T > &u, int n, bool sameElement=false) const |
Evaluate all derivatives upto order n,. More... | |
virtual void | evalAllDers_into (const gsMatrix< T > &u, int n, std::vector< gsMatrix< T > > &result, bool sameElement=false) const |
Evaluate the nonzero functions and their derivatives up to order n at points u into result. More... | |
const gsFunction< T > & | function (const index_t k) const |
Helper which casts and returns the k-th piece of this function set as a gsFunction. | |
virtual index_t | nPieces () const |
Number of pieces in the domain of definition. | |
virtual const gsFunctionSet & | piece (const index_t) const |
Returns the piece(s) of the function(s) at subdomain k. | |
virtual std::ostream & | print (std::ostream &os) const |
Prints the object as a string. | |
virtual index_t | size () const |
size More... | |
virtual short_t | targetDim () const |
Dimension of the target space. More... | |
Returns the indices of active (nonzero) functions at points u, as a list of indices.
Indices of active (non-zero) function(s) for each point.
The columns are sorted in increasing order, if on a point there are less active then the number of rows in the result matrix (some other point has more actives) then the rest of the column is filled with 0s.
u | |
result |
Reimplemented in gsHTensorBasis< d, T >, gsBasis< T >, gsBasis< Scalar >, gsBasis< real_t >, gsTensorBSplineBasis< 1, T >, gsTensorBSplineBasis< d, T >, gsTensorBSplineBasis< domainDim+1, T >, gsTensorBasis< d, T >, gsTensorBasis< 1, T >, gsTHBSplineBasis< d, T >, gsLagrangeBasis< T >, gsMappedSingleBasis< d, T >, gsFunction< T >, gsConstantBasis< T >, gsLegendreBasis< T >, gsMonomialBasis< T >, and gsMvLegendreBasis< T >.
|
virtual |
Computes function data.
This function evaluates the functions and their derivatives at the points in and writes them in the corresponding fields of out. Which field to write (and what to compute) is controlled by the out.flags (see also gsFuncData).
The input points in are expected to be compatible with the implementation/representation of the function, i.e. they should be points inside the domain of definitition of the function
[in] | in | |
[out] | out |
Reimplemented in gsGeometry< T >, and gsConstantFunction< T >.
Evaluate the derivatives,.
Evaluates the second derivatives of active (i.e., non-zero) functions at points u.
See documentation for deriv2_into() (the one without input parameter coefs) for details.
[in] | u | Evaluation points in columns. |
Second derivatives.
For scalar valued functions \(f_1, \ldots, f_S\) from \(\mathbb{R}^n\rightarrow\mathbb{R}\) format is:
\[ \left[ \begin{array}{ccccc} \partial_{1}\partial_{1}f_1(p_1) & \partial_{1}\partial_{1}f_1(p_2) & \ldots & \partial_{1}\partial_{1}f_1(p_N)\\ \partial_{2}\partial_{2}f_1(p_1) & \partial_{2}\partial_{2}f_1(p_2) & \ldots & \partial_{2}\partial_{2}f_1(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k}\partial_{k}f_1(p_1) & \partial_{k}\partial_{k}f_1(p_2) & \ldots & \partial_{k}\partial_{k}f_1(p_N)\\ \partial_{1}\partial_{2}f_1(p_1) & \partial_{1}\partial_{2}f_1(p_2) & \ldots & \partial_{1}\partial_{2}f_1(p_N)\\ \partial_{1}\partial_{3}f_1(p_1) & \partial_{1}\partial_{3}f_1(p_2) & \ldots & \partial_{1}\partial_{3}f_1(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{1}\partial_{k}f_1(p_1) & \partial_{1}\partial_{k}f_1(p_2) & \ldots & \partial_{1}\partial_{k}f_1(p_N)\\ \partial_{2}\partial_{3}f_1(p_1) & \partial_{2}\partial_{3}f_1(p_2) & \ldots & \partial_{2}\partial_{3}f_1(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{2}\partial_{k}f_1(p_1) & \partial_{2}\partial_{k}f_1(p_2) & \ldots & \partial_{2}\partial_{k}f_1(p_N)\\ \partial_{3}\partial_{4}f_1(p_1) & \partial_{3}\partial_{4}f_1(p_2) & \ldots & \partial_{3}\partial_{4}f_1(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k-1}\partial_{k}f_1(p_1) & \partial_{k-1}\partial_{k}f_1(p_2) & \ldots & \partial_{k-1}\partial_{k}f_1(p_N)\\ \partial_{1}\partial_{1}f_2(p_1) & \partial_{1}\partial_{1}f_2(p_2) & \ldots & \partial_{1}\partial_{1}f_2(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k-1}\partial_{k}f_S(p_1) & \partial_{k-1}\partial_{k}f_S(p_2) & \ldots & \partial_{k-1}\partial_{k}f_S(p_N)\\ \end{array} \right] \]
For vector valued functions function \(f_1, \ldots, f_S\) from \(\mathbb{R}^n\rightarrow\mathbb{R}^{m}\) the format is:
\[ \left[ \begin{array}{ccccc} \partial_{1}\partial_{1}f_1^{(1)}(p_1) & \partial_{1}\partial_{1}f_1^{(1)}(p_2) & \ldots & \partial_{1}\partial_{1}f_1^{(1)}(p_N)\\ \partial_{2}\partial_{2}f_1^{(1)}(p_1) & \partial_{2}\partial_{2}f_1^{(1)}(p_2) & \ldots & \partial_{2}\partial_{2}f_1^{(1)}(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k}\partial_{k}f_1^{(1)}(p_1) & \partial_{k}\partial_{k}f_1^{(1)}(p_2) & \ldots & \partial_{k}\partial_{k}f_1^{(1)}(p_N)\\ \partial_{1}\partial_{2}f_1^{(1)}(p_1) & \partial_{1}\partial_{2}f_1^{(1)}(p_2) & \ldots & \partial_{1}\partial_{2}f_1^{(1)}(p_N)\\ \partial_{1}\partial_{3}f_1^{(1)}(p_1) & \partial_{1}\partial_{3}f_1^{(1)}(p_2) & \ldots & \partial_{1}\partial_{3}f_1^{(1)}(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k-1}\partial_{k}f_1{(1)}(p_1) & \partial_{k-1}\partial_{k}f_1^{(1)}(p_2) & \ldots & \partial_{k-1}\partial_{k}f_1^{(1)}(p_N)\\ \partial_{1}\partial_{1}f_1^{(2)}(p_1) & \partial_{1}\partial_{1}f_1^{(2)}(p_2) & \ldots & \partial_{1}\partial_{1}f_1^{(2)}(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k-1}\partial_{k}f_1^{(m)}(p_1) & \partial_{k-1}\partial_{k}f_1^{(m)}(p_2) & \ldots & \partial_{k-1}\partial_{k}f_1^{(m)}(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k-1}\partial_{k}f_S^{(m)}(p_1) & \partial_{k-1}\partial_{k}f_S^{(m)}(p_2) & \ldots & \partial_{k-1}\partial_{k}f_S^{(m)}(p_N)\\ \end{array} \right] \]
where \( f^{(i)}_j\) is the \(i\)-th component of function \( f_j\) of the set.
u | |
result |
Reimplemented in gsBasis< T >, gsBasis< Scalar >, gsBasis< real_t >, gsTensorBSplineBasis< 1, T >, gsTensorBasis< d, T >, gsTensorBasis< 1, T >, gsGeometry< T >, gsFunctionExpr< T >, gsFunction< T >, gsLagrangeBasis< T >, gsMappedSingleBasis< d, T >, gsConstantFunction< T >, gsTHBSplineBasis< d, T >, gsPreCICEFunction< T >, gsLegendreBasis< T >, gsAffineFunction< T >, gsHBSplineBasis< d, T >, gsMappedSingleSpline< d, T >, gsMonomialBasis< T >, gsMvLegendreBasis< T >, gsGeometryTransform< T >, and gsFuncCoordinate< T >.
First derivatives.
For scalar valued functions \(f_1, \ldots, f_S\) from \(\mathbb{R}^n\rightarrow\mathbb{R}\) format is:
\[ \left[ \begin{array}{ccccc} \partial_{1}f_1(p_1) & \partial_{1}f_1(p_2) & \ldots & \partial_{1}f_1(p_N)\\ \partial_{2}f_1(p_1) & \partial_{2}f_1(p_2) & \ldots & \partial_{2}f_1(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k}f_1(p_1) & \partial_{k}f_1(p_2) & \ldots & \partial_{k}f_1(p_N)\\ \partial_{1}f_2(p_1) & \partial_{1}f_2(p_2) & \ldots & \partial_{1}f_2(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k}f_S(p_1) & \partial_{k}f_S(p_2) & \ldots & \partial_{k}f_S(p_N)\\ \end{array} \right] \]
For vector valued functions function \(f_1, \ldots, f_S\) from \(\mathbb{R}^n\rightarrow\mathbb{R}^{m}\) the format is:
\[ \left[ \begin{array}{ccccc} \partial_{1}f_1^1(p_1) & \partial_{1}f_1^1(p_2) & \ldots & \partial_{1}f_1^1(p_N)\\ \partial_{2}f_1^1(p_1) & \partial_{1}f_2^1(p_2) & \ldots & \partial_{1}f_2^1(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k}f_1^1(p_1) & \partial_{k}f_1^1(p_2) & \ldots & \partial_{k}f_1^1(p_N)\\ \partial_{1}f_1^2(p_1) & \partial_{1}f_1^2(p_2) & \ldots & \partial_{1}f_1^2(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k}f_1^2(p_1) & \partial_{k}f_1^2(p_2) & \ldots & \partial_{k}f_1^2(p_N)\\ \partial_{1}f_2^1(p_1) & \partial_{1}f_2^1(p_2) & \ldots & \partial_{1}f_2^1(p_N)\\ \vdots & \vdots & & \vdots\\ \partial_{k}f_S^{(m)}(p_1) & \partial_{k}f_S^{(m)}(p_2) & \ldots & \partial_{k}f_S^{(m)}(p_N) \end{array} \right] \]
where \(f^{(i)}_j\) is the \(i\)-th component of function \(f_j\) of the set.
u | |
result |
Reimplemented in gsBasis< T >, gsBasis< Scalar >, gsBasis< real_t >, gsTensorBasis< d, T >, gsTensorBasis< 1, T >, gsTensorBSplineBasis< 1, T >, gsGeometry< T >, gsFunctionExpr< T >, gsLagrangeBasis< T >, gsMappedSingleBasis< d, T >, gsFunction< T >, gsConstantFunction< T >, gsTHBSplineBasis< d, T >, gsPreCICEFunction< T >, gsAffineFunction< T >, gsHBSplineBasis< d, T >, gsLegendreBasis< T >, gsMappedSingleSpline< d, T >, gsMvLegendreBasis< T >, gsMonomialBasis< T >, gsConstantBasis< T >, gsGeometrySlice< T >, gsSquaredDistance< T >, gsBasisFun< T >, gsGeometryTransform< T >, and gsFuncCoordinate< T >.
|
pure virtual |
Dimension of the (source) domain.
Implemented in gsTHBSplineBasis< d, T >, gsPatchIdField< T >, gsGeometry< T >, gsParamField< T >, gsNormalField< T >, gsMultiBasis< T >, gsMultiBasis< real_t >, gsMultiPatch< T >, gsMultiPatch< real_t >, gsFsiLoad< T >, gsMaterialMatrixEvalSingle< T, out >, gsFunctionExpr< T >, gsMaterialMatrixIntegrateSingle< T, out >, gsTensorBSplineBasis< 1, T >, gsJacDetField< T >, gsDetFunction< T >, gsConstantFunction< T >, gsRationalBasis< SrcT >, gsRationalBasis< gsBSplineTraits< d, T >::Basis >, gsRationalBasis< gsBSplineBasis< T > >, gsRemapInterface< T >, gsGradientField< T >, gsShellStressFunction< T >, gsPiecewiseFunction< T >, gsLagrangeBasis< T >, gsTensorBasis< d, T >, gsTensorBasis< 1, T >, gsCPPInterface< T >, gsSquaredDistance< T >, gsHBSplineBasis< d, T >, gsAffineFunction< T >, gsPreCICEFunction< T >, gsElementErrorPlotter< T >, gsMappedSingleBasis< d, T >, gsMappedSingleSpline< d, T >, gsConstantBasis< T >, gsSurface< T >, gsLegendreBasis< T >, gsBasisFun< T >, gsCurve< T >, gsAbsError< T >, gsVolume< T >, gsBulk< T >, gsGeometryTransform< T >, gsFuncCoordinate< T >, gsGeometrySlice< T >, gsCauchyStressFunction< T >, and gsMonomialBasis< T >.
Evaluate the function,.
Evaluates the function(s).
For scalar valued functions \(f_1, \ldots, f_S\) from \(\mathbb{R}^n\rightarrow\mathbb{R}\) format is:
\[ \left[ \begin{array}{ccccc} f_1(p_1) & f_1(p_2) & \ldots & f_1(p_N)\\ f_2(p_1) & f_2(p_2) & \ldots & f_2(p_N)\\ \vdots & \vdots & & \vdots\\ f_S(p_1) & f_S(p_2) & \ldots & f_S(p_N) \end{array} \right] \]
For vector valued functions function \(f_1, \ldots, f_S\) from \(\mathbb{R}^n\rightarrow\mathbb{R}^m\) the format is:
\[ \left[ \begin{array}{ccccc} f_1^1(p_1) & f_1^{(1)}(p_2) & \ldots & f_1^{(1)}(p_N)\\ f_1^2(p_1) & f_1^{(2)}(p_2) & \ldots & f_1^{(2)}(p_N)\\ \vdots & \vdots & & \vdots\\ f_1^{(m)}(p_1) & f_1^{(m)}(p_2) & \ldots & f_1^{(m)}(p_N)\\ f_2^{(1)}(p_1) & f_2^{(1)}(p_2) & \ldots & f_2^{(1)}(p_N)\\ \vdots & \vdots & & \vdots\\ f_S^{(m)}(p_1) & f_S^{(m)}(p_2) & \ldots & f_S^{(m)}(p_N) \end{array} \right] \]
where \(f^{(i)}_j\) is the \(i\)-th component of function \(f_j\) of the set.
u | |
result |
Reimplemented in gsBasis< T >, gsBasis< Scalar >, gsBasis< real_t >, gsTHBSplineBasis< d, T >, gsPatchIdField< T >, gsMaterialMatrixEvalSingle< T, out >, gsParamField< T >, gsTensorBSplineBasis< 1, T >, gsTensorBasis< d, T >, gsTensorBasis< 1, T >, gsShellStressFunction< T >, gsMaterialMatrixIntegrateSingle< T, out >, gsNormalField< T >, gsFsiLoad< T >, gsFunctionExpr< T >, gsGeometry< T >, gsLagrangeBasis< T >, gsMappedSingleBasis< d, T >, gsPiecewiseFunction< T >, gsConstantFunction< T >, gsJacDetField< T >, gsDetFunction< T >, gsRemapInterface< T >, gsFunction< T >, gsGradientField< T >, gsCPPInterface< T >, gsHBSplineBasis< d, T >, gsAffineFunction< T >, gsPreCICEFunction< T >, gsLegendreBasis< T >, gsMvLegendreBasis< T >, gsMappedSingleSpline< d, T >, gsConstantBasis< T >, gsMonomialBasis< T >, gsCauchyStressFunction< T >, gsBasisFun< T >, gsGeometrySlice< T >, gsGeometryTransform< T >, gsFuncCoordinate< T >, gsAbsError< T >, gsElementErrorPlotter< T >, and gsSquaredDistance< T >.
std::vector< gsMatrix< T > > evalAllDers | ( | const gsMatrix< T > & | u, |
int | n, | ||
bool | sameElement = false |
||
) | const |
Evaluate all derivatives upto order n,.
|
virtual |
Evaluate the nonzero functions and their derivatives up to order n at points u into result.
The derivatives (the 0-th derivative is the function value) are stored in a result. result is a std::vector, where result[i] is a gsMatrix which contains the i-th derivatives.
The entries in result[0], result[1], and result[2] are ordered as in eval_into(), deriv_into(), and deriv2_into(), respectively. For i > 2, the derivatives are stored in lexicographical order, e.g. for order i = 3 and dimension 2 the derivatives are stored as follows: \( \partial_{xxx}, \, \partial_{xxy}, \, \partial_{xyy}, \, \partial_{yyy}.\, \)
[in] | u | Evaluation points, each column corresponds to one evaluation point. |
[in] | n | All derivatives up to order n are computed and stored in result. |
[in,out] | result | See above for format. |
Reimplemented in gsTensorBSplineBasis< 1, T >, gsTensorBasis< d, T >, gsTensorBasis< 1, T >, gsGeometry< T >, gsMappedSingleBasis< d, T >, gsConstantFunction< T >, gsTHBSplineBasis< d, T >, gsMappedSingleSpline< d, T >, and gsSquaredDistance< T >.
|
inlinevirtual |
size
Reimplemented in gsHTensorBasis< d, T >, gsFunction< T >, gsMappedSingleBasis< d, T >, gsMultiBasis< T >, gsLagrangeBasis< T >, gsMultiBasis< real_t >, gsTensorBSplineBasis< 1, T >, gsRationalBasis< SrcT >, gsRationalBasis< gsBSplineTraits< d, T >::Basis >, gsRationalBasis< gsBSplineBasis< T > >, gsPiecewiseFunction< T >, gsMultiPatch< T >, gsMultiPatch< real_t >, gsLegendreBasis< T >, gsMvLegendreBasis< T >, gsTensorBasis< d, T >, gsTensorBasis< 1, T >, gsConstantBasis< T >, gsMonomialBasis< T >, and gsPatchwiseFunction< T >.
|
inlinevirtual |
Dimension of the target space.
Reimplemented in gsPatchIdField< T >, gsGeometry< T >, gsParamField< T >, gsNormalField< T >, gsMultiBasis< T >, gsMultiBasis< real_t >, gsMultiPatch< T >, gsMultiPatch< real_t >, gsFsiLoad< T >, gsMaterialMatrixEvalSingle< T, out >, gsMaterialMatrixIntegrateSingle< T, out >, gsFunctionExpr< T >, gsJacDetField< T >, gsDetFunction< T >, gsConstantFunction< T >, gsShellStressFunction< T >, gsGradientField< T >, gsPiecewiseFunction< T >, gsSquaredDistance< T >, gsAffineFunction< T >, gsPreCICEFunction< T >, gsMappedSingleSpline< d, T >, gsBasisFun< T >, gsAbsError< T >, gsGeometrySlice< T >, gsGeometryTransform< T >, gsFuncCoordinate< T >, and gsCauchyStressFunction< T >.