From 4eb653da2397efc57deea6dbcde95ce01d7e95b8 Mon Sep 17 00:00:00 2001 From: Magnus Ulimoen Date: Sat, 19 Sep 2020 15:46:02 +0200 Subject: [PATCH] Documentation in SBP --- sbp/src/integrate.rs | 38 ++++++++++++++++++++++++++++++++++++++ sbp/src/lib.rs | 7 +++++++ sbp/src/operators.rs | 8 ++++++++ sbp/src/utils.rs | 7 +++++++ 4 files changed, 60 insertions(+) diff --git a/sbp/src/integrate.rs b/sbp/src/integrate.rs index 0a14f00..fe33a37 100644 --- a/sbp/src/integrate.rs +++ b/sbp/src/integrate.rs @@ -1,8 +1,22 @@ +//! Integration of explicit PDEs using different Butcher Tableaus +//! +//! Integration can be performed on all systems that can be represented +//! as using a transform into an [`ndarray::ArrayView`] for both the state +//! and the state difference. +//! +//! The integration functions are memory efficient, and relies +//! on the `k` parameter to hold the system state differences. +//! This parameter is tied to the Butcher Tableau + use super::Float; use ndarray::{ArrayView, ArrayViewMut}; +/// The Butcher Tableau, with the state transitions described as +/// [on wikipedia](https://en.wikipedia.org/wiki/Runge%E2%80%93Kutta_methods#Explicit_Runge%E2%80%93Kutta_methods). pub trait ButcherTableau { + /// This bound should not be overridden const S: usize = Self::B.len(); + /// Only the lower triangle will be used (explicit integration) const A: &'static [&'static [Float]]; const B: &'static [Float]; const C: &'static [Float]; @@ -137,6 +151,16 @@ impl EmbeddedButcherTableau for BogackiShampine { } #[allow(clippy::too_many_arguments)] +/// Integrates using the [`ButcherTableau`] specified. `rhs` should be the result +/// of the right hand side of $u_t = rhs$ +/// +/// rhs takes the old state and the current time, and outputs the state difference +/// in the first parameter +/// +/// Should be called as +/// ```rust,ignore +/// integrate::(...) +/// ``` pub fn integrate( mut rhs: RHS, prev: &F, @@ -191,6 +215,10 @@ pub fn integrate( } #[allow(clippy::too_many_arguments)] +/// Integrate using an [`EmbeddedButcherTableau`], else similar to [`integrate`] +/// +/// This produces two results, the most accurate result in `fut`, and the less accurate +/// result in `fut2`. This can be used for convergence testing and adaptive timesteps. pub fn integrate_embedded_rk( rhs: RHS, prev: &F, @@ -217,6 +245,16 @@ pub fn integrate_embedded_rk( #[cfg(feature = "rayon")] #[allow(clippy::too_many_arguments)] +/// Integrates a multigrid problem, much the same as [`integrate`], +/// using a `rayon` threadpool for parallelisation. +/// +/// note that `rhs` accepts the full system state, and is responsible +/// for computing the full state difference. +/// `rhs` can be a mutable closure, so buffers can be used +/// and mutated inside the closure. +/// +/// This function requires the `rayon` feature, and is not callable in +/// a `wasm` context. pub fn integrate_multigrid( mut rhs: RHS, prev: &[F], diff --git a/sbp/src/lib.rs b/sbp/src/lib.rs index bebe4ba..c63bca9 100644 --- a/sbp/src/lib.rs +++ b/sbp/src/lib.rs @@ -1,11 +1,14 @@ #![feature(min_specialization)] #![feature(core_intrinsics)] +/// Type used for floats, configure with the `f32` feature #[cfg(feature = "f32")] pub type Float = f32; #[cfg(not(feature = "f32"))] +/// Type used for floats, configure with the `f32` feature pub type Float = f64; +/// Associated constants for [`Float`] pub mod consts { #[cfg(feature = "f32")] pub use std::f32::consts::*; @@ -13,7 +16,11 @@ pub mod consts { pub use std::f64::consts::*; } +/// Grid and grid metrics pub mod grid; +/// RK operators and methods for implicit integration pub mod integrate; +/// SBP and interpolation operators pub mod operators; +/// General utilities pub mod utils; diff --git a/sbp/src/operators.rs b/sbp/src/operators.rs index b9fc3c2..380ef25 100644 --- a/sbp/src/operators.rs +++ b/sbp/src/operators.rs @@ -4,10 +4,14 @@ use crate::Float; use ndarray::{ArrayView1, ArrayView2, ArrayViewMut1, ArrayViewMut2}; +/// One-dimensional Summation By Parts operator pub trait SbpOperator1d: Send + Sync { + /// Differentiate on unit grid fn diff(&self, prev: ArrayView1, fut: ArrayViewMut1); + /// block component of `H`, with diagonal component 1 fn h(&self) -> &'static [Float]; + /// Whether the operator acts on a `h2` spaced computational grid fn is_h2(&self) -> bool { false } @@ -18,6 +22,7 @@ pub trait SbpOperator1d: Send + Sync { } pub trait SbpOperator1d2: SbpOperator1d { + /// Double differentiation on unit grid fn diff2(&self, prev: ArrayView1, fut: ArrayViewMut1); /// Result of H^-1 * d1, without the 1/h scaling fn d1(&self) -> &[Float]; @@ -103,6 +108,7 @@ impl SbpOperator2d for SBP { } pub trait UpwindOperator1d: SbpOperator1d + Send + Sync { + /// Dissipation operator fn diss(&self, prev: ArrayView1, fut: ArrayViewMut1); fn as_sbp(&self) -> &dyn SbpOperator1d; @@ -162,7 +168,9 @@ impl UpwindOperator2d for UO { } pub trait InterpolationOperator: Send + Sync { + /// Interpolation from a grid with twice resolution fn fine2coarse(&self, fine: ArrayView1, coarse: ArrayViewMut1); + /// Interpolation from a grid with half resolution fn coarse2fine(&self, coarse: ArrayView1, fine: ArrayViewMut1); } diff --git a/sbp/src/utils.rs b/sbp/src/utils.rs index 426ab4c..dd59187 100644 --- a/sbp/src/utils.rs +++ b/sbp/src/utils.rs @@ -13,6 +13,7 @@ use serde::{Deserialize, Serialize}; #[cfg_attr(feature = "serde", derive(Serialize, Deserialize))] #[derive(Copy, Clone, Debug, Default)] +/// struct to hold output for four directions pub struct Direction { pub north: T, pub south: T, @@ -45,6 +46,7 @@ impl Direction { east: f(self.east), } } + /// Combines two [`Direction`] into one pub fn zip(self, other: Direction) -> Direction<(T, U)> { Direction { north: (self.north, other.north), @@ -56,6 +58,8 @@ impl Direction { } impl Direction> { + /// Partially unwrap individual components in `self` or replace with + /// `other` pub fn unwrap_or(self, other: Direction) -> Direction { Direction { north: self.north.unwrap_or(other.north), @@ -66,6 +70,7 @@ impl Direction> { } } +/// Methods for borrowing impl Direction { pub fn north(&self) -> &T { &self.north @@ -93,6 +98,8 @@ impl Direction { } } +/// Linearly spaced parameters, apart from the boundaries which +/// only have a distance of `h/2` from the boundary pub fn h2linspace(start: Float, end: Float, n: usize) -> ndarray::Array1 { let h = (end - start) / (n - 2) as Float; ndarray::Array1::from_shape_fn(n, |i| match i {