pub struct Basis<'a> { /* private fields */ }
Implementations§
Source§impl<'a> Basis<'a>
impl<'a> Basis<'a>
pub fn create_tensor_H1( ceed: &Ceed, dim: usize, ncomp: usize, P1d: usize, Q1d: usize, interp1d: &[Scalar], grad1d: &[Scalar], qref1d: &[Scalar], qweight1d: &[Scalar], ) -> Result<Self>
pub fn create_tensor_H1_Lagrange( ceed: &Ceed, dim: usize, ncomp: usize, P: usize, Q: usize, qmode: QuadMode, ) -> Result<Self>
pub fn create_H1( ceed: &Ceed, topo: ElemTopology, ncomp: usize, nnodes: usize, nqpts: usize, interp: &[Scalar], grad: &[Scalar], qref: &[Scalar], qweight: &[Scalar], ) -> Result<Self>
pub fn create_Hdiv( ceed: &Ceed, topo: ElemTopology, ncomp: usize, nnodes: usize, nqpts: usize, interp: &[Scalar], div: &[Scalar], qref: &[Scalar], qweight: &[Scalar], ) -> Result<Self>
pub fn create_Hcurl( ceed: &Ceed, topo: ElemTopology, ncomp: usize, nnodes: usize, nqpts: usize, interp: &[Scalar], curl: &[Scalar], qref: &[Scalar], qweight: &[Scalar], ) -> Result<Self>
Sourcepub fn apply(
&self,
nelem: usize,
tmode: TransposeMode,
emode: EvalMode,
u: &Vector<'_>,
v: &mut Vector<'_>,
) -> Result<i32>
pub fn apply( &self, nelem: usize, tmode: TransposeMode, emode: EvalMode, u: &Vector<'_>, v: &mut Vector<'_>, ) -> Result<i32>
Apply basis evaluation from nodes to quadrature points or vice versa
nelem
- The number of elements to apply the basis evaluation totmode
-TrasposeMode::NoTranspose
to evaluate from nodes to quadrature points,TransposeMode::Transpose
to apply the transpose, mapping from quadrature points to nodesemode
-EvalMode::None
to use values directly,EvalMode::Interp
to use interpolated values,EvalMode::Grad
to use gradients,EvalMode::Weight
to use quadrature weightsu
- Input Vectorv
- Output Vector
const Q: usize = 6;
let bu = ceed.basis_tensor_H1_Lagrange(1, 1, Q, Q, QuadMode::GaussLobatto)?;
let bx = ceed.basis_tensor_H1_Lagrange(1, 1, 2, Q, QuadMode::Gauss)?;
let x_corners = ceed.vector_from_slice(&[-1., 1.])?;
let mut x_qpts = ceed.vector(Q)?;
let mut x_nodes = ceed.vector(Q)?;
bx.apply(
1,
TransposeMode::NoTranspose,
EvalMode::Interp,
&x_corners,
&mut x_nodes,
)?;
bu.apply(
1,
TransposeMode::NoTranspose,
EvalMode::Interp,
&x_nodes,
&mut x_qpts,
)?;
// Create function x^3 + 1 on Gauss Lobatto points
let mut u_arr = [0.; Q];
u_arr
.iter_mut()
.zip(x_nodes.view()?.iter())
.for_each(|(u, x)| *u = x * x * x + 1.);
let u = ceed.vector_from_slice(&u_arr)?;
// Map function to Gauss points
let mut v = ceed.vector(Q)?;
v.set_value(0.);
bu.apply(1, TransposeMode::NoTranspose, EvalMode::Interp, &u, &mut v)?;
// Verify results
v.view()?
.iter()
.zip(x_qpts.view()?.iter())
.for_each(|(v, x)| {
let true_value = x * x * x + 1.;
assert!(
(*v - true_value).abs() < 10.0 * libceed::EPSILON,
"Incorrect basis application"
);
});
Sourcepub fn dimension(&self) -> usize
pub fn dimension(&self) -> usize
Returns the dimension for given Basis
let dim = 2;
let b = ceed.basis_tensor_H1_Lagrange(dim, 1, 3, 4, QuadMode::Gauss)?;
let d = b.dimension();
assert_eq!(d, dim, "Incorrect dimension");
Sourcepub fn num_components(&self) -> usize
pub fn num_components(&self) -> usize
Returns number of components for given Basis
let ncomp = 2;
let b = ceed.basis_tensor_H1_Lagrange(1, ncomp, 3, 4, QuadMode::Gauss)?;
let n = b.num_components();
assert_eq!(n, ncomp, "Incorrect number of components");
Sourcepub fn num_nodes(&self) -> usize
pub fn num_nodes(&self) -> usize
Returns total number of nodes (in dim dimensions) of a Basis
let p = 3;
let b = ceed.basis_tensor_H1_Lagrange(2, 1, p, 4, QuadMode::Gauss)?;
let nnodes = b.num_nodes();
assert_eq!(nnodes, p * p, "Incorrect number of nodes");
Sourcepub fn num_quadrature_points(&self) -> usize
pub fn num_quadrature_points(&self) -> usize
Returns total number of quadrature points (in dim dimensions) of a Basis
let q = 4;
let b = ceed.basis_tensor_H1_Lagrange(2, 1, 3, q, QuadMode::Gauss)?;
let nqpts = b.num_quadrature_points();
assert_eq!(nqpts, q * q, "Incorrect number of quadrature points");
Sourcepub fn create_projection(&self, to: &Self) -> Result<Self>
pub fn create_projection(&self, to: &Self) -> Result<Self>
Create projection from self to specified Basis.
Both bases must have the same quadrature space. The input bases need not
be nested as function spaces; this interface solves a least squares
problem to find a representation in the to
basis that agrees at
quadrature points with the origin basis. Since the bases need not be
Lagrange, the resulting projection “basis” will have empty quadrature
points and weights.
let coarse = ceed.basis_tensor_H1_Lagrange(1, 1, 2, 3, QuadMode::Gauss)?;
let fine = ceed.basis_tensor_H1_Lagrange(1, 1, 3, 3, QuadMode::Gauss)?;
let proj = coarse.create_projection(&fine)?;
let u = ceed.vector_from_slice(&[1., 2.])?;
let mut v = ceed.vector(3)?;
proj.apply(1, TransposeMode::NoTranspose, EvalMode::Interp, &u, &mut v)?;
let expected = [1., 1.5, 2.];
for (a, b) in v.view()?.iter().zip(expected) {
assert!(
(a - b).abs() < 10.0 * libceed::EPSILON,
"Incorrect projection of linear Lagrange to quadratic Lagrange"
);
}
Trait Implementations§
Auto Trait Implementations§
impl<'a> Freeze for Basis<'a>
impl<'a> RefUnwindSafe for Basis<'a>
impl<'a> !Send for Basis<'a>
impl<'a> !Sync for Basis<'a>
impl<'a> Unpin for Basis<'a>
impl<'a> UnwindSafe for Basis<'a>
Blanket Implementations§
Source§impl<T> BorrowMut<T> for Twhere
T: ?Sized,
impl<T> BorrowMut<T> for Twhere
T: ?Sized,
Source§fn borrow_mut(&mut self) -> &mut T
fn borrow_mut(&mut self) -> &mut T
Mutably borrows from an owned value. Read more