Isogeometric analysis utilities.
The functions compute_bezier_extraction_1d() and eval_nurbs_basis_tp() implement the algorithms described in [1].
For a nD B-spline parametric domain, combine the 1D element extraction operators in each parametric dimension into a single operator for each nD element.
Parameters: | cs : list of lists of 2D arrays
|
---|---|
Returns: | ccs : list of 2D arrays
|
Compute the control points and weights of the Bezier mesh.
Parameters: | control_points : array
weights : array
ccs : list of 2D arrays
conn : array
bconn : array
|
---|---|
Returns: | bezier_control_points : array
bezier_weights : array
|
Compute local (element) Bezier extraction operators for a nD B-spline parametric domain.
Parameters: | knots : sequence of array or array
degrees : sequence of ints or int
|
---|---|
Returns: | cs : list of lists of 2D arrays
|
Compute local (element) Bezier extraction operators for a 1D B-spline parametric domain.
Parameters: | knots : array
degree : int
|
---|---|
Returns: | cs : array of 2D arrays (3D array)
|
Create boundary quadrature points from the surface quadrature points.
Uses the Bezier element tensor product structure.
Parameters: | coors : array, shape (n_qp, d)
dim : int
|
---|---|
Returns: | bcoors : array, shape (n_qp, d + 1)
|
Create connectivity arrays of nD Bezier elements.
Parameters: | n_els : sequence of ints
knots : sequence of array or array
degrees : sequence of ints or int
|
---|---|
Returns: | conn : array
bconn : array
|
Create connectivity arrays of 1D Bezier elements.
Parameters: | n_el : int
knots : array
degree : int
|
---|---|
Returns: | conn : array
bconn : array
|
Evaluate the Bernstein polynomial basis of the given degree, and its derivatives, in a point x in [0, 1].
Parameters: | x : float
degree : int
|
---|---|
Returns: | funs : array
ders : array
|
Evaluate data required for the isogeometric domain reference mapping in the given quadrature points. The quadrature points are the same for all Bezier elements and should correspond to the Bernstein basis degree.
Parameters: | qps : array
control_points : array
weights : array
degrees : sequence of ints or int
cs : list of lists of 2D arrays
conn : array
cells : array, optional
|
---|---|
Returns: | bfs : array
bfgs : array
dets : array
|
Evaluate the tensor-product NURBS shape functions in a quadrature point for a given Bezier element.
Parameters: | qp : array
ie : int
control_points : array
weights : array
degrees : sequence of ints or int
cs : list of lists of 2D arrays
conn : array
|
---|---|
Returns: | R : array
dR_dx : array
det : array
|
Evaluate a field variable in the given quadrature points. The quadrature points are the same for all Bezier elements and should correspond to the Bernstein basis degree. The field variable is defined by its DOFs - the coefficients of the NURBS basis.
Parameters: | variable : array
qps : array
control_points : array
weights : array
degrees : sequence of ints or int
cs : list of lists of 2D arrays
conn : array
cells : array, optional
|
---|---|
Returns: | coors : array
vals : array
dets : array
|
Get faces and edges of a Bezier mesh element in terms of indices into the element’s connectivity (reference Bezier element entities).
Parameters: | degrees : sequence of ints or int
|
---|---|
Returns: | faces : list of arrays
edges : list of arrays
vertices : list of arrays
|
Notes
The ordering of faces and edges has to be the same as in sfepy.discrete.fem.geometry_element.geometry_data.
Get a topology connectivity corresponding to the Bezier mesh connectivity.
In the referenced Bezier control points the Bezier mesh is interpolatory.
Parameters: | bconn : array
degrees : sequence of ints or int
|
---|---|
Returns: | tconn : array
|
For each reference Bezier element facet return the facet axes followed by the remaining (perpendicular) axis, as well as the remaining axis coordinate of the facet.
Parameters: | dim : int
|
---|---|
Returns: | axes : array
coors : array
|
Get box regions of Bezier topological mesh in terms of element corner vertices of Bezier mesh.
Parameters: | n_els : sequence of ints
degrees : sequence of ints or int
|
---|---|
Returns: | regions : dict
|
Get a global raveled index corresponding to nD indices into an array of the given shape.
Get degrees of the NURBS patch surfaces.
Parameters: | degrees : sequence of ints or int
|
---|---|
Returns: | sdegrees : list of arrays
|