From 382bc568e4eda16064c9b397f799640ed57ae5e2 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Wed, 3 Dec 2025 07:34:05 +0000 Subject: [PATCH 1/7] tidy docs and patch crates.io for local nd libraries --- Cargo.toml | 4 ++++ README.md | 19 ++++++++++++++++--- ndelement/README.md | 3 ++- ndelement/src/lib.rs | 2 +- ndfunctionspace/Cargo.toml | 37 +++++++++++++++++++++++++++++++++++++ ndfunctionspace/LICENSE | 29 +++++++++++++++++++++++++++++ ndfunctionspace/README.md | 31 +++++++++++++++++++++++++++++++ ndfunctionspace/src/lib.rs | 7 +++++++ ndgrid/Cargo.toml | 7 ++----- ndgrid/README.md | 2 +- ndgrid/src/lib.rs | 2 +- 11 files changed, 131 insertions(+), 12 deletions(-) create mode 100644 ndfunctionspace/Cargo.toml create mode 100644 ndfunctionspace/LICENSE create mode 100644 ndfunctionspace/README.md create mode 100644 ndfunctionspace/src/lib.rs diff --git a/Cargo.toml b/Cargo.toml index 1d9d20c..03b260b 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -2,8 +2,12 @@ members = [ "ndelement", "ndgrid", + "ndfunctionspace", ] resolver = "2" [patch.crates-io] scotch = { git = "https://github.com/mscroggs/scotch-rs.git", branch = "mscroggs/runtime", optional = true } +ndelement = { path = "ndelement" } +ndgrid = { path = "ndgrid" } +ndfunctionspace = { path = "ndfunctionspace" } diff --git a/README.md b/README.md index 82d35dc..dff7b36 100644 --- a/README.md +++ b/README.md @@ -12,7 +12,8 @@ method. [![docs.rs](https://img.shields.io/docsrs/ndelement?label=docs.rs)](https://docs.rs/ndelement/latest/ndelement/) [![PyPI](https://img.shields.io/pypi/v/ndelement?color=blue&label=PyPI&logo=pypi&logoColor=white)](https://pypi.org/project/ndelement/) -ndelement is an open-source library written in Rust that can be used to create finite elements on 1D, 2D, or 3D reference cells. +ndelement is an open-source Rust library for defining and tabulating finite elements on 1D, 2D, or +3D reference cells. #### Using ndelement ##### Rust @@ -37,15 +38,27 @@ python -m pytest ndelement/python/test ### [ndgrid](ndgrid/) [![crates.io](https://img.shields.io/crates/v/ndgrid?color=blue)](https://crates.io/crates/ndgrid) -ndgrid is an open-source library written in Rust for handling finite element grids/meshes. +ndgrid is an open-source Rust library for handling finite element grids/meshes. -### Using ndgrid +#### Using ndgrid You can use the latest release of ndgrid by adding the following to `[dependencies]` section of your Cargo.toml file: ```toml ndgrid = "0.1.5" ``` +### [ndfunctionspace](ndfunctionspace/) +[![crates.io](https://img.shields.io/crates/v/ndfunctionspace?color=blue)](https://crates.io/crates/ndfunctionspace) + +ndfunctionspace is an open-source Rust library for creating finite element DOF maps and function spaces. + +#### Using ndfunctionspace +You can use the latest version of ndfunctionspace by adding the following to `[dependencies]` section of your Cargo.toml file: + +```toml +ndfunctionspace = { git = "https://github.com/bempp/nd.git" } +``` + ## Documentation The latest documentation of the crates in this repo is available at [bempp.github.io/nd/](https://bempp.github.io/nd/). diff --git a/ndelement/README.md b/ndelement/README.md index 0ca17c6..9dbdadb 100644 --- a/ndelement/README.md +++ b/ndelement/README.md @@ -5,7 +5,8 @@ [![docs.rs](https://img.shields.io/docsrs/ndelement?label=docs.rs)](https://docs.rs/ndelement/latest/ndelement/) [![PyPI](https://img.shields.io/pypi/v/ndelement?color=blue&label=PyPI&logo=pypi&logoColor=white)](https://pypi.org/project/ndelement/) -ndelement is an open-source library written in Rust that can be used to create finite elements on 1D, 2D, or 3D reference cells. +ndelement is an open-source Rust library for defining and tabulating finite elements on 1D, 2D, or +3D reference cells. ## Using ndelement ### Rust diff --git a/ndelement/src/lib.rs b/ndelement/src/lib.rs index 994ddcf..80680fc 100644 --- a/ndelement/src/lib.rs +++ b/ndelement/src/lib.rs @@ -1,4 +1,4 @@ -//! A library for the definition and tabulation of finite elements in Rust. +//! A library for defining and tabulating finite elements on 1D, 2D, or 3D reference cells. //! //! `ndelement` provides definition of many frequently used low and high order finite elements //! and provides routines for the tabulation of their values. diff --git a/ndfunctionspace/Cargo.toml b/ndfunctionspace/Cargo.toml new file mode 100644 index 0000000..01439ce --- /dev/null +++ b/ndfunctionspace/Cargo.toml @@ -0,0 +1,37 @@ +[features] +strict = [] + +[package] +name = "ndfunctionspace" +version = "0.0.1-dev" +edition = "2024" +authors = [ + "Matthew Scroggs ", +] +description = "Finite element function spaces." +license = "BSD-3-Clause" +homepage = "https://github.com/bempp/nd" +repository = "https://github.com/bempp/nd" +readme = "README.md" +keywords = ["numerics"] +categories = ["mathematics", "science"] + +[lib] +name = "ndfunctionspace" +crate-type = ["lib", "cdylib"] + +[dependencies] +ndelement = "0.3.0-dev" +ndgrid = "0.1.5-dev" +rlst = { git = "https://github.com/linalg-rs/rlst.git" } + +[dev-dependencies] +approx = "0.5" +paste = "1.*" + +[package.metadata.docs.rs] +cargo-args = ["-Zunstable-options", "-Zrustdoc-scrape-examples"] +rustdoc-args = [ "--html-in-header", "../katex-header.html" ] + +[lints.clippy] +wildcard_imports = "forbid" diff --git a/ndfunctionspace/LICENSE b/ndfunctionspace/LICENSE new file mode 100644 index 0000000..3e7ead8 --- /dev/null +++ b/ndfunctionspace/LICENSE @@ -0,0 +1,29 @@ +BSD 3-Clause License + +Copyright (c) 2025-, Timo Betcke, Matthew Scroggs +All rights reserved. + +Redistribution and use in source and binary forms, with or without +modification, are permitted provided that the following conditions are met: + +1. Redistributions of source code must retain the above copyright notice, this + list of conditions and the following disclaimer. + +2. Redistributions in binary form must reproduce the above copyright notice, + this list of conditions and the following disclaimer in the documentation + and/or other materials provided with the distribution. + +3. Neither the name of the copyright holder nor the names of its + contributors may be used to endorse or promote products derived from + this software without specific prior written permission. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE +DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE +FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL +DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR +SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER +CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, +OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. diff --git a/ndfunctionspace/README.md b/ndfunctionspace/README.md new file mode 100644 index 0000000..b519c04 --- /dev/null +++ b/ndfunctionspace/README.md @@ -0,0 +1,31 @@ +# ndfunctionspace +[![crates.io](https://img.shields.io/crates/v/ndfunctionspace?color=blue)](https://crates.io/crates/ndfunctionspace) + +ndfunctionspace is an open-source Rust library for creating finite element DOF maps and function spaces. + +## Using ndfunctionspace +You can use the latest release of ndfunctionspace by adding the following to `[dependencies]` section of your Cargo.toml file: + +```toml +ndfunctionspace = { git = "https://github.com/bempp/nd.git" } +``` + +## Documentation +The latest documentation of nfunctionspace is available at [bempp.github.io/nd/rust/ndfunctionspace](https://bempp.github.io/nd/rust/ndfunctionspace/). + +## Testing +The Rust functionality of the library can be tested by running: +```bash +cargo test +``` + +## Examples +Examples of use can be found in the [examples](examples/) folder. + +## Getting help +Errors in should be added to the [nd GitHub issue tracker](https://github.com/bempp/nd/issues). + +Questions about use can be asked on the [Bempp Discourse](https://bempp.discourse.group). + +## Licence +ndfunctionspace is licensed under a BSD 3-Clause licence. diff --git a/ndfunctionspace/src/lib.rs b/ndfunctionspace/src/lib.rs new file mode 100644 index 0000000..f432716 --- /dev/null +++ b/ndfunctionspace/src/lib.rs @@ -0,0 +1,7 @@ +//! A library for creating finite element DOF maps and function spaces. + +#![cfg_attr(feature = "strict", deny(warnings), deny(unused_crate_dependencies))] +#![warn(missing_docs)] + +use ndgrid as _; +use ndelement as _; diff --git a/ndgrid/Cargo.toml b/ndgrid/Cargo.toml index 31ef714..08f8734 100644 --- a/ndgrid/Cargo.toml +++ b/ndgrid/Cargo.toml @@ -25,13 +25,13 @@ categories = ["mathematics", "science"] [lib] name = "ndgrid" -crate-type = ["lib", "cdylib"] +crate-type = ["lib"] [dependencies] coupe = { git = "https://github.com/LIHPC-Computational-Geometry/coupe.git", optional = true } itertools = "0.14.*" mpi = { git = "https://github.com/rsmpi/rsmpi.git", optional = true } -ndelement = { path = "../ndelement" } +ndelement = "0.3.0-dev" num = "0.4" rlst = { git = "https://github.com/linalg-rs/rlst.git" } serde = { version = "1", features = ["derive"], optional = true } @@ -42,9 +42,6 @@ scotch = { version = "0.2.1", optional = true } approx = "0.5" paste = "1.*" -[build-dependencies] -cbindgen = "0.29.2" - [package.metadata.docs.rs] cargo-args = ["-Zunstable-options", "-Zrustdoc-scrape-examples"] rustdoc-args = [ "--html-in-header", "../katex-header.html" ] diff --git a/ndgrid/README.md b/ndgrid/README.md index 4948167..94d1696 100644 --- a/ndgrid/README.md +++ b/ndgrid/README.md @@ -1,7 +1,7 @@ # ndgrid [![crates.io](https://img.shields.io/crates/v/ndgrid?color=blue)](https://crates.io/crates/ndgrid) -ndgrid is an open-source library written in Rust for handling finite element grids/meshes. +ndgrid is an open-source Rust library for handling finite element grids/meshes. ## Using ndgrid You can use the latest release of ndgrid by adding the following to `[dependencies]` section of your Cargo.toml file: diff --git a/ndgrid/src/lib.rs b/ndgrid/src/lib.rs index f58418a..14213c7 100644 --- a/ndgrid/src/lib.rs +++ b/ndgrid/src/lib.rs @@ -1,4 +1,4 @@ -//! Data structures to work with n-dimensional grids. +//! A library for handling Finite element grids/meshes. //! //! `ndgrid` builds upen [ndelement] to provide data structures for grids on either a single node or distributed via MPI. //! From 1262df37eaf0a815c8fea7079a4db8b37c9f1939 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Wed, 3 Dec 2025 07:36:31 +0000 Subject: [PATCH 2/7] remove badge that doesn't exist yet --- README.md | 1 - ndfunctionspace/README.md | 1 - 2 files changed, 2 deletions(-) diff --git a/README.md b/README.md index dff7b36..d0c45ba 100644 --- a/README.md +++ b/README.md @@ -48,7 +48,6 @@ ndgrid = "0.1.5" ``` ### [ndfunctionspace](ndfunctionspace/) -[![crates.io](https://img.shields.io/crates/v/ndfunctionspace?color=blue)](https://crates.io/crates/ndfunctionspace) ndfunctionspace is an open-source Rust library for creating finite element DOF maps and function spaces. diff --git a/ndfunctionspace/README.md b/ndfunctionspace/README.md index b519c04..a6c6535 100644 --- a/ndfunctionspace/README.md +++ b/ndfunctionspace/README.md @@ -1,5 +1,4 @@ # ndfunctionspace -[![crates.io](https://img.shields.io/crates/v/ndfunctionspace?color=blue)](https://crates.io/crates/ndfunctionspace) ndfunctionspace is an open-source Rust library for creating finite element DOF maps and function spaces. From 453a25992c1856b6f3b93ddb2c7044754fea72af Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Mon, 5 Jan 2026 11:11:25 +0000 Subject: [PATCH 3/7] start work on functionspace --- ndfunctionspace/src/function.rs | 445 ++++++++++++++++++++++++++ ndfunctionspace/src/function_space.rs | 59 ++++ ndfunctionspace/src/lib.rs | 6 +- ndfunctionspace/src/traits.rs | 39 +++ 4 files changed, 548 insertions(+), 1 deletion(-) create mode 100644 ndfunctionspace/src/function.rs create mode 100644 ndfunctionspace/src/function_space.rs create mode 100644 ndfunctionspace/src/traits.rs diff --git a/ndfunctionspace/src/function.rs b/ndfunctionspace/src/function.rs new file mode 100644 index 0000000..c725505 --- /dev/null +++ b/ndfunctionspace/src/function.rs @@ -0,0 +1,445 @@ +//! Functions and function spaces + +//mod function_space; + +use mpi::request::WaitGuard; +use mpi::traits::{Communicator, Destination, Source}; +use ndelement::ciarlet::CiarletElement; +use ndelement::traits::ElementFamily; +use ndelement::{traits::FiniteElement, types::ReferenceCellType}; +use ndgrid::traits::ParallelGrid; +use ndgrid::traits::{Entity, Topology}; +use ndgrid::{traits::Grid, types::Ownership}; +use rlst::{MatrixInverse, RlstScalar}; +use std::collections::HashMap; +use std::marker::PhantomData; + +type DofList = Vec>; +type OwnerData = Vec<(usize, usize, usize, usize)>; + +/// A function space +pub trait FunctionSpaceTrait { + /// Communicator + type C: Communicator; + /// Scalar type + type T: RlstScalar; + /// The grid type + type Grid: Grid::Real, EntityDescriptor = ReferenceCellType> + + ParallelGrid; + /// The finite element type + type FiniteElement: FiniteElement + Sync; + + /// Get the communicator + fn comm(&self) -> &Self::C; + + /// Get the grid that the element is defined on + fn grid(&self) -> &Self::Grid; + + /// Get the finite element used to define this function space + fn element(&self, cell_type: ReferenceCellType) -> &Self::FiniteElement; + + /// Check if the function space is stored in serial + fn is_serial(&self) -> bool; + + /// Get the DOF numbers on the local process associated with the given entity + fn get_local_dof_numbers(&self, entity_dim: usize, entity_number: usize) -> &[usize]; + + /// Get the number of DOFs associated with the local process + fn local_size(&self) -> usize; + + /// Get the number of DOFs on all processes + fn global_size(&self) -> usize; + + /// Get the local DOF numbers associated with a cell + fn cell_dofs(&self, cell: usize) -> Option<&[usize]>; + + /// Get the local DOF numbers associated with a cell + /// + /// # Safety + /// The function uses unchecked array access + unsafe fn cell_dofs_unchecked(&self, cell: usize) -> &[usize]; + + /// Compute a colouring of the cells so that no two cells that share an entity with DOFs associated with it are assigned the same colour + fn cell_colouring(&self) -> HashMap>>; + + /// Get the global DOF index associated with a local DOF index + fn global_dof_index(&self, local_dof_index: usize) -> usize; + + /// Get ownership of a local DOF + fn ownership(&self, local_dof_index: usize) -> Ownership; +} + +/// Implementation of a general function space. +pub struct FunctionSpace< + 'a, + C: Communicator, + T: RlstScalar + MatrixInverse, + GridImpl: ParallelGrid + Grid, +> { + grid: &'a GridImpl, + elements: HashMap>, + entity_dofs: [Vec>; 4], + cell_dofs: Vec>, + local_size: usize, + global_size: usize, + global_dof_numbers: Vec, + ownership: Vec, + _marker: std::marker::PhantomData, +} + +unsafe impl< + C: Communicator, + T: RlstScalar + MatrixInverse, + GridImpl: ParallelGrid + Grid, + > Sync for FunctionSpace<'_, C, T, GridImpl> +{ +} + +impl< + 'a, + C: Communicator, + T: RlstScalar + MatrixInverse, + GridImpl: ParallelGrid + Grid, + > FunctionSpace<'a, C, T, GridImpl> +{ + /// Create new function space + pub fn new( + grid: &'a GridImpl, + e_family: &impl ElementFamily< + T = T, + FiniteElement = CiarletElement, + CellType = ReferenceCellType, + >, + ) -> Self { + let comm = grid.comm(); + let rank = comm.rank(); + let size = comm.size(); + + // Create local space on current process + let (cell_dofs, entity_dofs, dofmap_size, owner_data) = + assign_dofs(rank as usize, grid.local_grid(), e_family); + + let mut elements = HashMap::new(); + for cell in grid.entity_types(grid.topology_dim()) { + elements.insert(*cell, e_family.element(*cell)); + } + + // Assign global DOF numbers + let mut global_dof_numbers = vec![0; dofmap_size]; + let mut ghost_indices = vec![vec![]; size as usize]; + let mut ghost_dims = vec![vec![]; size as usize]; + let mut ghost_entities = vec![vec![]; size as usize]; + let mut ghost_entity_dofs = vec![vec![]; size as usize]; + + let local_offset = if rank == 0 { + 0 + } else { + let (value, _status) = comm.process_at_rank(rank - 1).receive::(); + value + }; + let mut dof_n = local_offset; + for (i, ownership) in owner_data.iter().enumerate() { + if ownership.0 == rank as usize { + global_dof_numbers[i] = dof_n; + dof_n += 1; + } else { + ghost_indices[ownership.0].push(i); + ghost_dims[ownership.0].push(ownership.1); + ghost_entities[ownership.0].push(ownership.2); + ghost_entity_dofs[ownership.0].push(ownership.3); + } + } + if rank < size - 1 { + mpi::request::scope(|scope| { + let _ = + WaitGuard::from(comm.process_at_rank(rank + 1).immediate_send(scope, &dof_n)); + }); + } + + let global_size = if rank == size - 1 { + for p in 0..rank { + mpi::request::scope(|scope| { + let _ = WaitGuard::from(comm.process_at_rank(p).immediate_send(scope, &dof_n)); + }); + } + dof_n + } else { + let (gs, _status) = comm.process_at_rank(size - 1).receive::(); + gs + }; + + // Communicate information about ghosts + // send requests for ghost info + for p in 0..size { + if p != rank { + mpi::request::scope(|scope| { + let process = comm.process_at_rank(p); + let _ = WaitGuard::from(process.immediate_send(scope, &ghost_dims[p as usize])); + let _ = + WaitGuard::from(process.immediate_send(scope, &ghost_entities[p as usize])); + let _ = WaitGuard::from( + process.immediate_send(scope, &ghost_entity_dofs[p as usize]), + ); + }); + } + } + // accept requests and send ghost info + for p in 0..size { + if p != rank { + let process = comm.process_at_rank(p); + let (gdims, _status) = process.receive_vec::(); + let (gentities, _status) = process.receive_vec::(); + let (gentity_dofs, _status) = process.receive_vec::(); + let local_ghost_dofs = gdims + .iter() + .zip(gentities.iter().zip(&gentity_dofs)) + .map(|(c, (e, d))| entity_dofs[*c][*e][*d]) + .collect::>(); + let global_ghost_dofs = local_ghost_dofs + .iter() + .map(|i| global_dof_numbers[*i]) + .collect::>(); + mpi::request::scope(|scope| { + let _ = WaitGuard::from(process.immediate_send(scope, &local_ghost_dofs)); + let _ = WaitGuard::from(process.immediate_send(scope, &global_ghost_dofs)); + }); + } + } + + // receive ghost info + let mut ownership = vec![Ownership::Owned; dofmap_size]; + for p in 0..size { + if p != rank { + let process = comm.process_at_rank(p); + let (local_ghost_dofs, _status) = process.receive_vec::(); + let (global_ghost_dofs, _status) = process.receive_vec::(); + for (i, (l, g)) in ghost_indices[p as usize] + .iter() + .zip(local_ghost_dofs.iter().zip(&global_ghost_dofs)) + { + global_dof_numbers[*i] = *g; + ownership[*i] = Ownership::Ghost(p as usize, *l); + } + } + } + + Self { + grid, + elements, + entity_dofs, + cell_dofs, + local_size: dofmap_size, + global_size, + global_dof_numbers, + ownership, + _marker: PhantomData, + } + } +} + +impl< + C: Communicator, + T: RlstScalar + MatrixInverse, + GridImpl: ParallelGrid + Grid, + > FunctionSpaceTrait for FunctionSpace<'_, C, T, GridImpl> +{ + type T = T; + type Grid = GridImpl; + type FiniteElement = CiarletElement; + type C = C; + + fn grid(&self) -> &Self::Grid { + self.grid + } + fn element(&self, cell_type: ReferenceCellType) -> &CiarletElement { + &self.elements[&cell_type] + } + + fn get_local_dof_numbers(&self, entity_dim: usize, entity_number: usize) -> &[usize] { + &self.entity_dofs[entity_dim][entity_number] + } + fn is_serial(&self) -> bool { + self.grid.comm().size() == 1 + } + fn local_size(&self) -> usize { + self.local_size + } + + fn global_size(&self) -> usize { + self.global_size + } + unsafe fn cell_dofs_unchecked(&self, cell: usize) -> &[usize] { + self.cell_dofs.get_unchecked(cell) + } + fn cell_dofs(&self, cell: usize) -> Option<&[usize]> { + if cell < self.cell_dofs.len() { + Some(unsafe { self.cell_dofs_unchecked(cell) }) + } else { + None + } + } + fn cell_colouring(&self) -> HashMap>> { + let mut colouring = HashMap::new(); + //: HashMap>> + for cell in self.grid.entity_types(2) { + colouring.insert(*cell, vec![]); + } + let mut edim = 0; + while self.elements[&self.grid.entity_types(2)[0]] + .entity_dofs(edim, 0) + .unwrap() + .is_empty() + { + edim += 1; + } + + let mut entity_colours = vec![ + vec![]; + if edim == 0 { + self.grid.entity_count(ReferenceCellType::Point) + } else if edim == 1 { + self.grid.entity_count(ReferenceCellType::Interval) + } else if edim == 2 && self.grid.topology_dim() == 2 { + self.grid + .entity_types(2) + .iter() + .map(|&i| self.grid.entity_count(i)) + .sum::() + } else { + unimplemented!(); + } + ]; + + for cell in self.grid.entity_iter(2) { + let cell_type = cell.entity_type(); + let indices = cell.topology().sub_entity_iter(edim).collect::>(); + + let c = { + let mut c = 0; + while c < colouring[&cell_type].len() { + let mut found = false; + for v in &indices { + if entity_colours[*v].contains(&c) { + found = true; + break; + } + } + + if !found { + break; + } + c += 1; + } + c + }; + if c == colouring[&cell_type].len() { + for ct in self.grid.entity_types(2) { + colouring.get_mut(ct).unwrap().push(if *ct == cell_type { + vec![cell.local_index()] + } else { + vec![] + }); + } + } else { + colouring.get_mut(&cell_type).unwrap()[c].push(cell.local_index()); + } + for v in &indices { + entity_colours[*v].push(c); + } + } + colouring + } + fn global_dof_index(&self, local_dof_index: usize) -> usize { + self.global_dof_numbers[local_dof_index] + } + fn ownership(&self, local_dof_index: usize) -> Ownership { + self.ownership[local_dof_index] + } + + fn comm(&self) -> &C { + self.grid.comm() + } +} + +/// Assign DOFs to entities. +pub fn assign_dofs< + T: RlstScalar + MatrixInverse, + GridImpl: Grid + Sync, +>( + rank: usize, + grid: &GridImpl, + e_family: &impl ElementFamily< + T = T, + FiniteElement = CiarletElement, + CellType = ReferenceCellType, + >, +) -> (DofList, [DofList; 4], usize, OwnerData) { + let mut size = 0; + let mut entity_dofs: [Vec>; 4] = [vec![], vec![], vec![], vec![]]; + let mut owner_data = vec![]; + let tdim = grid.topology_dim(); + + let mut elements = HashMap::new(); + let mut element_dims = HashMap::new(); + for cell in grid.entity_types(2) { + elements.insert(*cell, e_family.element(*cell)); + element_dims.insert(*cell, elements[cell].dim()); + } + + let entity_counts = (0..tdim + 1) + .map(|d| { + grid.entity_types(d) + .iter() + .map(|&i| grid.entity_count(i)) + .sum::() + }) + .collect::>(); + if tdim > 2 { + unimplemented!("DOF maps not implemented for cells with tdim > 2."); + } + + for d in 0..tdim + 1 { + entity_dofs[d] = vec![vec![]; entity_counts[d]]; + } + let mut cell_dofs = vec![vec![]; entity_counts[tdim]]; + + let mut max_rank = rank; + for cell in grid.entity_iter(tdim) { + if let Ownership::Ghost(process, _index) = cell.ownership() { + if process > max_rank { + max_rank = process; + } + } + } + for cell in grid.entity_iter(tdim) { + cell_dofs[cell.local_index()] = vec![0; element_dims[&cell.entity_type()]]; + let element = &elements[&cell.entity_type()]; + let topology = cell.topology(); + + // Assign DOFs to entities + for (d, edofs_d) in entity_dofs.iter_mut().take(tdim + 1).enumerate() { + for (i, e) in topology.sub_entity_iter(d).enumerate() { + let e_dofs = element.entity_dofs(d, i).unwrap(); + if !e_dofs.is_empty() { + if edofs_d[e].is_empty() { + for (dof_i, _d) in e_dofs.iter().enumerate() { + edofs_d[e].push(size); + if let Ownership::Ghost(process, index) = + grid.entity(d, e).unwrap().ownership() + { + owner_data.push((process, d, index, dof_i)); + } else { + owner_data.push((rank, d, e, dof_i)); + } + size += 1; + } + } + for (local_dof, dof) in e_dofs.iter().zip(&edofs_d[e]) { + cell_dofs[cell.local_index()][*local_dof] = *dof; + } + } + } + } + } + (cell_dofs, entity_dofs, size, owner_data) +} diff --git a/ndfunctionspace/src/function_space.rs b/ndfunctionspace/src/function_space.rs new file mode 100644 index 0000000..be881ed --- /dev/null +++ b/ndfunctionspace/src/function_space.rs @@ -0,0 +1,59 @@ +//! Function space traits +use crate::traits::FunctionSpace as FunctionSpaceTrait; +use ndelement::traits::FiniteElement; +use ndgrid::{traits::Grid, types::Ownership}; +use std::collections::HashMap; + +/// Function space. +pub struct FunctionSpace<'a, G: Grid, F: FiniteElement> { + grid: &'a G, + elements: Vec, + entity_dofs: HashMap>>, +} + +impl<'a, G: Grid, F: FiniteElement> FunctionSpaceTrait for FunctionSpace<'a, G, F> { + type Grid = G; + type FiniteElement = F; + + fn grid(&self) -> &G { + self.grid + } + + fn elements(&self) -> &[F] { + &self.elements + } + + fn entities_by_element(&self, element_index: usize) -> Option<&[usize]> { + todo!(); + } + + fn entity_dofs( + &self, + entity_type: ::EntityDescriptor, + entity_number: usize, + ) -> Option<&[usize]> { + if let Some(i) = self.entity_dofs.get(&entity_type) { + if let Some(j) = i.get(entity_number) { + Some(j) + } else { + None + } + } else { None } + } + + fn local_size(&self) -> usize { + todo!(); + } + + fn global_size(&self) -> usize { + todo!(); + } + + fn global_dof_index(&self, local_dof_index: usize) -> usize { + todo!(); + } + + fn ownership(&self, local_dof_index: usize) -> Ownership { + todo!(); + } +} diff --git a/ndfunctionspace/src/lib.rs b/ndfunctionspace/src/lib.rs index f432716..98d7ad2 100644 --- a/ndfunctionspace/src/lib.rs +++ b/ndfunctionspace/src/lib.rs @@ -3,5 +3,9 @@ #![cfg_attr(feature = "strict", deny(warnings), deny(unused_crate_dependencies))] #![warn(missing_docs)] -use ndgrid as _; +mod function_space; +pub mod traits; +pub use function_space::FunctionSpace; + use ndelement as _; +use ndgrid as _; diff --git a/ndfunctionspace/src/traits.rs b/ndfunctionspace/src/traits.rs new file mode 100644 index 0000000..eaebdc1 --- /dev/null +++ b/ndfunctionspace/src/traits.rs @@ -0,0 +1,39 @@ +//! Function space traits +use ndelement::traits::FiniteElement; +use ndgrid::{traits::Grid, types::Ownership}; + +/// Function space. +pub trait FunctionSpace { + /// The type for the grid this function space is defined on + type Grid: Grid; + /// The type for the finite element this grid is defined by + type FiniteElement: FiniteElement; + + /// The grid that this function space is defined on + fn grid(&self) -> &Self::Grid; + + /// The finite elements used in this function space + fn elements(&self) -> &[Self::FiniteElement]; + + /// A list of cell indices that use the element with the given index + fn cells_by_element(&self, element_index: usize) -> Option<&[usize]>; + + /// Get the local DOFs numbers associated with the given entity + fn entity_dofs( + &self, + entity_type: ::EntityDescriptor, + entity_number: usize, + ) -> Option<&[usize]>; + + /// Get the number of DOFs associated with the local process + fn local_size(&self) -> usize; + + /// Get the number of DOFs on all processes + fn global_size(&self) -> usize; + + /// Get the global DOF index associated with a local DOF index + fn global_dof_index(&self, local_dof_index: usize) -> usize; + + /// Get ownership of a local DOF + fn ownership(&self, local_dof_index: usize) -> Ownership; +} From c427f7066d837ad862fa2f28f552f21346324b51 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Mon, 5 Jan 2026 14:58:34 +0000 Subject: [PATCH 4/7] implement serial function space --- ndelement/src/reference_cell.rs | 36 ++ ndfunctionspace/Cargo.toml | 1 + ndfunctionspace/src/function_space.rs | 596 +++++++++++++++++++++++++- ndfunctionspace/src/lib.rs | 2 +- ndfunctionspace/src/traits.rs | 21 +- 5 files changed, 633 insertions(+), 23 deletions(-) diff --git a/ndelement/src/reference_cell.rs b/ndelement/src/reference_cell.rs index 063e3b0..6626ad0 100644 --- a/ndelement/src/reference_cell.rs +++ b/ndelement/src/reference_cell.rs @@ -375,6 +375,42 @@ pub fn entity_types(cell: ReferenceCellType) -> Vec> { } } +/// The local indices of the subentities of the reference cell, with each entity type indexed separately +pub fn entity_indices(cell: ReferenceCellType) -> Vec> { + match cell { + ReferenceCellType::Point => vec![vec![0], vec![], vec![], vec![]], + ReferenceCellType::Interval => vec![vec![0, 1], vec![0], vec![], vec![]], + ReferenceCellType::Triangle => vec![vec![0, 1, 2], vec![0, 1, 2], vec![0], vec![]], + ReferenceCellType::Quadrilateral => { + vec![vec![0, 1, 2, 3], vec![0, 1, 2, 3], vec![0], vec![]] + } + ReferenceCellType::Tetrahedron => vec![ + vec![0, 1, 2, 3], + vec![0, 1, 2, 3, 4, 5], + vec![0, 1, 2, 3], + vec![0], + ], + ReferenceCellType::Hexahedron => vec![ + vec![0, 1, 2, 3, 4, 5, 6, 7], + vec![0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11], + vec![0, 1, 2, 3, 4, 5], + vec![0], + ], + ReferenceCellType::Prism => vec![ + vec![0, 1, 2, 3, 4, 5], + vec![0, 1, 2, 3, 4, 5, 6, 7, 8], + vec![0, 0, 1, 2, 1], + vec![0], + ], + ReferenceCellType::Pyramid => vec![ + vec![0, 1, 2, 3, 4], + vec![0, 1, 2, 3, 4, 5, 6, 7], + vec![0, 0, 1, 2, 3], + vec![0], + ], + } +} + /// The number of subentities of each dimension pub fn entity_counts(cell: ReferenceCellType) -> Vec { match cell { diff --git a/ndfunctionspace/Cargo.toml b/ndfunctionspace/Cargo.toml index 7695694..45cf5ec 100644 --- a/ndfunctionspace/Cargo.toml +++ b/ndfunctionspace/Cargo.toml @@ -25,6 +25,7 @@ ndelement = { path = "../ndelement" } ndgrid = { path = "../ndgrid" } num = "0.4" rlst = "0.4" +itertools = "0.14.*" [dev-dependencies] approx = "0.5" diff --git a/ndfunctionspace/src/function_space.rs b/ndfunctionspace/src/function_space.rs index be881ed..8a7263f 100644 --- a/ndfunctionspace/src/function_space.rs +++ b/ndfunctionspace/src/function_space.rs @@ -1,17 +1,118 @@ //! Function space traits use crate::traits::FunctionSpace as FunctionSpaceTrait; -use ndelement::traits::FiniteElement; -use ndgrid::{traits::Grid, types::Ownership}; +use itertools::izip; +use ndelement::{ + reference_cell, + traits::{ElementFamily, MappedFiniteElement}, + types::ReferenceCellType, +}; +use ndgrid::{ + traits::{Entity, Grid, Topology}, + types::Ownership, +}; use std::collections::HashMap; +use std::fmt::Debug; +use std::hash::Hash; /// Function space. -pub struct FunctionSpace<'a, G: Grid, F: FiniteElement> { +pub struct SerialFunctionSpace< + 'a, + E: Debug + PartialEq + Eq + Clone + Copy + Hash, + G: Grid, + F: MappedFiniteElement, +> { grid: &'a G, elements: Vec, - entity_dofs: HashMap>>, + entity_dofs: HashMap>>, + entity_closure_dofs: HashMap>>, + ndofs: usize, + entities_by_element: HashMap>, } -impl<'a, G: Grid, F: FiniteElement> FunctionSpaceTrait for FunctionSpace<'a, G, F> { +impl< + 'a, + G: Grid, + F: MappedFiniteElement, +> SerialFunctionSpace<'a, ReferenceCellType, G, F> +{ + /// Create a new serial function space + pub fn new>( + grid: &'a G, + family: &EF, + ) -> Self { + let elements = grid + .entity_types(grid.topology_dim()) + .iter() + .map(|e| family.element(*e)) + .collect::>(); + let mut entity_dofs = HashMap::new(); + let mut entity_closure_dofs = HashMap::new(); + let mut entities_by_element = HashMap::new(); + + for d in 0..=grid.topology_dim() { + for e in grid.entity_types(d) { + entity_dofs.insert(*e, vec![vec![]; grid.entity_count(*e)]); + entity_closure_dofs.insert(*e, vec![vec![]; grid.entity_count(*e)]); + } + } + + let mut ndofs = 0; + + for (cell_type, element) in izip!(grid.entity_types(grid.topology_dim()), &elements) { + let sub_entity_types = reference_cell::entity_types(*cell_type); + let sub_entity_indices = reference_cell::entity_indices(*cell_type); + let mut cell_indices = vec![]; + for cell in grid.entity_iter(*cell_type) { + cell_indices.push(cell.local_index()); + let mut cell_dofs = vec![]; // &mut entity_closure_dofs.get_mut(cell_type).unwrap()[cell.local_index()]; + for (d, (types, indices)) in + izip!(&sub_entity_types, &sub_entity_indices).enumerate() + { + for (t, i) in izip!(types, indices) { + let reference_entity_dofs = element.entity_dofs(d, *i).unwrap(); + if reference_entity_dofs.len() > 0 && !reference_entity_dofs.is_empty() { + let ed = &mut entity_dofs.get_mut(t).unwrap() + [cell.topology().sub_entity(*t, *i)]; + while ed.len() < reference_entity_dofs.len() { + ed.push(ndofs); + ndofs += 1; + } + for dof in ed { + cell_dofs.push(*dof); + } + } + entity_closure_dofs.get_mut(t).unwrap() + [cell.topology().sub_entity(*t, *i)] = element + .entity_closure_dofs(d, *i) + .unwrap() + .iter() + .map(|e| cell_dofs[*e]) + .collect::>(); + } + } + } + entities_by_element.insert(*cell_type, cell_indices); + } + + Self { + grid, + elements, + entity_dofs, + entity_closure_dofs, + ndofs, + entities_by_element, + } + } +} + +impl< + 'a, + E: Debug + PartialEq + Eq + Clone + Copy + Hash, + G: Grid, + F: MappedFiniteElement, +> FunctionSpaceTrait for SerialFunctionSpace<'a, E, G, F> +{ + type EntityDescriptor = E; type Grid = G; type FiniteElement = F; @@ -24,36 +125,497 @@ impl<'a, G: Grid, F: FiniteElement> FunctionSpaceTrait for FunctionSpace<'a, G, } fn entities_by_element(&self, element_index: usize) -> Option<&[usize]> { - todo!(); + self.entities_by_element + .get(&self.elements[element_index].cell_type()) + .map(|v| &**v) } - fn entity_dofs( - &self, - entity_type: ::EntityDescriptor, - entity_number: usize, - ) -> Option<&[usize]> { + fn entity_dofs(&self, entity_type: E, entity_number: usize) -> Option<&[usize]> { if let Some(i) = self.entity_dofs.get(&entity_type) { if let Some(j) = i.get(entity_number) { Some(j) } else { None } - } else { None } + } else { + None + } + } + + fn entity_closure_dofs(&self, entity_type: E, entity_number: usize) -> Option<&[usize]> { + if let Some(i) = self.entity_closure_dofs.get(&entity_type) { + if let Some(j) = i.get(entity_number) { + Some(j) + } else { + None + } + } else { + None + } } fn local_size(&self) -> usize { - todo!(); + self.ndofs } fn global_size(&self) -> usize { - todo!(); + self.ndofs } fn global_dof_index(&self, local_dof_index: usize) -> usize { - todo!(); + local_dof_index + } + + fn ownership(&self, _local_dof_index: usize) -> Ownership { + Ownership::Owned } +} + +#[cfg(test)] +mod test { + use super::*; + use ndelement::{ + ciarlet::LagrangeElementFamily, + types::{Continuity, ReferenceCellType}, + }; + use ndgrid::shapes::unit_cube_boundary; + + #[test] + fn test_dp0() { + let grid = unit_cube_boundary::(2, 2, 2, ReferenceCellType::Triangle); + let family = LagrangeElementFamily::::new(0, Continuity::Discontinuous); + + let space = SerialFunctionSpace::new(&grid, &family); + + assert_eq!( + space.local_size(), + grid.entity_count(ReferenceCellType::Triangle) + ); + + for cell in grid.entity_iter(ReferenceCellType::Point) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 0 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 0 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Interval) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 0 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 0 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Triangle) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Triangle, cell.local_index()) + .unwrap() + .len(), + 1 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Triangle, cell.local_index()) + .unwrap() + .len(), + 1 + ); + } + } + + #[test] + fn test_p1() { + let grid = unit_cube_boundary::(2, 2, 2, ReferenceCellType::Triangle); + let family = LagrangeElementFamily::::new(1, Continuity::Standard); + + let space = SerialFunctionSpace::new(&grid, &family); + + assert_eq!( + space.local_size(), + grid.entity_count(ReferenceCellType::Point) + ); + + for cell in grid.entity_iter(ReferenceCellType::Point) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 1 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 1 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Interval) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 0 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 2 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Triangle) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Triangle, cell.local_index()) + .unwrap() + .len(), + 0 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Triangle, cell.local_index()) + .unwrap() + .len(), + 3 + ); + } + } + + #[test] + fn test_dp1() { + let grid = unit_cube_boundary::(2, 2, 2, ReferenceCellType::Triangle); + let family = LagrangeElementFamily::::new(1, Continuity::Discontinuous); + + let space = SerialFunctionSpace::new(&grid, &family); + + assert_eq!( + space.local_size(), + 3 * grid.entity_count(ReferenceCellType::Triangle) + ); + + for cell in grid.entity_iter(ReferenceCellType::Point) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 0 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 0 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Interval) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 0 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 0 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Triangle) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Triangle, cell.local_index()) + .unwrap() + .len(), + 3 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Triangle, cell.local_index()) + .unwrap() + .len(), + 3 + ); + } + } + + #[test] + fn test_p2() { + let grid = unit_cube_boundary::(2, 2, 2, ReferenceCellType::Triangle); + let family = LagrangeElementFamily::::new(2, Continuity::Standard); + + let space = SerialFunctionSpace::new(&grid, &family); + + assert_eq!( + space.local_size(), + grid.entity_count(ReferenceCellType::Point) + + grid.entity_count(ReferenceCellType::Interval) + ); + + for cell in grid.entity_iter(ReferenceCellType::Point) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 1 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 1 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Interval) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 1 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 3 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Triangle) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Triangle, cell.local_index()) + .unwrap() + .len(), + 0 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Triangle, cell.local_index()) + .unwrap() + .len(), + 6 + ); + } + } + + #[test] + fn test_p3() { + let grid = unit_cube_boundary::(2, 2, 2, ReferenceCellType::Triangle); + let family = LagrangeElementFamily::::new(3, Continuity::Standard); + + let space = SerialFunctionSpace::new(&grid, &family); + + assert_eq!( + space.local_size(), + grid.entity_count(ReferenceCellType::Point) + + 2 * grid.entity_count(ReferenceCellType::Interval) + + grid.entity_count(ReferenceCellType::Triangle) + ); + + for cell in grid.entity_iter(ReferenceCellType::Point) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 1 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 1 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Interval) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 2 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 4 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Triangle) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Triangle, cell.local_index()) + .unwrap() + .len(), + 1 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Triangle, cell.local_index()) + .unwrap() + .len(), + 10 + ); + } + } + + #[test] + fn test_p1_quad() { + let grid = unit_cube_boundary::(2, 2, 2, ReferenceCellType::Quadrilateral); + let family = LagrangeElementFamily::::new(1, Continuity::Standard); + + let space = SerialFunctionSpace::new(&grid, &family); + + assert_eq!( + space.local_size(), + grid.entity_count(ReferenceCellType::Point) + ); + + for cell in grid.entity_iter(ReferenceCellType::Point) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 1 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 1 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Interval) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 0 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 2 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Quadrilateral) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Quadrilateral, cell.local_index()) + .unwrap() + .len(), + 0 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Quadrilateral, cell.local_index()) + .unwrap() + .len(), + 4 + ); + } + } + + #[test] + fn test_p2_quad() { + let grid = unit_cube_boundary::(2, 2, 2, ReferenceCellType::Quadrilateral); + let family = LagrangeElementFamily::::new(2, Continuity::Standard); + + let space = SerialFunctionSpace::new(&grid, &family); + + assert_eq!( + space.local_size(), + grid.entity_count(ReferenceCellType::Point) + + grid.entity_count(ReferenceCellType::Interval) + + grid.entity_count(ReferenceCellType::Quadrilateral) + ); - fn ownership(&self, local_dof_index: usize) -> Ownership { - todo!(); + for cell in grid.entity_iter(ReferenceCellType::Point) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 1 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Point, cell.local_index()) + .unwrap() + .len(), + 1 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Interval) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 1 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Interval, cell.local_index()) + .unwrap() + .len(), + 3 + ); + } + for cell in grid.entity_iter(ReferenceCellType::Quadrilateral) { + assert_eq!( + space + .entity_dofs(ReferenceCellType::Quadrilateral, cell.local_index()) + .unwrap() + .len(), + 1 + ); + assert_eq!( + space + .entity_closure_dofs(ReferenceCellType::Quadrilateral, cell.local_index()) + .unwrap() + .len(), + 9 + ); + } } } diff --git a/ndfunctionspace/src/lib.rs b/ndfunctionspace/src/lib.rs index 98d7ad2..d9c490b 100644 --- a/ndfunctionspace/src/lib.rs +++ b/ndfunctionspace/src/lib.rs @@ -5,7 +5,7 @@ mod function_space; pub mod traits; -pub use function_space::FunctionSpace; +pub use function_space::SerialFunctionSpace; use ndelement as _; use ndgrid as _; diff --git a/ndfunctionspace/src/traits.rs b/ndfunctionspace/src/traits.rs index eaebdc1..3777a0a 100644 --- a/ndfunctionspace/src/traits.rs +++ b/ndfunctionspace/src/traits.rs @@ -1,13 +1,17 @@ //! Function space traits use ndelement::traits::FiniteElement; use ndgrid::{traits::Grid, types::Ownership}; +use std::fmt::Debug; +use std::hash::Hash; /// Function space. pub trait FunctionSpace { + /// Type used as identifier of different entity types + type EntityDescriptor: Debug + PartialEq + Eq + Clone + Copy + Hash; /// The type for the grid this function space is defined on - type Grid: Grid; + type Grid: Grid; /// The type for the finite element this grid is defined by - type FiniteElement: FiniteElement; + type FiniteElement: FiniteElement; /// The grid that this function space is defined on fn grid(&self) -> &Self::Grid; @@ -15,13 +19,20 @@ pub trait FunctionSpace { /// The finite elements used in this function space fn elements(&self) -> &[Self::FiniteElement]; - /// A list of cell indices that use the element with the given index - fn cells_by_element(&self, element_index: usize) -> Option<&[usize]>; + /// A list of entity indices that use the element with the given index + fn entities_by_element(&self, element_index: usize) -> Option<&[usize]>; /// Get the local DOFs numbers associated with the given entity fn entity_dofs( &self, - entity_type: ::EntityDescriptor, + entity_type: Self::EntityDescriptor, + entity_number: usize, + ) -> Option<&[usize]>; + + /// Get the local DOFs numbers associated with the closure of the given entity + fn entity_closure_dofs( + &self, + entity_type: Self::EntityDescriptor, entity_number: usize, ) -> Option<&[usize]>; From 917a46fce5c41f66d195063939d282fd5c8a7335 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Mon, 5 Jan 2026 15:00:24 +0000 Subject: [PATCH 5/7] clippy --- ndfunctionspace/Cargo.toml | 2 -- ndfunctionspace/src/function_space.rs | 2 +- 2 files changed, 1 insertion(+), 3 deletions(-) diff --git a/ndfunctionspace/Cargo.toml b/ndfunctionspace/Cargo.toml index 45cf5ec..6b89f8a 100644 --- a/ndfunctionspace/Cargo.toml +++ b/ndfunctionspace/Cargo.toml @@ -23,8 +23,6 @@ crate-type = ["lib", "cdylib"] [dependencies] ndelement = { path = "../ndelement" } ndgrid = { path = "../ndgrid" } -num = "0.4" -rlst = "0.4" itertools = "0.14.*" [dev-dependencies] diff --git a/ndfunctionspace/src/function_space.rs b/ndfunctionspace/src/function_space.rs index 8a7263f..6370930 100644 --- a/ndfunctionspace/src/function_space.rs +++ b/ndfunctionspace/src/function_space.rs @@ -70,7 +70,7 @@ impl< { for (t, i) in izip!(types, indices) { let reference_entity_dofs = element.entity_dofs(d, *i).unwrap(); - if reference_entity_dofs.len() > 0 && !reference_entity_dofs.is_empty() { + if !reference_entity_dofs.is_empty() && !reference_entity_dofs.is_empty() { let ed = &mut entity_dofs.get_mut(t).unwrap() [cell.topology().sub_entity(*t, *i)]; while ed.len() < reference_entity_dofs.len() { From 9c86a6d29824ed10ffdbfc21dbd8d26cdc287bf8 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Mon, 5 Jan 2026 15:13:47 +0000 Subject: [PATCH 6/7] unused deps --- ndfunctionspace/Cargo.toml | 4 ---- 1 file changed, 4 deletions(-) diff --git a/ndfunctionspace/Cargo.toml b/ndfunctionspace/Cargo.toml index 6b89f8a..6aebd03 100644 --- a/ndfunctionspace/Cargo.toml +++ b/ndfunctionspace/Cargo.toml @@ -25,10 +25,6 @@ ndelement = { path = "../ndelement" } ndgrid = { path = "../ndgrid" } itertools = "0.14.*" -[dev-dependencies] -approx = "0.5" -paste = "1.*" - [package.metadata.docs.rs] cargo-args = ["-Zunstable-options", "-Zrustdoc-scrape-examples"] rustdoc-args = [ "--html-in-header", "../katex-header.html" ] From f19f66ccc76eb36ff39fa9daa28f8a3e571baba5 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Mon, 5 Jan 2026 15:15:07 +0000 Subject: [PATCH 7/7] rlst version ++ --- ndelement/Cargo.toml | 2 +- ndgrid/Cargo.toml | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/ndelement/Cargo.toml b/ndelement/Cargo.toml index b1f52cf..f1075c5 100644 --- a/ndelement/Cargo.toml +++ b/ndelement/Cargo.toml @@ -26,7 +26,7 @@ bempp-quadrature = { version = "0.1.0" } itertools = "0.14.*" mpi = { version = "0.8.0", optional = true } num = "0.4" -rlst = "0.4" +rlst = "0.5" serde = { version = "1", features = ["derive"], optional = true } strum = "0.27" strum_macros = "0.27" diff --git a/ndgrid/Cargo.toml b/ndgrid/Cargo.toml index dea7027..efba65c 100644 --- a/ndgrid/Cargo.toml +++ b/ndgrid/Cargo.toml @@ -33,7 +33,7 @@ itertools = "0.14.*" mpi = { version = "0.8.0", optional = true } ndelement = { path = "../ndelement" } num = "0.4" -rlst = "0.4" +rlst = "0.5" serde = { version = "1", features = ["derive"], optional = true } ron = { version = "0.12", optional = true } scotch = { version = "0.2.1", optional = true }