Skip to content
This repository was archived by the owner on Dec 2, 2025. It is now read-only.
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
14 changes: 7 additions & 7 deletions src/bindings.rs
Original file line number Diff line number Diff line change
Expand Up @@ -254,7 +254,7 @@ pub mod ciarlet {
ciarlet::CiarletElement,
map::{ContravariantPiolaMap, CovariantPiolaMap, IdentityMap},
reference_cell,
traits::{ElementFamily, FiniteElement, Map},
traits::{ElementFamily, FiniteElement, Map, MappedFiniteElement},
types::{Continuity, ReferenceCellType},
};
use c_api_tools::{DType, DTypeIdentifier, cfuncs, concretise_types};
Expand Down Expand Up @@ -495,7 +495,7 @@ pub mod ciarlet {
gen_type(name = "maptype", replace_with = ["IdentityMap", "CovariantPiolaMap", "ContravariantPiolaMap"]),
field(arg = 0, name = "element", wrapper = "CiarletElementT", replace_with = ["CiarletElement<{{dtype}}, {{maptype}}>"])
)]
pub fn ciarlet_element_physical_value_size<E: FiniteElement>(
pub fn ciarlet_element_physical_value_size<E: MappedFiniteElement>(
element: &E,
gdim: usize,
) -> usize {
Expand All @@ -507,7 +507,7 @@ pub mod ciarlet {
gen_type(name = "maptype", replace_with = ["IdentityMap", "CovariantPiolaMap", "ContravariantPiolaMap"]),
field(arg = 0, name = "element", wrapper = "CiarletElementT", replace_with = ["CiarletElement<{{dtype}}, {{maptype}}>"])
)]
pub fn ciarlet_element_physical_value_rank<E: FiniteElement>(
pub fn ciarlet_element_physical_value_rank<E: MappedFiniteElement>(
element: &E,
gdim: usize,
) -> usize {
Expand All @@ -519,7 +519,7 @@ pub mod ciarlet {
gen_type(name = "maptype", replace_with = ["IdentityMap", "CovariantPiolaMap", "ContravariantPiolaMap"]),
field(arg = 0, name = "element", wrapper = "CiarletElementT", replace_with = ["CiarletElement<{{dtype}}, {{maptype}}>"])
)]
pub fn ciarlet_element_physical_value_shape<E: FiniteElement>(
pub fn ciarlet_element_physical_value_shape<E: MappedFiniteElement>(
element: &E,
gdim: usize,
shape: *mut usize,
Expand All @@ -537,7 +537,7 @@ pub mod ciarlet {
gen_type(name = "maptype", replace_with = ["IdentityMap", "CovariantPiolaMap", "ContravariantPiolaMap"]),
field(arg = 0, name = "element", wrapper = "CiarletElementT", replace_with = ["CiarletElement<{{dtype}}, {{maptype}}>"])
)]
pub fn ciarlet_element_push_forward<E: FiniteElement<CellType = ReferenceCellType>>(
pub fn ciarlet_element_push_forward<E: MappedFiniteElement<CellType = ReferenceCellType>>(
element: &E,
npoints: usize,
nfunctions: usize,
Expand Down Expand Up @@ -606,7 +606,7 @@ pub mod ciarlet {
gen_type(name = "maptype", replace_with = ["IdentityMap", "CovariantPiolaMap", "ContravariantPiolaMap"]),
field(arg = 0, name = "element", wrapper = "CiarletElementT", replace_with = ["CiarletElement<{{dtype}}, {{maptype}}>"])
)]
pub fn ciarlet_element_pull_back<E: FiniteElement<CellType = ReferenceCellType>>(
pub fn ciarlet_element_pull_back<E: MappedFiniteElement<CellType = ReferenceCellType>>(
element: &E,
npoints: usize,
nfunctions: usize,
Expand Down Expand Up @@ -685,7 +685,7 @@ pub mod ciarlet {
gen_type(name = "maptype", replace_with = ["IdentityMap", "CovariantPiolaMap", "ContravariantPiolaMap"]),
field(arg = 0, name = "element", wrapper = "CiarletElementT", replace_with = ["CiarletElement<{{dtype}}, {{maptype}}>"])
)]
pub fn ciarlet_element_embedded_superdegree<E: FiniteElement>(element: &E) -> usize {
pub fn ciarlet_element_embedded_superdegree<E: MappedFiniteElement>(element: &E) -> usize {
element.lagrange_superdegree()
}

Expand Down
38 changes: 27 additions & 11 deletions src/ciarlet.rs
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,7 @@ extern crate lapack_src;
use crate::math;
use crate::polynomials::{legendre_shape, polynomial_count, tabulate_legendre_polynomials};
use crate::reference_cell;
use crate::traits::{FiniteElement, Map};
use crate::traits::{FiniteElement, Map, MappedFiniteElement};
use crate::types::{Continuity, DofTransformation, ReferenceCellType, Transformation};
use itertools::izip;
use num::{One, Zero};
Expand Down Expand Up @@ -621,26 +621,27 @@ where
&self.interpolation_weights
}
}

impl<T: RlstScalar, M: Map> FiniteElement for CiarletElement<T, M> {
type CellType = ReferenceCellType;
type TransformationType = Transformation;
type T = T;

fn value_shape(&self) -> &[usize] {
&self.value_shape
}

fn value_size(&self) -> usize {
self.value_size
}

fn cell_type(&self) -> ReferenceCellType {
self.cell_type
}
fn lagrange_superdegree(&self) -> usize {
self.embedded_superdegree
}

fn dim(&self) -> usize {
self.dim
}

fn tabulate<
Array2Impl: ValueArrayImpl<<Self::T as RlstScalar>::Real, 2>,
Array4MutImpl: MutableArrayImpl<Self::T, 4>,
Expand Down Expand Up @@ -681,27 +682,39 @@ impl<T: RlstScalar, M: Map> FiniteElement for CiarletElement<T, M> {
}
}
}

fn tabulate_array_shape(&self, nderivs: usize, npoints: usize) -> [usize; 4] {
let deriv_count = compute_derivative_count(nderivs, self.cell_type());
let point_count = npoints;
let basis_count = self.dim();
let value_size = self.value_size();
[deriv_count, point_count, basis_count, value_size]
}

fn entity_dofs(&self, entity_dim: usize, entity_number: usize) -> Option<&[usize]> {
if entity_dim < 4 && entity_number < self.entity_dofs[entity_dim].len() {
Some(&self.entity_dofs[entity_dim][entity_number])
} else {
None
}
}

fn entity_closure_dofs(&self, entity_dim: usize, entity_number: usize) -> Option<&[usize]> {
if entity_dim < 4 && entity_number < self.entity_closure_dofs[entity_dim].len() {
Some(&self.entity_closure_dofs[entity_dim][entity_number])
} else {
None
}
}
fn tabulate_array_shape(&self, nderivs: usize, npoints: usize) -> [usize; 4] {
let deriv_count = compute_derivative_count(nderivs, self.cell_type());
let point_count = npoints;
let basis_count = self.dim();
let value_size = self.value_size();
[deriv_count, point_count, basis_count, value_size]
}

impl<T: RlstScalar, M: Map> MappedFiniteElement for CiarletElement<T, M> {
type TransformationType = Transformation;

fn lagrange_superdegree(&self) -> usize {
self.embedded_superdegree
}

fn push_forward<
Array3RealImpl: ValueArrayImpl<<Self::T as RlstScalar>::Real, 3>,
Array4Impl: ValueArrayImpl<Self::T, 4>,
Expand All @@ -724,6 +737,7 @@ impl<T: RlstScalar, M: Map> FiniteElement for CiarletElement<T, M> {
physical_values,
)
}

fn pull_back<
Array3RealImpl: ValueArrayImpl<<Self::T as RlstScalar>::Real, 3>,
Array4Impl: ValueArrayImpl<Self::T, 4>,
Expand All @@ -746,9 +760,11 @@ impl<T: RlstScalar, M: Map> FiniteElement for CiarletElement<T, M> {
reference_values,
)
}

fn physical_value_shape(&self, gdim: usize) -> Vec<usize> {
self.map.physical_value_shape(gdim)
}

fn dof_transformation(
&self,
entity: ReferenceCellType,
Expand Down
37 changes: 20 additions & 17 deletions src/traits.rs
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@ use rlst::{Array, MutableArrayImpl, RlstScalar, ValueArrayImpl};
use std::fmt::Debug;
use std::hash::Hash;

/// This trait provides the definition of a finite element.
/// A finite element.
pub trait FiniteElement {
/// The scalar type
type T: RlstScalar;
Expand All @@ -14,21 +14,9 @@ pub trait FiniteElement {
/// The cell type is typically defined through [ReferenceCellType](crate::types::ReferenceCellType).
type CellType: Debug + PartialEq + Eq + Clone + Copy + Hash;

/// Transformation type
///
/// The Transformation type specifies possible transformations of the dofs on the reference element.
/// In most cases these will be rotations and reflections as defined in [Transformation](crate::types::Transformation).
type TransformationType: Debug + PartialEq + Eq + Clone + Copy + Hash;

/// The reference cell type, eg one of `Point`, `Interval`, `Triangle`, etc.
fn cell_type(&self) -> Self::CellType;

/// The smallest degree Lagrange space that contains all possible polynomials of the finite element's polynomial space.
///
/// Details on the definition of the degree of Lagrange spaces of finite elements are
/// given [here](https://defelement.org/ciarlet.html#The+degree+of+a+finite+element).
fn lagrange_superdegree(&self) -> usize;

/// The number of basis functions.
fn dim(&self) -> usize;

Expand Down Expand Up @@ -75,6 +63,9 @@ pub trait FiniteElement {
data: &mut Array<Array4MutImpl, 4>,
);

/// Get the required shape for a tabulation array.
fn tabulate_array_shape(&self, nderivs: usize, npoints: usize) -> [usize; 4];

/// Return the dof indices that are associated with the subentity with index `entity_number` and dimension `entity_dim`.
///
/// - For `entity_dim = 0` this returns the degrees of freedom (dofs) associated with the corresponding point.
Expand All @@ -92,9 +83,21 @@ pub trait FiniteElement {
/// associated with the boundary of an entity. For an edge (for example) it returns the dofs associated
/// with the vertices at the boundary of the edge (as well as the dofs associated with the edge itself).
fn entity_closure_dofs(&self, entity_dim: usize, entity_number: usize) -> Option<&[usize]>;
}

/// Get the required shape for a tabulation array.
fn tabulate_array_shape(&self, nderivs: usize, npoints: usize) -> [usize; 4];
/// A finite element that is mapped from a reference cell.
pub trait MappedFiniteElement: FiniteElement {
/// Transformation type
///
/// The Transformation type specifies possible transformations of the dofs on the reference element.
/// In most cases these will be rotations and reflections as defined in [Transformation](crate::types::Transformation).
type TransformationType: Debug + PartialEq + Eq + Clone + Copy + Hash;

/// The smallest degree Lagrange space that contains all possible polynomials of the finite element's polynomial space.
///
/// Details on the definition of the degree of Lagrange spaces of finite elements are
/// given [here](https://defelement.org/ciarlet.html#The+degree+of+a+finite+element).
fn lagrange_superdegree(&self) -> usize;

/// Push function values forward to a physical cell.
///
Expand All @@ -117,7 +120,7 @@ pub trait FiniteElement {
/// is the topological dimension, and the third dimension is the geometric dimension. If the Jacobian is rectangular then the
/// inverse Jacobian is the pseudo-inverse of the Jacobian, ie the matrix $J^\dagger$ such that $J^\dagger J = I$.
/// - `physical_values`: The output array of the push operation. This shape of this array is the same as the `reference_values`
/// input, with the [FiniteElement::physical_value_size] used instead of the reference value size.
/// input, with the [MappedFiniteElement::physical_value_size] used instead of the reference value size.
fn push_forward<
Array3RealImpl: ValueArrayImpl<<Self::T as RlstScalar>::Real, 3>,
Array4Impl: ValueArrayImpl<Self::T, 4>,
Expand All @@ -134,7 +137,7 @@ pub trait FiniteElement {

/// Pull function values back to the reference cell.
///
/// This is the inverse operation to [FiniteElement::push_forward].
/// This is the inverse operation to [MappedFiniteElement::push_forward].
fn pull_back<
Array3RealImpl: ValueArrayImpl<<Self::T as RlstScalar>::Real, 3>,
Array4Impl: ValueArrayImpl<Self::T, 4>,
Expand Down