From cb5006ff07ce00cb3e017ddb389cd5c362c00121 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Mon, 19 Aug 2024 14:26:08 +0100 Subject: [PATCH 1/5] update version numbers --- Cargo.toml | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/Cargo.toml b/Cargo.toml index d95a28de..e65e0b20 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -24,11 +24,11 @@ itertools = "0.13.*" mpi = { version = "0.8.*", optional = true } num = "0.4" lazy_static = "1.4" -ndelement = { git = "https://github.com/bempp/ndelement.git" } -ndgrid = { git = "https://github.com/bempp/ndgrid.git" } +ndelement = "0.1.0" +ndgrid = "0.1.0" rayon = "1.9" -rlst = { version = "0.2" } -green-kernels = { git = "https://github.com/bempp/green-kernels.git" } +rlst = "0.2" +green-kernels = "0.2.0" [dev-dependencies] approx = "0.5" From 9abf8287749e1ab0270ae266afe13fcdc44d2111 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Mon, 19 Aug 2024 15:07:57 +0100 Subject: [PATCH 2/5] restructure files --- examples/assembly.rs | 4 +- examples/test_parallel_assembly.rs | 6 +- src/assembly.rs | 9 +- src/assembly/batched.rs | 2 +- src/assembly/{batched => }/boundary.rs | 128 ++++++++++++++++-- .../boundary/adjoint_double_layer.rs | 24 ++-- .../{batched => }/boundary/double_layer.rs | 24 ++-- .../{batched => }/boundary/hypersingular.rs | 41 +++--- .../{batched => }/boundary/single_layer.rs | 24 ++-- src/assembly/common.rs | 5 +- src/assembly/{batched => }/potential.rs | 24 ++-- .../{batched => }/potential/double_layer.rs | 26 ++-- .../{batched => }/potential/single_layer.rs | 26 ++-- tests/compare_to_bempp_cl.rs | 28 ++-- tests/fmm.rs.bak | 8 +- 15 files changed, 242 insertions(+), 137 deletions(-) rename src/assembly/{batched => }/boundary.rs (92%) rename src/assembly/{batched => }/boundary/adjoint_double_layer.rs (87%) rename src/assembly/{batched => }/boundary/double_layer.rs (87%) rename src/assembly/{batched => }/boundary/hypersingular.rs (93%) rename src/assembly/{batched => }/boundary/single_layer.rs (82%) rename src/assembly/{batched => }/potential.rs (92%) rename src/assembly/{batched => }/potential/double_layer.rs (79%) rename src/assembly/{batched => }/potential/single_layer.rs (75%) diff --git a/examples/assembly.rs b/examples/assembly.rs index 5a80449b..73ac26e3 100644 --- a/examples/assembly.rs +++ b/examples/assembly.rs @@ -1,4 +1,4 @@ -use bempp::assembly::{batched, batched::BatchedAssembler}; +use bempp::assembly::{boundary, boundary::BoundaryAssembler}; use bempp::function::SerialFunctionSpace; use bempp::traits::FunctionSpace; use ndelement::ciarlet::LagrangeElementFamily; @@ -17,7 +17,7 @@ fn main() { let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); // Create an assembler for the Laplace single layer operator - let mut a = batched::LaplaceSingleLayerAssembler::::default(); + let mut a = boundary::LaplaceSingleLayerAssembler::::default(); // Adjust the quadrature degree for non-singular integrals on a triangle. // This makes the integrals use a quadrature rule with 16 points diff --git a/examples/test_parallel_assembly.rs b/examples/test_parallel_assembly.rs index 1d650544..9a881877 100644 --- a/examples/test_parallel_assembly.rs +++ b/examples/test_parallel_assembly.rs @@ -4,8 +4,8 @@ use approx::assert_relative_eq; #[cfg(feature = "mpi")] use bempp::{ - assembly::batched, - assembly::batched::BatchedAssembler, + assembly::boundary, + assembly::boundary::BoundaryAssembler, function::{ParallelFunctionSpace, SerialFunctionSpace}, traits::FunctionSpace, }; @@ -92,7 +92,7 @@ fn test_parallel_assembly_single_element_grid( let element = LagrangeElementFamily::::new(degree, cont); let space = ParallelFunctionSpace::new(&grid, &element); - let a = batched::LaplaceSingleLayerAssembler::::default(); + let a = boundary::LaplaceSingleLayerAssembler::::default(); let matrix = a.parallel_assemble_singular_into_csr(&space, &space); diff --git a/src/assembly.rs b/src/assembly.rs index b578a399..5e422b3a 100644 --- a/src/assembly.rs +++ b/src/assembly.rs @@ -1,12 +1,13 @@ //! Boundary operator assembly -pub mod batched; +pub mod boundary; +pub mod potential; pub(crate) mod common; pub mod fmm_tools; #[cfg(test)] mod test { - use super::batched::BatchedAssembler; use super::*; + use super::boundary::BoundaryAssembler; use crate::function::SerialFunctionSpace; use crate::traits::FunctionSpace; use cauchy::{c32, c64}; @@ -122,12 +123,12 @@ mod test { macro_rules! create_assembler { (Laplace, $operator:ident, $dtype:ident) => { paste! { - batched::[]::<[<$dtype>]>::default() + boundary::[]::<[<$dtype>]>::default() } }; (Helmholtz, $operator:ident, $dtype:ident) => { paste! { - batched::[]::<[<$dtype>]>::new(3.0) + boundary::[]::<[<$dtype>]>::new(3.0) } }; } diff --git a/src/assembly/batched.rs b/src/assembly/batched.rs index fecf40c3..7cf2b4c7 100644 --- a/src/assembly/batched.rs +++ b/src/assembly/batched.rs @@ -1,6 +1,6 @@ //! Batched dense assembly use crate::assembly::common::SparseMatrixData; -use green_kernels::types::GreenKernelEvalType; +use green_kernels::types::EvalType as GreenKernelEvalType; use rlst::{Array, BaseArray, VectorContainer}; mod boundary; diff --git a/src/assembly/batched/boundary.rs b/src/assembly/boundary.rs similarity index 92% rename from src/assembly/batched/boundary.rs rename to src/assembly/boundary.rs index 48e77b97..5406a2bf 100644 --- a/src/assembly/batched/boundary.rs +++ b/src/assembly/boundary.rs @@ -1,10 +1,19 @@ -//! Batched dense assembly of boundary operators +//! Assembly of boundary operators pub(crate) mod adjoint_double_layer; pub(crate) mod double_layer; pub(crate) mod hypersingular; pub(crate) mod single_layer; -use crate::assembly::common::{equal_grids, RawData2D, SparseMatrixData}; +pub use adjoint_double_layer::{LaplaceAdjointDoubleLayerAssembler, +HelmholtzAdjointDoubleLayerAssembler}; +pub use double_layer::{LaplaceDoubleLayerAssembler, +HelmholtzDoubleLayerAssembler}; +pub use single_layer::{LaplaceSingleLayerAssembler, +HelmholtzSingleLayerAssembler}; +pub use hypersingular::{LaplaceHypersingularAssembler, +HelmholtzHypersingularAssembler}; + +use crate::assembly::common::{equal_grids, RawData2D, SparseMatrixData, RlstArray}; use crate::quadrature::duffy::{ quadrilateral_duffy, quadrilateral_triangle_duffy, triangle_duffy, triangle_quadrilateral_duffy, }; @@ -27,8 +36,6 @@ use rlst::{ }; use std::collections::HashMap; -use super::RlstArray; - fn neighbours( test_grid: &TestGrid, trial_grid: &TrialGrid, @@ -106,7 +113,7 @@ fn assemble_batch_singular< TrialGrid: Grid, Element: FiniteElement + Sync, >( - assembler: &impl BatchedAssembler, + assembler: &impl BoundaryAssembler, deriv_size: usize, shape: [usize; 2], trial_cell_type: ReferenceCellType, @@ -222,7 +229,7 @@ fn assemble_batch_nonadjacent< TrialGrid: Grid, Element: FiniteElement + Sync, >( - assembler: &impl BatchedAssembler, + assembler: &impl BoundaryAssembler, deriv_size: usize, output: &RawData2D, trial_cell_type: ReferenceCellType, @@ -362,7 +369,7 @@ fn assemble_batch_singular_correction< TrialGrid: Grid, Element: FiniteElement + Sync, >( - assembler: &impl BatchedAssembler, + assembler: &impl BoundaryAssembler, deriv_size: usize, shape: [usize; 2], trial_cell_type: ReferenceCellType, @@ -500,8 +507,8 @@ fn get_pairs_if_smallest( Some(pairs) } -/// Options for a batched assembler -pub struct BatchedAssemblerOptions { +/// Options for a boundary assembler +pub struct BoundaryAssemblerOptions { /// Number of points used in quadrature for non-singular integrals quadrature_degrees: HashMap, /// Quadrature degrees to be used for singular integrals @@ -510,7 +517,7 @@ pub struct BatchedAssemblerOptions { batch_size: usize, } -impl Default for BatchedAssemblerOptions { +impl Default for BoundaryAssemblerOptions { fn default() -> Self { use ReferenceCellType::{Quadrilateral, Triangle}; Self { @@ -526,8 +533,8 @@ impl Default for BatchedAssemblerOptions { } } -pub trait BatchedAssembler: Sync + Sized { - //! Batched assembler +pub trait BoundaryAssembler: Sync + Sized { + //! Boundary assembler //! //! Assemble operators by processing batches of cells in parallel @@ -539,10 +546,10 @@ pub trait BatchedAssembler: Sync + Sized { const TABLE_DERIVS: usize; /// Get assembler options - fn options(&self) -> &BatchedAssemblerOptions; + fn options(&self) -> &BoundaryAssemblerOptions; /// Get mutable assembler options - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions; + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions; /// Set (non-singular) quadrature degree for a cell type fn quadrature_degree(&mut self, cell: ReferenceCellType, degree: usize) { @@ -1298,3 +1305,96 @@ pub trait BatchedAssembler: Sync + Sized { } } } + +#[cfg(test)] +mod test { + use super::*; + use crate::function::SerialFunctionSpace; + use crate::traits::FunctionSpace; + use approx::*; + use ndelement::ciarlet::LagrangeElementFamily; + use ndelement::types::Continuity; + use ndgrid::shapes::regular_sphere; + use rlst::rlst_dynamic_array2; + use rlst::RandomAccessByRef; + + #[test] + fn test_singular_dp0() { + let grid = regular_sphere::(0); + let element = LagrangeElementFamily::::new(0, Continuity::Discontinuous); + let space = SerialFunctionSpace::new(&grid, &element); + + let ndofs = space.global_size(); + + let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); + let assembler = LaplaceSingleLayerAssembler::::default(); + assembler.assemble_singular_into_dense(&mut matrix, &space, &space); + let csr = assembler.assemble_singular_into_csr(&space, &space); + + let indptr = csr.indptr(); + let indices = csr.indices(); + let data = csr.data(); + + let mut row = 0; + for (i, j) in indices.iter().enumerate() { + while i >= indptr[row + 1] { + row += 1; + } + assert_relative_eq!(*matrix.get([row, *j]).unwrap(), data[i], epsilon = 1e-8); + } + } + + #[test] + fn test_singular_p1() { + let grid = regular_sphere::(0); + let element = LagrangeElementFamily::::new(1, Continuity::Standard); + let space = SerialFunctionSpace::new(&grid, &element); + + let ndofs = space.global_size(); + + let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); + let assembler = LaplaceSingleLayerAssembler::::default(); + assembler.assemble_singular_into_dense(&mut matrix, &space, &space); + let csr = assembler.assemble_singular_into_csr(&space, &space); + + let indptr = csr.indptr(); + let indices = csr.indices(); + let data = csr.data(); + + let mut row = 0; + for (i, j) in indices.iter().enumerate() { + while i >= indptr[row + 1] { + row += 1; + } + assert_relative_eq!(*matrix.get([row, *j]).unwrap(), data[i], epsilon = 1e-8); + } + } + + #[test] + fn test_singular_dp0_p1() { + let grid = regular_sphere::(0); + let element0 = LagrangeElementFamily::::new(0, Continuity::Discontinuous); + let element1 = LagrangeElementFamily::::new(1, Continuity::Standard); + let space0 = SerialFunctionSpace::new(&grid, &element0); + let space1 = SerialFunctionSpace::new(&grid, &element1); + + let ndofs0 = space0.global_size(); + let ndofs1 = space1.global_size(); + + let mut matrix = rlst_dynamic_array2!(f64, [ndofs1, ndofs0]); + let assembler = LaplaceSingleLayerAssembler::::default(); + assembler.assemble_singular_into_dense(&mut matrix, &space0, &space1); + let csr = assembler.assemble_singular_into_csr(&space0, &space1); + let indptr = csr.indptr(); + let indices = csr.indices(); + let data = csr.data(); + + let mut row = 0; + for (i, j) in indices.iter().enumerate() { + while i >= indptr[row + 1] { + row += 1; + } + assert_relative_eq!(*matrix.get([row, *j]).unwrap(), data[i], epsilon = 1e-8); + } + } +} diff --git a/src/assembly/batched/boundary/adjoint_double_layer.rs b/src/assembly/boundary/adjoint_double_layer.rs similarity index 87% rename from src/assembly/batched/boundary/adjoint_double_layer.rs rename to src/assembly/boundary/adjoint_double_layer.rs index 504a7146..ee4074bf 100644 --- a/src/assembly/batched/boundary/adjoint_double_layer.rs +++ b/src/assembly/boundary/adjoint_double_layer.rs @@ -1,32 +1,32 @@ //! Adjoint double layer assemblers use super::{ - super::{GreenKernelEvalType, RlstArray}, - BatchedAssembler, BatchedAssemblerOptions, + BoundaryAssembler, BoundaryAssemblerOptions, }; +use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; /// Assembler for a Laplace adjoint double layer operator pub struct LaplaceAdjointDoubleLayerAssembler { kernel: Laplace3dKernel, - options: BatchedAssemblerOptions, + options: BoundaryAssemblerOptions, } impl Default for LaplaceAdjointDoubleLayerAssembler { fn default() -> Self { Self { kernel: Laplace3dKernel::::new(), - options: BatchedAssemblerOptions::default(), + options: BoundaryAssemblerOptions::default(), } } } -impl BatchedAssembler for LaplaceAdjointDoubleLayerAssembler { +impl BoundaryAssembler for LaplaceAdjointDoubleLayerAssembler { const DERIV_SIZE: usize = 4; const TABLE_DERIVS: usize = 0; type T = T; - fn options(&self) -> &BatchedAssemblerOptions { + fn options(&self) -> &BoundaryAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions { + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { &mut self.options } unsafe fn singular_kernel_value( @@ -76,27 +76,27 @@ impl BatchedAssembler for LaplaceAdjointDoubleLay /// Assembler for a Helmholtz adjoint double layer boundary operator pub struct HelmholtzAdjointDoubleLayerAssembler + MatrixInverse> { kernel: Helmholtz3dKernel, - options: BatchedAssemblerOptions, + options: BoundaryAssemblerOptions, } impl + MatrixInverse> HelmholtzAdjointDoubleLayerAssembler { /// Create a new assembler pub fn new(wavenumber: T::Real) -> Self { Self { kernel: Helmholtz3dKernel::::new(wavenumber), - options: BatchedAssemblerOptions::default(), + options: BoundaryAssemblerOptions::default(), } } } -impl + MatrixInverse> BatchedAssembler +impl + MatrixInverse> BoundaryAssembler for HelmholtzAdjointDoubleLayerAssembler { const DERIV_SIZE: usize = 4; const TABLE_DERIVS: usize = 0; type T = T; - fn options(&self) -> &BatchedAssemblerOptions { + fn options(&self) -> &BoundaryAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions { + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { &mut self.options } unsafe fn singular_kernel_value( diff --git a/src/assembly/batched/boundary/double_layer.rs b/src/assembly/boundary/double_layer.rs similarity index 87% rename from src/assembly/batched/boundary/double_layer.rs rename to src/assembly/boundary/double_layer.rs index de5c5d0a..68240487 100644 --- a/src/assembly/batched/boundary/double_layer.rs +++ b/src/assembly/boundary/double_layer.rs @@ -1,32 +1,32 @@ //! Double layer assemblers use super::{ - super::{GreenKernelEvalType, RlstArray}, - BatchedAssembler, BatchedAssemblerOptions, + BoundaryAssembler, BoundaryAssemblerOptions, }; +use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; /// Assembler for a Laplace double layer operator pub struct LaplaceDoubleLayerAssembler { kernel: Laplace3dKernel, - options: BatchedAssemblerOptions, + options: BoundaryAssemblerOptions, } impl Default for LaplaceDoubleLayerAssembler { fn default() -> Self { Self { kernel: Laplace3dKernel::::new(), - options: BatchedAssemblerOptions::default(), + options: BoundaryAssemblerOptions::default(), } } } -impl BatchedAssembler for LaplaceDoubleLayerAssembler { +impl BoundaryAssembler for LaplaceDoubleLayerAssembler { const DERIV_SIZE: usize = 4; const TABLE_DERIVS: usize = 0; type T = T; - fn options(&self) -> &BatchedAssemblerOptions { + fn options(&self) -> &BoundaryAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions { + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { &mut self.options } unsafe fn singular_kernel_value( @@ -76,27 +76,27 @@ impl BatchedAssembler for LaplaceDoubleLayerAssem /// Assembler for a Helmholtz double layer boundary operator pub struct HelmholtzDoubleLayerAssembler + MatrixInverse> { kernel: Helmholtz3dKernel, - options: BatchedAssemblerOptions, + options: BoundaryAssemblerOptions, } impl + MatrixInverse> HelmholtzDoubleLayerAssembler { /// Create a new assembler pub fn new(wavenumber: T::Real) -> Self { Self { kernel: Helmholtz3dKernel::::new(wavenumber), - options: BatchedAssemblerOptions::default(), + options: BoundaryAssemblerOptions::default(), } } } -impl + MatrixInverse> BatchedAssembler +impl + MatrixInverse> BoundaryAssembler for HelmholtzDoubleLayerAssembler { const DERIV_SIZE: usize = 4; const TABLE_DERIVS: usize = 0; type T = T; - fn options(&self) -> &BatchedAssemblerOptions { + fn options(&self) -> &BoundaryAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions { + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { &mut self.options } unsafe fn singular_kernel_value( diff --git a/src/assembly/batched/boundary/hypersingular.rs b/src/assembly/boundary/hypersingular.rs similarity index 93% rename from src/assembly/batched/boundary/hypersingular.rs rename to src/assembly/boundary/hypersingular.rs index f183292d..4a6d2533 100644 --- a/src/assembly/batched/boundary/hypersingular.rs +++ b/src/assembly/boundary/hypersingular.rs @@ -1,9 +1,8 @@ //! Hypersingular assemblers use super::{ - super::{GreenKernelEvalType, RlstArray, SparseMatrixData}, - BatchedAssembler, BatchedAssemblerOptions, + BoundaryAssembler, BoundaryAssemblerOptions, }; -use crate::assembly::common::equal_grids; +use crate::assembly::common::{GreenKernelEvalType, RlstArray, SparseMatrixData, equal_grids}; use crate::traits::FunctionSpace; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use ndelement::traits::FiniteElement; @@ -69,24 +68,24 @@ unsafe fn hyp_test_trial_product( /// Assembler for a Laplace hypersingular operator pub struct LaplaceHypersingularAssembler { kernel: Laplace3dKernel, - options: BatchedAssemblerOptions, + options: BoundaryAssemblerOptions, } impl Default for LaplaceHypersingularAssembler { fn default() -> Self { Self { kernel: Laplace3dKernel::::new(), - options: BatchedAssemblerOptions::default(), + options: BoundaryAssemblerOptions::default(), } } } -impl BatchedAssembler for LaplaceHypersingularAssembler { +impl BoundaryAssembler for LaplaceHypersingularAssembler { const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 1; type T = T; - fn options(&self) -> &BatchedAssemblerOptions { + fn options(&self) -> &BoundaryAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions { + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { &mut self.options } unsafe fn singular_kernel_value( @@ -153,27 +152,27 @@ impl BatchedAssembler for LaplaceHypersingularAss /// Assembler for curl-curl term of Helmholtz hypersingular operator struct HelmholtzHypersingularCurlCurlAssembler + MatrixInverse> { kernel: Helmholtz3dKernel, - options: BatchedAssemblerOptions, + options: BoundaryAssemblerOptions, } impl + MatrixInverse> HelmholtzHypersingularCurlCurlAssembler { /// Create a new assembler pub fn new(wavenumber: T::Real) -> Self { Self { kernel: Helmholtz3dKernel::::new(wavenumber), - options: BatchedAssemblerOptions::default(), + options: BoundaryAssemblerOptions::default(), } } } -impl + MatrixInverse> BatchedAssembler +impl + MatrixInverse> BoundaryAssembler for HelmholtzHypersingularCurlCurlAssembler { const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 1; type T = T; - fn options(&self) -> &BatchedAssemblerOptions { + fn options(&self) -> &BoundaryAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions { + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { &mut self.options } unsafe fn singular_kernel_value( @@ -241,7 +240,7 @@ impl + MatrixInverse> BatchedAssembler struct HelmholtzHypersingularNormalNormalAssembler + MatrixInverse> { kernel: Helmholtz3dKernel, wavenumber: T::Real, - options: BatchedAssemblerOptions, + options: BoundaryAssemblerOptions, } impl + MatrixInverse> HelmholtzHypersingularNormalNormalAssembler { /// Create a new assembler @@ -249,20 +248,20 @@ impl + MatrixInverse> HelmholtzHypersingularNormalNor Self { kernel: Helmholtz3dKernel::::new(wavenumber), wavenumber, - options: BatchedAssemblerOptions::default(), + options: BoundaryAssemblerOptions::default(), } } } -impl + MatrixInverse> BatchedAssembler +impl + MatrixInverse> BoundaryAssembler for HelmholtzHypersingularNormalNormalAssembler { const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 0; type T = T; - fn options(&self) -> &BatchedAssemblerOptions { + fn options(&self) -> &BoundaryAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions { + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { &mut self.options } unsafe fn singular_kernel_value( @@ -331,16 +330,16 @@ impl + MatrixInverse> HelmholtzHypersingularAssembler } } } -impl + MatrixInverse> BatchedAssembler +impl + MatrixInverse> BoundaryAssembler for HelmholtzHypersingularAssembler { const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 1; type T = T; - fn options(&self) -> &BatchedAssemblerOptions { + fn options(&self) -> &BoundaryAssemblerOptions { panic!("Cannot directly use HelmholtzHypersingularAssembler"); } - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions { + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { panic!("Cannot directly use HelmholtzHypersingularAssembler"); } fn quadrature_degree(&mut self, cell: ReferenceCellType, degree: usize) { diff --git a/src/assembly/batched/boundary/single_layer.rs b/src/assembly/boundary/single_layer.rs similarity index 82% rename from src/assembly/batched/boundary/single_layer.rs rename to src/assembly/boundary/single_layer.rs index 4976b491..9cbec60e 100644 --- a/src/assembly/batched/boundary/single_layer.rs +++ b/src/assembly/boundary/single_layer.rs @@ -1,32 +1,32 @@ //! Single layer assemblers use super::{ - super::{GreenKernelEvalType, RlstArray}, - BatchedAssembler, BatchedAssemblerOptions, + BoundaryAssembler, BoundaryAssemblerOptions, }; +use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; /// Assembler for a Laplace single layer operator pub struct LaplaceSingleLayerAssembler { kernel: Laplace3dKernel, - options: BatchedAssemblerOptions, + options: BoundaryAssemblerOptions, } impl Default for LaplaceSingleLayerAssembler { fn default() -> Self { Self { kernel: Laplace3dKernel::::new(), - options: BatchedAssemblerOptions::default(), + options: BoundaryAssemblerOptions::default(), } } } -impl BatchedAssembler for LaplaceSingleLayerAssembler { +impl BoundaryAssembler for LaplaceSingleLayerAssembler { const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 0; type T = T; - fn options(&self) -> &BatchedAssemblerOptions { + fn options(&self) -> &BoundaryAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions { + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { &mut self.options } unsafe fn singular_kernel_value( @@ -66,27 +66,27 @@ impl BatchedAssembler for LaplaceSingleLayerAssem /// Assembler for a Helmholtz single layer boundary operator pub struct HelmholtzSingleLayerAssembler + MatrixInverse> { kernel: Helmholtz3dKernel, - options: BatchedAssemblerOptions, + options: BoundaryAssemblerOptions, } impl + MatrixInverse> HelmholtzSingleLayerAssembler { /// Create a new assembler pub fn new(wavenumber: T::Real) -> Self { Self { kernel: Helmholtz3dKernel::::new(wavenumber), - options: BatchedAssemblerOptions::default(), + options: BoundaryAssemblerOptions::default(), } } } -impl + MatrixInverse> BatchedAssembler +impl + MatrixInverse> BoundaryAssembler for HelmholtzSingleLayerAssembler { const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 0; type T = T; - fn options(&self) -> &BatchedAssemblerOptions { + fn options(&self) -> &BoundaryAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedAssemblerOptions { + fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { &mut self.options } unsafe fn singular_kernel_value( diff --git a/src/assembly/common.rs b/src/assembly/common.rs index 2e1630eb..61ce04c1 100644 --- a/src/assembly/common.rs +++ b/src/assembly/common.rs @@ -1,6 +1,9 @@ //! Common utility functions use ndgrid::traits::Grid; -use rlst::{MatrixInverse, RlstScalar}; +use rlst::{MatrixInverse, RlstScalar, Array, BaseArray, VectorContainer}; +pub(crate) use green_kernels::types::EvalType as GreenKernelEvalType; + +pub(crate) type RlstArray = Array, DIM>, DIM>; pub(crate) fn equal_grids( test_grid: &TestGrid, diff --git a/src/assembly/batched/potential.rs b/src/assembly/potential.rs similarity index 92% rename from src/assembly/batched/potential.rs rename to src/assembly/potential.rs index 077f14ab..b7d47d7c 100644 --- a/src/assembly/batched/potential.rs +++ b/src/assembly/potential.rs @@ -1,7 +1,10 @@ -//! Batched dense assembly of potential operators +//! Assembly of potential operators pub(crate) mod double_layer; pub(crate) mod single_layer; +pub use single_layer::{HelmholtzSingleLayerPotentialAssembler, LaplaceSingleLayerPotentialAssembler}; +pub use double_layer::{HelmholtzDoubleLayerPotentialAssembler, LaplaceDoubleLayerPotentialAssembler}; + use crate::assembly::common::RawData2D; use crate::quadrature::simplex_rules::simplex_rule; use crate::traits::FunctionSpace; @@ -14,8 +17,7 @@ use rlst::{ RawAccess, RawAccessMut, RlstScalar, Shape, UnsafeRandomAccessByRef, }; use std::collections::HashMap; - -use super::RlstArray; +use crate::assembly::common::RlstArray; /// Assemble the contribution to the terms of a matrix for a batch of non-adjacent cells #[allow(clippy::too_many_arguments)] @@ -24,7 +26,7 @@ fn assemble_batch< G: Grid, Element: FiniteElement + Sync, >( - assembler: &impl BatchedPotentialAssembler, + assembler: &impl PotentialAssembler, deriv_size: usize, output: &RawData2D, cell_type: ReferenceCellType, @@ -87,15 +89,15 @@ fn assemble_batch< 1 } -/// Options for a batched assembler -pub struct BatchedPotentialAssemblerOptions { +/// Options for a potential assembler +pub struct PotentialAssemblerOptions { /// Number of points used in quadrature for non-singular integrals quadrature_degrees: HashMap, /// Maximum size of each batch of cells to send to an assembly function batch_size: usize, } -impl Default for BatchedPotentialAssemblerOptions { +impl Default for PotentialAssemblerOptions { fn default() -> Self { use ReferenceCellType::{Quadrilateral, Triangle}; Self { @@ -105,8 +107,8 @@ impl Default for BatchedPotentialAssemblerOptions { } } -pub trait BatchedPotentialAssembler: Sync + Sized { - //! Batched potential assembler +pub trait PotentialAssembler: Sync + Sized { + //! Potential assembler //! //! Assemble potential operators by processing batches of cells in parallel @@ -116,10 +118,10 @@ pub trait BatchedPotentialAssembler: Sync + Sized { const DERIV_SIZE: usize; /// Get assembler options - fn options(&self) -> &BatchedPotentialAssemblerOptions; + fn options(&self) -> &PotentialAssemblerOptions; /// Get mutable assembler options - fn options_mut(&mut self) -> &mut BatchedPotentialAssemblerOptions; + fn options_mut(&mut self) -> &mut PotentialAssemblerOptions; /// Set (non-singular) quadrature degree for a cell type fn quadrature_degree(&mut self, cell: ReferenceCellType, degree: usize) { diff --git a/src/assembly/batched/potential/double_layer.rs b/src/assembly/potential/double_layer.rs similarity index 79% rename from src/assembly/batched/potential/double_layer.rs rename to src/assembly/potential/double_layer.rs index 77cb6c4b..91184613 100644 --- a/src/assembly/batched/potential/double_layer.rs +++ b/src/assembly/potential/double_layer.rs @@ -1,35 +1,35 @@ -//! Batched dense assembly of boundary operators +//! Double layer potential assemblers use super::{ - super::{GreenKernelEvalType, RlstArray}, - BatchedPotentialAssembler, BatchedPotentialAssemblerOptions, + PotentialAssembler, PotentialAssemblerOptions, }; +use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; /// Assembler for a Laplace double layer potential operator pub struct LaplaceDoubleLayerPotentialAssembler { kernel: Laplace3dKernel, - options: BatchedPotentialAssemblerOptions, + options: PotentialAssemblerOptions, } impl Default for LaplaceDoubleLayerPotentialAssembler { fn default() -> Self { Self { kernel: Laplace3dKernel::::new(), - options: BatchedPotentialAssemblerOptions::default(), + options: PotentialAssemblerOptions::default(), } } } -impl BatchedPotentialAssembler +impl PotentialAssembler for LaplaceDoubleLayerPotentialAssembler { const DERIV_SIZE: usize = 4; type T = T; - fn options(&self) -> &BatchedPotentialAssemblerOptions { + fn options(&self) -> &PotentialAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedPotentialAssemblerOptions { + fn options_mut(&mut self) -> &mut PotentialAssemblerOptions { &mut self.options } @@ -62,28 +62,28 @@ impl BatchedPotentialAssembler /// Assembler for a Helmholtz double layer potential operator pub struct HelmholtzDoubleLayerPotentialAssembler + MatrixInverse> { kernel: Helmholtz3dKernel, - options: BatchedPotentialAssemblerOptions, + options: PotentialAssemblerOptions, } impl + MatrixInverse> HelmholtzDoubleLayerPotentialAssembler { /// Create a new assembler pub fn new(wavenumber: T::Real) -> Self { Self { kernel: Helmholtz3dKernel::::new(wavenumber), - options: BatchedPotentialAssemblerOptions::default(), + options: PotentialAssemblerOptions::default(), } } } -impl + MatrixInverse> BatchedPotentialAssembler +impl + MatrixInverse> PotentialAssembler for HelmholtzDoubleLayerPotentialAssembler { const DERIV_SIZE: usize = 4; type T = T; - fn options(&self) -> &BatchedPotentialAssemblerOptions { + fn options(&self) -> &PotentialAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedPotentialAssemblerOptions { + fn options_mut(&mut self) -> &mut PotentialAssemblerOptions { &mut self.options } diff --git a/src/assembly/batched/potential/single_layer.rs b/src/assembly/potential/single_layer.rs similarity index 75% rename from src/assembly/batched/potential/single_layer.rs rename to src/assembly/potential/single_layer.rs index d777aad5..26d19540 100644 --- a/src/assembly/batched/potential/single_layer.rs +++ b/src/assembly/potential/single_layer.rs @@ -1,35 +1,35 @@ -//! Batched dense assembly of boundary operators +//! Double layer potential assemblers use super::{ - super::{GreenKernelEvalType, RlstArray}, - BatchedPotentialAssembler, BatchedPotentialAssemblerOptions, + PotentialAssembler, PotentialAssemblerOptions, }; +use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; /// Assembler for a Laplace single layer potential operator pub struct LaplaceSingleLayerPotentialAssembler { kernel: Laplace3dKernel, - options: BatchedPotentialAssemblerOptions, + options: PotentialAssemblerOptions, } impl Default for LaplaceSingleLayerPotentialAssembler { fn default() -> Self { Self { kernel: Laplace3dKernel::::new(), - options: BatchedPotentialAssemblerOptions::default(), + options: PotentialAssemblerOptions::default(), } } } -impl BatchedPotentialAssembler +impl PotentialAssembler for LaplaceSingleLayerPotentialAssembler { const DERIV_SIZE: usize = 1; type T = T; - fn options(&self) -> &BatchedPotentialAssemblerOptions { + fn options(&self) -> &PotentialAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedPotentialAssemblerOptions { + fn options_mut(&mut self) -> &mut PotentialAssemblerOptions { &mut self.options } @@ -57,28 +57,28 @@ impl BatchedPotentialAssembler /// Assembler for a Helmholtz single layer potential operator pub struct HelmholtzSingleLayerPotentialAssembler + MatrixInverse> { kernel: Helmholtz3dKernel, - options: BatchedPotentialAssemblerOptions, + options: PotentialAssemblerOptions, } impl + MatrixInverse> HelmholtzSingleLayerPotentialAssembler { /// Create a new assembler pub fn new(wavenumber: T::Real) -> Self { Self { kernel: Helmholtz3dKernel::::new(wavenumber), - options: BatchedPotentialAssemblerOptions::default(), + options: PotentialAssemblerOptions::default(), } } } -impl + MatrixInverse> BatchedPotentialAssembler +impl + MatrixInverse> PotentialAssembler for HelmholtzSingleLayerPotentialAssembler { const DERIV_SIZE: usize = 1; type T = T; - fn options(&self) -> &BatchedPotentialAssemblerOptions { + fn options(&self) -> &PotentialAssemblerOptions { &self.options } - fn options_mut(&mut self) -> &mut BatchedPotentialAssemblerOptions { + fn options_mut(&mut self) -> &mut PotentialAssemblerOptions { &mut self.options } diff --git a/tests/compare_to_bempp_cl.rs b/tests/compare_to_bempp_cl.rs index de6d5cde..02e552ac 100644 --- a/tests/compare_to_bempp_cl.rs +++ b/tests/compare_to_bempp_cl.rs @@ -1,5 +1,5 @@ use approx::*; -use bempp::assembly::{batched, batched::BatchedAssembler, batched::BatchedPotentialAssembler}; +use bempp::assembly::{boundary, boundary::BoundaryAssembler, potential, potential::PotentialAssembler}; use bempp::function::SerialFunctionSpace; use bempp::traits::FunctionSpace; use cauchy::c64; @@ -18,7 +18,7 @@ fn test_laplace_single_layer_dp0_dp0() { let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let a = batched::LaplaceSingleLayerAssembler::::default(); + let a = boundary::LaplaceSingleLayerAssembler::::default(); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -41,7 +41,7 @@ fn test_laplace_double_layer_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let a = batched::LaplaceDoubleLayerAssembler::::default(); + let a = boundary::LaplaceDoubleLayerAssembler::::default(); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -63,7 +63,7 @@ fn test_laplace_adjoint_double_layer_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let a = batched::LaplaceAdjointDoubleLayerAssembler::::default(); + let a = boundary::LaplaceAdjointDoubleLayerAssembler::::default(); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -86,7 +86,7 @@ fn test_laplace_hypersingular_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let a = batched::LaplaceHypersingularAssembler::::default(); + let a = boundary::LaplaceHypersingularAssembler::::default(); a.assemble_into_dense(&mut matrix, &space, &space); for i in 0..ndofs { @@ -105,7 +105,7 @@ fn test_laplace_hypersingular_p1_p1() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let a = batched::LaplaceHypersingularAssembler::::default(); + let a = boundary::LaplaceHypersingularAssembler::::default(); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -134,7 +134,7 @@ fn test_helmholtz_single_layer_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(c64, [ndofs, ndofs]); - let a = batched::HelmholtzSingleLayerAssembler::::new(3.0); + let a = boundary::HelmholtzSingleLayerAssembler::::new(3.0); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -157,7 +157,7 @@ fn test_helmholtz_double_layer_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(c64, [ndofs, ndofs]); - let a = batched::HelmholtzDoubleLayerAssembler::::new(3.0); + let a = boundary::HelmholtzDoubleLayerAssembler::::new(3.0); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -179,7 +179,7 @@ fn test_helmholtz_adjoint_double_layer_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(c64, [ndofs, ndofs]); - let a = batched::HelmholtzAdjointDoubleLayerAssembler::::new(3.0); + let a = boundary::HelmholtzAdjointDoubleLayerAssembler::::new(3.0); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -202,7 +202,7 @@ fn test_helmholtz_hypersingular_p1_p1() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(c64, [ndofs, ndofs]); - let a = batched::HelmholtzHypersingularAssembler::::new(3.0); + let a = boundary::HelmholtzHypersingularAssembler::::new(3.0); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -237,7 +237,7 @@ fn test_laplace_single_layer_potential_dp0() { *points.get_mut([1, 1]).unwrap() = 2.0; *points.get_mut([2, 2]).unwrap() = 2.0; - let a = batched::LaplaceSingleLayerPotentialAssembler::::default(); + let a = potential::LaplaceSingleLayerPotentialAssembler::::default(); a.assemble_into_dense(&mut matrix, &space, &points); // Compare to result from bempp-cl @@ -265,7 +265,7 @@ fn test_helmholtz_single_layer_potential_dp0() { *points.get_mut([1, 1]).unwrap() = 2.0; *points.get_mut([2, 2]).unwrap() = 2.0; - let a = batched::HelmholtzSingleLayerPotentialAssembler::::new(3.0); + let a = potential::HelmholtzSingleLayerPotentialAssembler::::new(3.0); a.assemble_into_dense(&mut matrix, &space, &points); // Compare to result from bempp-cl @@ -293,7 +293,7 @@ fn test_laplace_double_layer_potential_dp0() { *points.get_mut([1, 1]).unwrap() = 2.0; *points.get_mut([2, 2]).unwrap() = 2.0; - let a = batched::LaplaceDoubleLayerPotentialAssembler::::default(); + let a = potential::LaplaceDoubleLayerPotentialAssembler::::default(); a.assemble_into_dense(&mut matrix, &space, &points); // Compare to result from bempp-cl @@ -321,7 +321,7 @@ fn test_helmholtz_double_layer_potential_dp0() { *points.get_mut([1, 1]).unwrap() = 2.0; *points.get_mut([2, 2]).unwrap() = 2.0; - let a = batched::HelmholtzDoubleLayerPotentialAssembler::::new(3.0); + let a = potential::HelmholtzDoubleLayerPotentialAssembler::::new(3.0); a.assemble_into_dense(&mut matrix, &space, &points); // Compare to result from bempp-cl diff --git a/tests/fmm.rs.bak b/tests/fmm.rs.bak index 77a90f64..38ea78f5 100644 --- a/tests/fmm.rs.bak +++ b/tests/fmm.rs.bak @@ -1,6 +1,6 @@ use approx::*; -use bempp::assembly::batched::BatchedAssembler; -use bempp::assembly::{batched, fmm_tools}; +use bempp::assembly::boundary::BoundaryAssembler; +use bempp::assembly::{boundary, fmm_tools}; use bempp::function::SerialFunctionSpace; use ndgrid::{shapes::regular_sphere, traits::Grid}; use bempp::traits::FunctionSpace; @@ -40,7 +40,7 @@ fn fmm_prototype + // Compute dense let mut matrix = rlst_dynamic_array2!(f64, [test_ndofs, trial_ndofs]); - let mut a = batched::LaplaceSingleLayerAssembler::::default(); + let mut a = boundary::LaplaceSingleLayerAssembler::::default(); a.quadrature_degree(ReferenceCellType::Triangle, npts); a.assemble_into_dense(&mut matrix, trial_space, test_space); @@ -115,7 +115,7 @@ fn fmm_matvec + S // Compute dense let mut matrix = rlst_dynamic_array2!(f64, [test_ndofs, trial_ndofs]); - let mut a = batched::LaplaceSingleLayerAssembler::::default(); + let mut a = boundary::LaplaceSingleLayerAssembler::::default(); a.quadrature_degree(ReferenceCellType::Triangle, npts); a.assemble_into_dense(&mut matrix, trial_space, test_space); From 586d40df396a160a32984a9b9bf67fbe6d595210 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Mon, 19 Aug 2024 15:08:14 +0100 Subject: [PATCH 3/5] fmt --- src/assembly.rs | 4 ++-- src/assembly/boundary.rs | 18 ++++++++---------- src/assembly/boundary/adjoint_double_layer.rs | 4 +--- src/assembly/boundary/double_layer.rs | 4 +--- src/assembly/boundary/hypersingular.rs | 6 ++---- src/assembly/boundary/single_layer.rs | 4 +--- src/assembly/common.rs | 7 ++++--- src/assembly/potential.rs | 10 +++++++--- src/assembly/potential/double_layer.rs | 8 ++------ src/assembly/potential/single_layer.rs | 8 ++------ tests/compare_to_bempp_cl.rs | 4 +++- 11 files changed, 33 insertions(+), 44 deletions(-) diff --git a/src/assembly.rs b/src/assembly.rs index 5e422b3a..e173aeb3 100644 --- a/src/assembly.rs +++ b/src/assembly.rs @@ -1,13 +1,13 @@ //! Boundary operator assembly pub mod boundary; -pub mod potential; pub(crate) mod common; pub mod fmm_tools; +pub mod potential; #[cfg(test)] mod test { - use super::*; use super::boundary::BoundaryAssembler; + use super::*; use crate::function::SerialFunctionSpace; use crate::traits::FunctionSpace; use cauchy::{c32, c64}; diff --git a/src/assembly/boundary.rs b/src/assembly/boundary.rs index 5406a2bf..95c7338b 100644 --- a/src/assembly/boundary.rs +++ b/src/assembly/boundary.rs @@ -4,16 +4,14 @@ pub(crate) mod double_layer; pub(crate) mod hypersingular; pub(crate) mod single_layer; -pub use adjoint_double_layer::{LaplaceAdjointDoubleLayerAssembler, -HelmholtzAdjointDoubleLayerAssembler}; -pub use double_layer::{LaplaceDoubleLayerAssembler, -HelmholtzDoubleLayerAssembler}; -pub use single_layer::{LaplaceSingleLayerAssembler, -HelmholtzSingleLayerAssembler}; -pub use hypersingular::{LaplaceHypersingularAssembler, -HelmholtzHypersingularAssembler}; - -use crate::assembly::common::{equal_grids, RawData2D, SparseMatrixData, RlstArray}; +pub use adjoint_double_layer::{ + HelmholtzAdjointDoubleLayerAssembler, LaplaceAdjointDoubleLayerAssembler, +}; +pub use double_layer::{HelmholtzDoubleLayerAssembler, LaplaceDoubleLayerAssembler}; +pub use hypersingular::{HelmholtzHypersingularAssembler, LaplaceHypersingularAssembler}; +pub use single_layer::{HelmholtzSingleLayerAssembler, LaplaceSingleLayerAssembler}; + +use crate::assembly::common::{equal_grids, RawData2D, RlstArray, SparseMatrixData}; use crate::quadrature::duffy::{ quadrilateral_duffy, quadrilateral_triangle_duffy, triangle_duffy, triangle_quadrilateral_duffy, }; diff --git a/src/assembly/boundary/adjoint_double_layer.rs b/src/assembly/boundary/adjoint_double_layer.rs index ee4074bf..0331b2bb 100644 --- a/src/assembly/boundary/adjoint_double_layer.rs +++ b/src/assembly/boundary/adjoint_double_layer.rs @@ -1,7 +1,5 @@ //! Adjoint double layer assemblers -use super::{ - BoundaryAssembler, BoundaryAssemblerOptions, -}; +use super::{BoundaryAssembler, BoundaryAssemblerOptions}; use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; diff --git a/src/assembly/boundary/double_layer.rs b/src/assembly/boundary/double_layer.rs index 68240487..3c65e1b8 100644 --- a/src/assembly/boundary/double_layer.rs +++ b/src/assembly/boundary/double_layer.rs @@ -1,7 +1,5 @@ //! Double layer assemblers -use super::{ - BoundaryAssembler, BoundaryAssemblerOptions, -}; +use super::{BoundaryAssembler, BoundaryAssemblerOptions}; use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; diff --git a/src/assembly/boundary/hypersingular.rs b/src/assembly/boundary/hypersingular.rs index 4a6d2533..479df259 100644 --- a/src/assembly/boundary/hypersingular.rs +++ b/src/assembly/boundary/hypersingular.rs @@ -1,8 +1,6 @@ //! Hypersingular assemblers -use super::{ - BoundaryAssembler, BoundaryAssemblerOptions, -}; -use crate::assembly::common::{GreenKernelEvalType, RlstArray, SparseMatrixData, equal_grids}; +use super::{BoundaryAssembler, BoundaryAssemblerOptions}; +use crate::assembly::common::{equal_grids, GreenKernelEvalType, RlstArray, SparseMatrixData}; use crate::traits::FunctionSpace; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use ndelement::traits::FiniteElement; diff --git a/src/assembly/boundary/single_layer.rs b/src/assembly/boundary/single_layer.rs index 9cbec60e..8d816280 100644 --- a/src/assembly/boundary/single_layer.rs +++ b/src/assembly/boundary/single_layer.rs @@ -1,7 +1,5 @@ //! Single layer assemblers -use super::{ - BoundaryAssembler, BoundaryAssemblerOptions, -}; +use super::{BoundaryAssembler, BoundaryAssemblerOptions}; use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; diff --git a/src/assembly/common.rs b/src/assembly/common.rs index 61ce04c1..fb3db393 100644 --- a/src/assembly/common.rs +++ b/src/assembly/common.rs @@ -1,9 +1,10 @@ //! Common utility functions -use ndgrid::traits::Grid; -use rlst::{MatrixInverse, RlstScalar, Array, BaseArray, VectorContainer}; pub(crate) use green_kernels::types::EvalType as GreenKernelEvalType; +use ndgrid::traits::Grid; +use rlst::{Array, BaseArray, MatrixInverse, RlstScalar, VectorContainer}; -pub(crate) type RlstArray = Array, DIM>, DIM>; +pub(crate) type RlstArray = + Array, DIM>, DIM>; pub(crate) fn equal_grids( test_grid: &TestGrid, diff --git a/src/assembly/potential.rs b/src/assembly/potential.rs index b7d47d7c..d0d1eaa1 100644 --- a/src/assembly/potential.rs +++ b/src/assembly/potential.rs @@ -2,10 +2,15 @@ pub(crate) mod double_layer; pub(crate) mod single_layer; -pub use single_layer::{HelmholtzSingleLayerPotentialAssembler, LaplaceSingleLayerPotentialAssembler}; -pub use double_layer::{HelmholtzDoubleLayerPotentialAssembler, LaplaceDoubleLayerPotentialAssembler}; +pub use double_layer::{ + HelmholtzDoubleLayerPotentialAssembler, LaplaceDoubleLayerPotentialAssembler, +}; +pub use single_layer::{ + HelmholtzSingleLayerPotentialAssembler, LaplaceSingleLayerPotentialAssembler, +}; use crate::assembly::common::RawData2D; +use crate::assembly::common::RlstArray; use crate::quadrature::simplex_rules::simplex_rule; use crate::traits::FunctionSpace; use ndelement::traits::FiniteElement; @@ -17,7 +22,6 @@ use rlst::{ RawAccess, RawAccessMut, RlstScalar, Shape, UnsafeRandomAccessByRef, }; use std::collections::HashMap; -use crate::assembly::common::RlstArray; /// Assemble the contribution to the terms of a matrix for a batch of non-adjacent cells #[allow(clippy::too_many_arguments)] diff --git a/src/assembly/potential/double_layer.rs b/src/assembly/potential/double_layer.rs index 91184613..cbe51102 100644 --- a/src/assembly/potential/double_layer.rs +++ b/src/assembly/potential/double_layer.rs @@ -1,7 +1,5 @@ //! Double layer potential assemblers -use super::{ - PotentialAssembler, PotentialAssemblerOptions, -}; +use super::{PotentialAssembler, PotentialAssemblerOptions}; use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; @@ -20,9 +18,7 @@ impl Default for LaplaceDoubleLayerPotentialAssem } } -impl PotentialAssembler - for LaplaceDoubleLayerPotentialAssembler -{ +impl PotentialAssembler for LaplaceDoubleLayerPotentialAssembler { const DERIV_SIZE: usize = 4; type T = T; diff --git a/src/assembly/potential/single_layer.rs b/src/assembly/potential/single_layer.rs index 26d19540..03d0ff06 100644 --- a/src/assembly/potential/single_layer.rs +++ b/src/assembly/potential/single_layer.rs @@ -1,7 +1,5 @@ //! Double layer potential assemblers -use super::{ - PotentialAssembler, PotentialAssemblerOptions, -}; +use super::{PotentialAssembler, PotentialAssemblerOptions}; use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; @@ -20,9 +18,7 @@ impl Default for LaplaceSingleLayerPotentialAssem } } -impl PotentialAssembler - for LaplaceSingleLayerPotentialAssembler -{ +impl PotentialAssembler for LaplaceSingleLayerPotentialAssembler { const DERIV_SIZE: usize = 1; type T = T; diff --git a/tests/compare_to_bempp_cl.rs b/tests/compare_to_bempp_cl.rs index 02e552ac..fd709216 100644 --- a/tests/compare_to_bempp_cl.rs +++ b/tests/compare_to_bempp_cl.rs @@ -1,5 +1,7 @@ use approx::*; -use bempp::assembly::{boundary, boundary::BoundaryAssembler, potential, potential::PotentialAssembler}; +use bempp::assembly::{ + boundary, boundary::BoundaryAssembler, potential, potential::PotentialAssembler, +}; use bempp::function::SerialFunctionSpace; use bempp::traits::FunctionSpace; use cauchy::c64; From 76af9cbb694d6cfdf38b5221a0df9b207eb7dcb9 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Mon, 19 Aug 2024 16:49:47 +0100 Subject: [PATCH 4/5] reduce code duplication --- benches/assembly_benchmark.rs | 4 +- examples/assembly.rs | 2 +- examples/test_parallel_assembly.rs | 2 +- src/assembly.rs | 4 +- src/assembly/boundary.rs | 16 +- src/assembly/boundary/adjoint_double_layer.rs | 92 ++------ src/assembly/boundary/double_layer.rs | 90 ++----- src/assembly/boundary/hypersingular.rs | 221 +++++++++--------- src/assembly/boundary/single_layer.rs | 80 ++----- src/assembly/potential.rs | 8 +- src/assembly/potential/double_layer.rs | 75 ++---- src/assembly/potential/single_layer.rs | 72 ++---- tests/compare_to_bempp_cl.rs | 26 +-- 13 files changed, 232 insertions(+), 460 deletions(-) diff --git a/benches/assembly_benchmark.rs b/benches/assembly_benchmark.rs index ba46e374..1306d765 100644 --- a/benches/assembly_benchmark.rs +++ b/benches/assembly_benchmark.rs @@ -1,4 +1,4 @@ -use bempp::assembly::{batched, batched::BatchedAssembler}; +use bempp::assembly::{boundary, boundary::BoundaryAssembler}; use bempp::function::SerialFunctionSpace; use bempp::traits::FunctionSpace; use criterion::{criterion_group, criterion_main, Criterion}; @@ -19,7 +19,7 @@ pub fn assembly_parts_benchmark(c: &mut Criterion) { let mut matrix = rlst_dynamic_array2!(f64, [space.global_size(), space.global_size()]); let colouring = space.cell_colouring(); - let mut a = batched::LaplaceSingleLayerAssembler::::default(); + let mut a = boundary::SingleLayerAssembler::::new_laplace(); a.quadrature_degree(ReferenceCellType::Triangle, 16); a.singular_quadrature_degree( (ReferenceCellType::Triangle, ReferenceCellType::Triangle), diff --git a/examples/assembly.rs b/examples/assembly.rs index 73ac26e3..4dd25635 100644 --- a/examples/assembly.rs +++ b/examples/assembly.rs @@ -17,7 +17,7 @@ fn main() { let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); // Create an assembler for the Laplace single layer operator - let mut a = boundary::LaplaceSingleLayerAssembler::::default(); + let mut a = boundary::SingleLayerAssembler::::new_laplace(); // Adjust the quadrature degree for non-singular integrals on a triangle. // This makes the integrals use a quadrature rule with 16 points diff --git a/examples/test_parallel_assembly.rs b/examples/test_parallel_assembly.rs index 9a881877..7f15660a 100644 --- a/examples/test_parallel_assembly.rs +++ b/examples/test_parallel_assembly.rs @@ -92,7 +92,7 @@ fn test_parallel_assembly_single_element_grid( let element = LagrangeElementFamily::::new(degree, cont); let space = ParallelFunctionSpace::new(&grid, &element); - let a = boundary::LaplaceSingleLayerAssembler::::default(); + let a = boundary::SingleLayerAssembler::::new_laplace(); let matrix = a.parallel_assemble_singular_into_csr(&space, &space); diff --git a/src/assembly.rs b/src/assembly.rs index e173aeb3..7398c597 100644 --- a/src/assembly.rs +++ b/src/assembly.rs @@ -123,12 +123,12 @@ mod test { macro_rules! create_assembler { (Laplace, $operator:ident, $dtype:ident) => { paste! { - boundary::[]::<[<$dtype>]>::default() + boundary::[<$operator Assembler>]::<[<$dtype>], _>::new_laplace() } }; (Helmholtz, $operator:ident, $dtype:ident) => { paste! { - boundary::[]::<[<$dtype>]>::new(3.0) + boundary::[<$operator Assembler>]::<[<$dtype>], _>::new_helmholtz(3.0) } }; } diff --git a/src/assembly/boundary.rs b/src/assembly/boundary.rs index 95c7338b..48afc285 100644 --- a/src/assembly/boundary.rs +++ b/src/assembly/boundary.rs @@ -4,12 +4,10 @@ pub(crate) mod double_layer; pub(crate) mod hypersingular; pub(crate) mod single_layer; -pub use adjoint_double_layer::{ - HelmholtzAdjointDoubleLayerAssembler, LaplaceAdjointDoubleLayerAssembler, -}; -pub use double_layer::{HelmholtzDoubleLayerAssembler, LaplaceDoubleLayerAssembler}; -pub use hypersingular::{HelmholtzHypersingularAssembler, LaplaceHypersingularAssembler}; -pub use single_layer::{HelmholtzSingleLayerAssembler, LaplaceSingleLayerAssembler}; +pub use adjoint_double_layer::AdjointDoubleLayerAssembler; +pub use double_layer::DoubleLayerAssembler; +pub use hypersingular::HypersingularAssembler; +pub use single_layer::SingleLayerAssembler; use crate::assembly::common::{equal_grids, RawData2D, RlstArray, SparseMatrixData}; use crate::quadrature::duffy::{ @@ -1325,7 +1323,7 @@ mod test { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let assembler = LaplaceSingleLayerAssembler::::default(); + let assembler = SingleLayerAssembler::::new_laplace(); assembler.assemble_singular_into_dense(&mut matrix, &space, &space); let csr = assembler.assemble_singular_into_csr(&space, &space); @@ -1351,7 +1349,7 @@ mod test { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let assembler = LaplaceSingleLayerAssembler::::default(); + let assembler = SingleLayerAssembler::::new_laplace(); assembler.assemble_singular_into_dense(&mut matrix, &space, &space); let csr = assembler.assemble_singular_into_csr(&space, &space); @@ -1380,7 +1378,7 @@ mod test { let ndofs1 = space1.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs1, ndofs0]); - let assembler = LaplaceSingleLayerAssembler::::default(); + let assembler = SingleLayerAssembler::::new_laplace(); assembler.assemble_singular_into_dense(&mut matrix, &space0, &space1); let csr = assembler.assemble_singular_into_csr(&space0, &space1); let indptr = csr.indptr(); diff --git a/src/assembly/boundary/adjoint_double_layer.rs b/src/assembly/boundary/adjoint_double_layer.rs index 0331b2bb..36f3d000 100644 --- a/src/assembly/boundary/adjoint_double_layer.rs +++ b/src/assembly/boundary/adjoint_double_layer.rs @@ -4,89 +4,37 @@ use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; -/// Assembler for a Laplace adjoint double layer operator -pub struct LaplaceAdjointDoubleLayerAssembler { - kernel: Laplace3dKernel, +/// Assembler for a adjoint double layer operator +pub struct AdjointDoubleLayerAssembler> { + kernel: K, options: BoundaryAssemblerOptions, } -impl Default for LaplaceAdjointDoubleLayerAssembler { - fn default() -> Self { +impl> AdjointDoubleLayerAssembler { + /// Create a new adjoint double layer assembler + pub fn new(kernel: K) -> Self { Self { - kernel: Laplace3dKernel::::new(), + kernel, options: BoundaryAssemblerOptions::default(), } } } -impl BoundaryAssembler for LaplaceAdjointDoubleLayerAssembler { - const DERIV_SIZE: usize = 4; - const TABLE_DERIVS: usize = 0; - type T = T; - fn options(&self) -> &BoundaryAssemblerOptions { - &self.options - } - fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { - &mut self.options - } - unsafe fn singular_kernel_value( - &self, - k: &RlstArray, - test_normals: &RlstArray, - _trial_normals: &RlstArray, - index: usize, - ) -> T { - -*k.get_unchecked([1, index]) - * num::cast::(*test_normals.get_unchecked([0, index])).unwrap() - - *k.get_unchecked([2, index]) - * num::cast::(*test_normals.get_unchecked([1, index])).unwrap() - - *k.get_unchecked([3, index]) - * num::cast::(*test_normals.get_unchecked([2, index])).unwrap() - } - unsafe fn nonsingular_kernel_value( - &self, - k: &RlstArray, - test_normals: &RlstArray, - _trial_normals: &RlstArray, - test_index: usize, - trial_index: usize, - ) -> T { - -*k.get_unchecked([1, test_index, trial_index]) - * num::cast::(*test_normals.get_unchecked([0, test_index])).unwrap() - - *k.get_unchecked([2, test_index, trial_index]) - * num::cast::(*test_normals.get_unchecked([1, test_index])).unwrap() - - *k.get_unchecked([3, test_index, trial_index]) - * num::cast::(*test_normals.get_unchecked([2, test_index])).unwrap() +impl AdjointDoubleLayerAssembler> { + /// Create a new Laplace adjoint double layer assembler + pub fn new_laplace() -> Self { + Self::new(Laplace3dKernel::::new()) } - fn kernel_assemble_pairwise_st( - &self, - sources: &[T::Real], - targets: &[T::Real], - result: &mut [T], - ) { - self.kernel - .assemble_pairwise_st(GreenKernelEvalType::ValueDeriv, sources, targets, result); - } - fn kernel_assemble_st(&self, sources: &[T::Real], targets: &[T::Real], result: &mut [T]) { - self.kernel - .assemble_st(GreenKernelEvalType::ValueDeriv, sources, targets, result); - } -} - -/// Assembler for a Helmholtz adjoint double layer boundary operator -pub struct HelmholtzAdjointDoubleLayerAssembler + MatrixInverse> { - kernel: Helmholtz3dKernel, - options: BoundaryAssemblerOptions, } -impl + MatrixInverse> HelmholtzAdjointDoubleLayerAssembler { - /// Create a new assembler - pub fn new(wavenumber: T::Real) -> Self { - Self { - kernel: Helmholtz3dKernel::::new(wavenumber), - options: BoundaryAssemblerOptions::default(), - } +impl + MatrixInverse> + AdjointDoubleLayerAssembler> +{ + /// Create a new Helmholtz adjoint double layer assembler + pub fn new_helmholtz(wavenumber: T::Real) -> Self { + Self::new(Helmholtz3dKernel::::new(wavenumber)) } } -impl + MatrixInverse> BoundaryAssembler - for HelmholtzAdjointDoubleLayerAssembler + +impl> BoundaryAssembler + for AdjointDoubleLayerAssembler { const DERIV_SIZE: usize = 4; const TABLE_DERIVS: usize = 0; diff --git a/src/assembly/boundary/double_layer.rs b/src/assembly/boundary/double_layer.rs index 3c65e1b8..51083622 100644 --- a/src/assembly/boundary/double_layer.rs +++ b/src/assembly/boundary/double_layer.rs @@ -4,89 +4,35 @@ use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; -/// Assembler for a Laplace double layer operator -pub struct LaplaceDoubleLayerAssembler { - kernel: Laplace3dKernel, +/// Assembler for a double layer operator +pub struct DoubleLayerAssembler> { + kernel: K, options: BoundaryAssemblerOptions, } -impl Default for LaplaceDoubleLayerAssembler { - fn default() -> Self { +impl> DoubleLayerAssembler { + /// Create a new double layer assembler + pub fn new(kernel: K) -> Self { Self { - kernel: Laplace3dKernel::::new(), + kernel, options: BoundaryAssemblerOptions::default(), } } } -impl BoundaryAssembler for LaplaceDoubleLayerAssembler { - const DERIV_SIZE: usize = 4; - const TABLE_DERIVS: usize = 0; - type T = T; - fn options(&self) -> &BoundaryAssemblerOptions { - &self.options - } - fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { - &mut self.options - } - unsafe fn singular_kernel_value( - &self, - k: &RlstArray, - _test_normals: &RlstArray, - trial_normals: &RlstArray, - index: usize, - ) -> T { - *k.get_unchecked([1, index]) - * num::cast::(*trial_normals.get_unchecked([0, index])).unwrap() - + *k.get_unchecked([2, index]) - * num::cast::(*trial_normals.get_unchecked([1, index])).unwrap() - + *k.get_unchecked([3, index]) - * num::cast::(*trial_normals.get_unchecked([2, index])).unwrap() - } - unsafe fn nonsingular_kernel_value( - &self, - k: &RlstArray, - _test_normals: &RlstArray, - trial_normals: &RlstArray, - test_index: usize, - trial_index: usize, - ) -> T { - *k.get_unchecked([1, test_index, trial_index]) - * num::cast::(*trial_normals.get_unchecked([0, trial_index])).unwrap() - + *k.get_unchecked([2, test_index, trial_index]) - * num::cast::(*trial_normals.get_unchecked([1, trial_index])).unwrap() - + *k.get_unchecked([3, test_index, trial_index]) - * num::cast::(*trial_normals.get_unchecked([2, trial_index])).unwrap() - } - fn kernel_assemble_pairwise_st( - &self, - sources: &[T::Real], - targets: &[T::Real], - result: &mut [T], - ) { - self.kernel - .assemble_pairwise_st(GreenKernelEvalType::ValueDeriv, sources, targets, result); - } - fn kernel_assemble_st(&self, sources: &[T::Real], targets: &[T::Real], result: &mut [T]) { - self.kernel - .assemble_st(GreenKernelEvalType::ValueDeriv, sources, targets, result); +impl DoubleLayerAssembler> { + /// Create a new Laplace double layer assembler + pub fn new_laplace() -> Self { + Self::new(Laplace3dKernel::::new()) } } - -/// Assembler for a Helmholtz double layer boundary operator -pub struct HelmholtzDoubleLayerAssembler + MatrixInverse> { - kernel: Helmholtz3dKernel, - options: BoundaryAssemblerOptions, -} -impl + MatrixInverse> HelmholtzDoubleLayerAssembler { - /// Create a new assembler - pub fn new(wavenumber: T::Real) -> Self { - Self { - kernel: Helmholtz3dKernel::::new(wavenumber), - options: BoundaryAssemblerOptions::default(), - } +impl + MatrixInverse> DoubleLayerAssembler> { + /// Create a new Helmholtz double layer assembler + pub fn new_helmholtz(wavenumber: T::Real) -> Self { + Self::new(Helmholtz3dKernel::::new(wavenumber)) } } -impl + MatrixInverse> BoundaryAssembler - for HelmholtzDoubleLayerAssembler + +impl> BoundaryAssembler + for DoubleLayerAssembler { const DERIV_SIZE: usize = 4; const TABLE_DERIVS: usize = 0; diff --git a/src/assembly/boundary/hypersingular.rs b/src/assembly/boundary/hypersingular.rs index 479df259..468cb9aa 100644 --- a/src/assembly/boundary/hypersingular.rs +++ b/src/assembly/boundary/hypersingular.rs @@ -63,20 +63,36 @@ unsafe fn hyp_test_trial_product( .unwrap() } -/// Assembler for a Laplace hypersingular operator -pub struct LaplaceHypersingularAssembler { - kernel: Laplace3dKernel, +/// Assembler for a hypersingular operator +pub struct HypersingularAssembler> { + kernel: K, options: BoundaryAssemblerOptions, } -impl Default for LaplaceHypersingularAssembler { - fn default() -> Self { +impl> HypersingularAssembler { + /// Create a new hypersingular assembler + pub fn new(kernel: K) -> Self { Self { - kernel: Laplace3dKernel::::new(), + kernel, options: BoundaryAssemblerOptions::default(), } } } -impl BoundaryAssembler for LaplaceHypersingularAssembler { +impl HypersingularAssembler> { + /// Create a new Laplace hypersingular assembler + pub fn new_laplace() -> Self { + Self::new(Laplace3dKernel::::new()) + } +} +impl + MatrixInverse> HypersingularAssembler> { + /// Create a new Helmholtz hypersingular assembler + pub fn new_helmholtz(wavenumber: T::Real) -> Self { + Self::new(Helmholtz3dKernel::::new(wavenumber)) + } +} + +impl BoundaryAssembler + for HypersingularAssembler> +{ const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 1; type T = T; @@ -148,21 +164,20 @@ impl BoundaryAssembler for LaplaceHypersingularAs } /// Assembler for curl-curl term of Helmholtz hypersingular operator -struct HelmholtzHypersingularCurlCurlAssembler + MatrixInverse> { - kernel: Helmholtz3dKernel, - options: BoundaryAssemblerOptions, +struct HelmholtzHypersingularCurlCurlAssembler<'a, T: RlstScalar + MatrixInverse> { + kernel: &'a Helmholtz3dKernel, + options: &'a BoundaryAssemblerOptions, } -impl + MatrixInverse> HelmholtzHypersingularCurlCurlAssembler { +impl<'a, T: RlstScalar + MatrixInverse> + HelmholtzHypersingularCurlCurlAssembler<'a, T> +{ /// Create a new assembler - pub fn new(wavenumber: T::Real) -> Self { - Self { - kernel: Helmholtz3dKernel::::new(wavenumber), - options: BoundaryAssemblerOptions::default(), - } + pub fn new(kernel: &'a Helmholtz3dKernel, options: &'a BoundaryAssemblerOptions) -> Self { + Self { kernel, options } } } -impl + MatrixInverse> BoundaryAssembler - for HelmholtzHypersingularCurlCurlAssembler +impl<'a, T: RlstScalar + MatrixInverse> BoundaryAssembler + for HelmholtzHypersingularCurlCurlAssembler<'a, T> { const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 1; @@ -171,7 +186,7 @@ impl + MatrixInverse> BoundaryAssembler &self.options } fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { - &mut self.options + panic!("Cannot get mutable options") } unsafe fn singular_kernel_value( &self, @@ -235,23 +250,20 @@ impl + MatrixInverse> BoundaryAssembler } /// Assembler for normal normal term of Helmholtz hypersingular boundary operator -struct HelmholtzHypersingularNormalNormalAssembler + MatrixInverse> { - kernel: Helmholtz3dKernel, - wavenumber: T::Real, - options: BoundaryAssemblerOptions, +struct HelmholtzHypersingularNormalNormalAssembler<'a, T: RlstScalar + MatrixInverse> { + kernel: &'a Helmholtz3dKernel, + options: &'a BoundaryAssemblerOptions, } -impl + MatrixInverse> HelmholtzHypersingularNormalNormalAssembler { +impl<'a, T: RlstScalar + MatrixInverse> + HelmholtzHypersingularNormalNormalAssembler<'a, T> +{ /// Create a new assembler - pub fn new(wavenumber: T::Real) -> Self { - Self { - kernel: Helmholtz3dKernel::::new(wavenumber), - wavenumber, - options: BoundaryAssemblerOptions::default(), - } + pub fn new(kernel: &'a Helmholtz3dKernel, options: &'a BoundaryAssemblerOptions) -> Self { + Self { kernel, options } } } -impl + MatrixInverse> BoundaryAssembler - for HelmholtzHypersingularNormalNormalAssembler +impl<'a, T: RlstScalar + MatrixInverse> BoundaryAssembler + for HelmholtzHypersingularNormalNormalAssembler<'a, T> { const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 0; @@ -260,7 +272,7 @@ impl + MatrixInverse> BoundaryAssembler &self.options } fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { - &mut self.options + panic!("Cannot get mutable options") } unsafe fn singular_kernel_value( &self, @@ -269,14 +281,17 @@ impl + MatrixInverse> BoundaryAssembler trial_normals: &RlstArray, index: usize, ) -> T { - -num::cast::(self.wavenumber.powi(2)).unwrap() - * *k.get_unchecked([0, index]) - * (num::cast::(*trial_normals.get_unchecked([0, index])).unwrap() - * num::cast::(*test_normals.get_unchecked([0, index])).unwrap() - + num::cast::(*trial_normals.get_unchecked([1, index])).unwrap() - * num::cast::(*test_normals.get_unchecked([1, index])).unwrap() - + num::cast::(*trial_normals.get_unchecked([2, index])).unwrap() - * num::cast::(*test_normals.get_unchecked([2, index])).unwrap()) + -*k.get_unchecked([0, index]) + * num::cast::( + self.kernel.wavenumber.powi(2) + * (*trial_normals.get_unchecked([0, index]) + * *test_normals.get_unchecked([0, index]) + + *trial_normals.get_unchecked([1, index]) + * *test_normals.get_unchecked([1, index]) + + *trial_normals.get_unchecked([2, index]) + * *test_normals.get_unchecked([2, index])), + ) + .unwrap() } unsafe fn nonsingular_kernel_value( &self, @@ -286,7 +301,7 @@ impl + MatrixInverse> BoundaryAssembler test_index: usize, trial_index: usize, ) -> T { - -num::cast::(self.wavenumber.powi(2)).unwrap() + -num::cast::(self.kernel.wavenumber.powi(2)).unwrap() * *k.get_unchecked([0, test_index, trial_index]) * (num::cast::(*trial_normals.get_unchecked([0, trial_index])).unwrap() * num::cast::(*test_normals.get_unchecked([0, test_index])).unwrap() @@ -312,51 +327,17 @@ impl + MatrixInverse> BoundaryAssembler } } -/// Assembler for curl-curl term of Helmholtz hypersingular operator -pub struct HelmholtzHypersingularAssembler + MatrixInverse> { - curl_curl_assembler: HelmholtzHypersingularCurlCurlAssembler, - normal_normal_assembler: HelmholtzHypersingularNormalNormalAssembler, -} -impl + MatrixInverse> HelmholtzHypersingularAssembler { - /// Create a new assembler - pub fn new(wavenumber: T::Real) -> Self { - Self { - curl_curl_assembler: HelmholtzHypersingularCurlCurlAssembler::::new(wavenumber), - normal_normal_assembler: HelmholtzHypersingularNormalNormalAssembler::::new( - wavenumber, - ), - } - } -} impl + MatrixInverse> BoundaryAssembler - for HelmholtzHypersingularAssembler + for HypersingularAssembler> { const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 1; type T = T; fn options(&self) -> &BoundaryAssemblerOptions { - panic!("Cannot directly use HelmholtzHypersingularAssembler"); + &self.options } fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { - panic!("Cannot directly use HelmholtzHypersingularAssembler"); - } - fn quadrature_degree(&mut self, cell: ReferenceCellType, degree: usize) { - self.curl_curl_assembler.quadrature_degree(cell, degree); - self.normal_normal_assembler.quadrature_degree(cell, degree); - } - fn singular_quadrature_degree( - &mut self, - cells: (ReferenceCellType, ReferenceCellType), - degree: usize, - ) { - self.curl_curl_assembler - .singular_quadrature_degree(cells, degree); - self.normal_normal_assembler - .singular_quadrature_degree(cells, degree); - } - fn batch_size(&mut self, size: usize) { - self.curl_curl_assembler.batch_size(size); - self.normal_normal_assembler.batch_size(size); + &mut self.options } unsafe fn singular_kernel_value( @@ -366,7 +347,7 @@ impl + MatrixInverse> BoundaryAssembler _trial_normals: &RlstArray, _index: usize, ) -> T { - panic!("Cannot directly use HelmholtzHypersingularAssembler"); + panic!("Cannot directly use HypersingularAssembler for Helmholtz"); } unsafe fn nonsingular_kernel_value( &self, @@ -376,7 +357,7 @@ impl + MatrixInverse> BoundaryAssembler _test_index: usize, _trial_index: usize, ) -> T { - panic!("Cannot directly use HelmholtzHypersingularAssembler"); + panic!("Cannot directly use HypersingularAssembler for Helmholtz"); } fn kernel_assemble_pairwise_st( &self, @@ -384,10 +365,10 @@ impl + MatrixInverse> BoundaryAssembler _targets: &[T::Real], _result: &mut [T], ) { - panic!("Cannot directly use HelmholtzHypersingularAssembler"); + panic!("Cannot directly use HypersingularAssembler for Helmholtz"); } fn kernel_assemble_st(&self, _sources: &[T::Real], _targets: &[T::Real], _result: &mut [T]) { - panic!("Cannot directly use HelmholtzHypersingularAssembler"); + panic!("Cannot directly use HypersingularAssembler for Helmholtz"); } fn assemble_singular< @@ -400,12 +381,22 @@ impl + MatrixInverse> BoundaryAssembler trial_space: &(impl FunctionSpace + Sync), test_space: &(impl FunctionSpace + Sync), ) -> SparseMatrixData { - let mut curlcurl = self - .curl_curl_assembler - .assemble_singular::(shape, trial_space, test_space); + let curl_curl_assembler = + HelmholtzHypersingularCurlCurlAssembler::::new(&self.kernel, &self.options); + let normal_normal_assembler = + HelmholtzHypersingularNormalNormalAssembler::::new(&self.kernel, &self.options); + + let mut curlcurl = curl_curl_assembler.assemble_singular::( + shape, + trial_space, + test_space, + ); curlcurl.add( - self.normal_normal_assembler - .assemble_singular::(shape, trial_space, test_space), + normal_normal_assembler.assemble_singular::( + shape, + trial_space, + test_space, + ), ); curlcurl } @@ -425,20 +416,23 @@ impl + MatrixInverse> BoundaryAssembler return SparseMatrixData::new(shape); } - let mut curlcurl = self - .curl_curl_assembler + let curl_curl_assembler = + HelmholtzHypersingularCurlCurlAssembler::::new(&self.kernel, &self.options); + let normal_normal_assembler = + HelmholtzHypersingularNormalNormalAssembler::::new(&self.kernel, &self.options); + + let mut curlcurl = curl_curl_assembler .assemble_singular_correction::( shape, trial_space, test_space, ); curlcurl.add( - self.normal_normal_assembler - .assemble_singular_correction::( - shape, - trial_space, - test_space, - ), + normal_normal_assembler.assemble_singular_correction::( + shape, + trial_space, + test_space, + ), ); curlcurl } @@ -465,21 +459,24 @@ impl + MatrixInverse> BoundaryAssembler panic!("Matrix has wrong shape"); } - self.curl_curl_assembler - .assemble_nonsingular_into_dense::( - output, - trial_space, - test_space, - trial_colouring, - test_colouring, - ); - self.normal_normal_assembler - .assemble_nonsingular_into_dense::( - output, - trial_space, - test_space, - trial_colouring, - test_colouring, - ); + let curl_curl_assembler = + HelmholtzHypersingularCurlCurlAssembler::::new(&self.kernel, &self.options); + let normal_normal_assembler = + HelmholtzHypersingularNormalNormalAssembler::::new(&self.kernel, &self.options); + + curl_curl_assembler.assemble_nonsingular_into_dense::( + output, + trial_space, + test_space, + trial_colouring, + test_colouring, + ); + normal_normal_assembler.assemble_nonsingular_into_dense::( + output, + trial_space, + test_space, + trial_colouring, + test_colouring, + ); } } diff --git a/src/assembly/boundary/single_layer.rs b/src/assembly/boundary/single_layer.rs index 8d816280..e3f812c8 100644 --- a/src/assembly/boundary/single_layer.rs +++ b/src/assembly/boundary/single_layer.rs @@ -4,79 +4,35 @@ use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; -/// Assembler for a Laplace single layer operator -pub struct LaplaceSingleLayerAssembler { - kernel: Laplace3dKernel, +/// Assembler for a single layer operator +pub struct SingleLayerAssembler> { + kernel: K, options: BoundaryAssemblerOptions, } -impl Default for LaplaceSingleLayerAssembler { - fn default() -> Self { +impl> SingleLayerAssembler { + /// Create a new single layer assembler + pub fn new(kernel: K) -> Self { Self { - kernel: Laplace3dKernel::::new(), + kernel, options: BoundaryAssemblerOptions::default(), } } } -impl BoundaryAssembler for LaplaceSingleLayerAssembler { - const DERIV_SIZE: usize = 1; - const TABLE_DERIVS: usize = 0; - type T = T; - fn options(&self) -> &BoundaryAssemblerOptions { - &self.options - } - fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { - &mut self.options - } - unsafe fn singular_kernel_value( - &self, - k: &RlstArray, - _test_normals: &RlstArray, - _trial_normals: &RlstArray, - index: usize, - ) -> T { - *k.get_unchecked([0, index]) - } - unsafe fn nonsingular_kernel_value( - &self, - k: &RlstArray, - _test_normals: &RlstArray, - _trial_normals: &RlstArray, - test_index: usize, - trial_index: usize, - ) -> T { - *k.get_unchecked([0, test_index, trial_index]) - } - fn kernel_assemble_pairwise_st( - &self, - sources: &[T::Real], - targets: &[T::Real], - result: &mut [T], - ) { - self.kernel - .assemble_pairwise_st(GreenKernelEvalType::Value, sources, targets, result); - } - fn kernel_assemble_st(&self, sources: &[T::Real], targets: &[T::Real], result: &mut [T]) { - self.kernel - .assemble_st(GreenKernelEvalType::Value, sources, targets, result); +impl SingleLayerAssembler> { + /// Create a new Laplace single layer assembler + pub fn new_laplace() -> Self { + Self::new(Laplace3dKernel::::new()) } } - -/// Assembler for a Helmholtz single layer boundary operator -pub struct HelmholtzSingleLayerAssembler + MatrixInverse> { - kernel: Helmholtz3dKernel, - options: BoundaryAssemblerOptions, -} -impl + MatrixInverse> HelmholtzSingleLayerAssembler { - /// Create a new assembler - pub fn new(wavenumber: T::Real) -> Self { - Self { - kernel: Helmholtz3dKernel::::new(wavenumber), - options: BoundaryAssemblerOptions::default(), - } +impl + MatrixInverse> SingleLayerAssembler> { + /// Create a new Helmholtz single layer assembler + pub fn new_helmholtz(wavenumber: T::Real) -> Self { + Self::new(Helmholtz3dKernel::::new(wavenumber)) } } -impl + MatrixInverse> BoundaryAssembler - for HelmholtzSingleLayerAssembler + +impl> BoundaryAssembler + for SingleLayerAssembler { const DERIV_SIZE: usize = 1; const TABLE_DERIVS: usize = 0; diff --git a/src/assembly/potential.rs b/src/assembly/potential.rs index d0d1eaa1..207c2117 100644 --- a/src/assembly/potential.rs +++ b/src/assembly/potential.rs @@ -2,12 +2,8 @@ pub(crate) mod double_layer; pub(crate) mod single_layer; -pub use double_layer::{ - HelmholtzDoubleLayerPotentialAssembler, LaplaceDoubleLayerPotentialAssembler, -}; -pub use single_layer::{ - HelmholtzSingleLayerPotentialAssembler, LaplaceSingleLayerPotentialAssembler, -}; +pub use double_layer::DoubleLayerPotentialAssembler; +pub use single_layer::SingleLayerPotentialAssembler; use crate::assembly::common::RawData2D; use crate::assembly::common::RlstArray; diff --git a/src/assembly/potential/double_layer.rs b/src/assembly/potential/double_layer.rs index cbe51102..0ae223da 100644 --- a/src/assembly/potential/double_layer.rs +++ b/src/assembly/potential/double_layer.rs @@ -4,74 +4,37 @@ use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; -/// Assembler for a Laplace double layer potential operator -pub struct LaplaceDoubleLayerPotentialAssembler { - kernel: Laplace3dKernel, +/// Assembler for a double layer potential operator +pub struct DoubleLayerPotentialAssembler> { + kernel: K, options: PotentialAssemblerOptions, } -impl Default for LaplaceDoubleLayerPotentialAssembler { - fn default() -> Self { +impl> DoubleLayerPotentialAssembler { + /// Create a new double layer potential assembler + pub fn new(kernel: K) -> Self { Self { - kernel: Laplace3dKernel::::new(), + kernel, options: PotentialAssemblerOptions::default(), } } } - -impl PotentialAssembler for LaplaceDoubleLayerPotentialAssembler { - const DERIV_SIZE: usize = 4; - type T = T; - - fn options(&self) -> &PotentialAssemblerOptions { - &self.options +impl DoubleLayerPotentialAssembler> { + /// Create a new Laplace double layer potential assembler + pub fn new_laplace() -> Self { + Self::new(Laplace3dKernel::::new()) } - fn options_mut(&mut self) -> &mut PotentialAssemblerOptions { - &mut self.options - } - - unsafe fn kernel_value( - &self, - k: &RlstArray, - normals: &RlstArray, - index: usize, - point_index: usize, - ) -> T { - -*k.get_unchecked([1, index, point_index]) - * num::cast::(*normals.get_unchecked([0, index])).unwrap() - - *k.get_unchecked([2, index, point_index]) - * num::cast::(*normals.get_unchecked([1, index])).unwrap() - - *k.get_unchecked([3, index, point_index]) - * num::cast::(*normals.get_unchecked([2, index])).unwrap() - } - - fn kernel_assemble_st( - &self, - sources: &[::Real], - targets: &[::Real], - result: &mut [Self::T], - ) { - self.kernel - .assemble_st(GreenKernelEvalType::ValueDeriv, sources, targets, result); - } -} - -/// Assembler for a Helmholtz double layer potential operator -pub struct HelmholtzDoubleLayerPotentialAssembler + MatrixInverse> { - kernel: Helmholtz3dKernel, - options: PotentialAssemblerOptions, } -impl + MatrixInverse> HelmholtzDoubleLayerPotentialAssembler { - /// Create a new assembler - pub fn new(wavenumber: T::Real) -> Self { - Self { - kernel: Helmholtz3dKernel::::new(wavenumber), - options: PotentialAssemblerOptions::default(), - } +impl + MatrixInverse> + DoubleLayerPotentialAssembler> +{ + /// Create a new Helmholtz double layer potential assembler + pub fn new_helmholtz(wavenumber: T::Real) -> Self { + Self::new(Helmholtz3dKernel::::new(wavenumber)) } } -impl + MatrixInverse> PotentialAssembler - for HelmholtzDoubleLayerPotentialAssembler +impl> PotentialAssembler + for DoubleLayerPotentialAssembler { const DERIV_SIZE: usize = 4; type T = T; diff --git a/src/assembly/potential/single_layer.rs b/src/assembly/potential/single_layer.rs index 03d0ff06..657f04b3 100644 --- a/src/assembly/potential/single_layer.rs +++ b/src/assembly/potential/single_layer.rs @@ -1,72 +1,40 @@ -//! Double layer potential assemblers +//! Single layer potential assemblers use super::{PotentialAssembler, PotentialAssemblerOptions}; use crate::assembly::common::{GreenKernelEvalType, RlstArray}; use green_kernels::{helmholtz_3d::Helmholtz3dKernel, laplace_3d::Laplace3dKernel, traits::Kernel}; use rlst::{MatrixInverse, RlstScalar, UnsafeRandomAccessByRef}; -/// Assembler for a Laplace single layer potential operator -pub struct LaplaceSingleLayerPotentialAssembler { - kernel: Laplace3dKernel, +/// Assembler for a single layer potential operator +pub struct SingleLayerPotentialAssembler> { + kernel: K, options: PotentialAssemblerOptions, } -impl Default for LaplaceSingleLayerPotentialAssembler { - fn default() -> Self { +impl> SingleLayerPotentialAssembler { + /// Create a new single layer potential assembler + pub fn new(kernel: K) -> Self { Self { - kernel: Laplace3dKernel::::new(), + kernel, options: PotentialAssemblerOptions::default(), } } } - -impl PotentialAssembler for LaplaceSingleLayerPotentialAssembler { - const DERIV_SIZE: usize = 1; - type T = T; - - fn options(&self) -> &PotentialAssemblerOptions { - &self.options - } - fn options_mut(&mut self) -> &mut PotentialAssemblerOptions { - &mut self.options - } - - unsafe fn kernel_value( - &self, - k: &RlstArray, - _normals: &RlstArray, - index: usize, - point_index: usize, - ) -> T { - *k.get_unchecked([0, index, point_index]) - } - - fn kernel_assemble_st( - &self, - sources: &[::Real], - targets: &[::Real], - result: &mut [Self::T], - ) { - self.kernel - .assemble_st(GreenKernelEvalType::Value, sources, targets, result); +impl SingleLayerPotentialAssembler> { + /// Create a new Laplace single layer potential assembler + pub fn new_laplace() -> Self { + Self::new(Laplace3dKernel::::new()) } } - -/// Assembler for a Helmholtz single layer potential operator -pub struct HelmholtzSingleLayerPotentialAssembler + MatrixInverse> { - kernel: Helmholtz3dKernel, - options: PotentialAssemblerOptions, -} -impl + MatrixInverse> HelmholtzSingleLayerPotentialAssembler { - /// Create a new assembler - pub fn new(wavenumber: T::Real) -> Self { - Self { - kernel: Helmholtz3dKernel::::new(wavenumber), - options: PotentialAssemblerOptions::default(), - } +impl + MatrixInverse> + SingleLayerPotentialAssembler> +{ + /// Create a new Helmholtz single layer potential assembler + pub fn new_helmholtz(wavenumber: T::Real) -> Self { + Self::new(Helmholtz3dKernel::::new(wavenumber)) } } -impl + MatrixInverse> PotentialAssembler - for HelmholtzSingleLayerPotentialAssembler +impl> PotentialAssembler + for SingleLayerPotentialAssembler { const DERIV_SIZE: usize = 1; type T = T; diff --git a/tests/compare_to_bempp_cl.rs b/tests/compare_to_bempp_cl.rs index fd709216..308ec67d 100644 --- a/tests/compare_to_bempp_cl.rs +++ b/tests/compare_to_bempp_cl.rs @@ -20,7 +20,7 @@ fn test_laplace_single_layer_dp0_dp0() { let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let a = boundary::LaplaceSingleLayerAssembler::::default(); + let a = boundary::SingleLayerAssembler::::new_laplace(); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -43,7 +43,7 @@ fn test_laplace_double_layer_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let a = boundary::LaplaceDoubleLayerAssembler::::default(); + let a = boundary::DoubleLayerAssembler::::new_laplace(); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -65,7 +65,7 @@ fn test_laplace_adjoint_double_layer_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let a = boundary::LaplaceAdjointDoubleLayerAssembler::::default(); + let a = boundary::AdjointDoubleLayerAssembler::::new_laplace(); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -88,7 +88,7 @@ fn test_laplace_hypersingular_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let a = boundary::LaplaceHypersingularAssembler::::default(); + let a = boundary::HypersingularAssembler::::new_laplace(); a.assemble_into_dense(&mut matrix, &space, &space); for i in 0..ndofs { @@ -107,7 +107,7 @@ fn test_laplace_hypersingular_p1_p1() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(f64, [ndofs, ndofs]); - let a = boundary::LaplaceHypersingularAssembler::::default(); + let a = boundary::HypersingularAssembler::::new_laplace(); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -136,7 +136,7 @@ fn test_helmholtz_single_layer_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(c64, [ndofs, ndofs]); - let a = boundary::HelmholtzSingleLayerAssembler::::new(3.0); + let a = boundary::SingleLayerAssembler::::new_helmholtz(3.0); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -159,7 +159,7 @@ fn test_helmholtz_double_layer_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(c64, [ndofs, ndofs]); - let a = boundary::HelmholtzDoubleLayerAssembler::::new(3.0); + let a = boundary::DoubleLayerAssembler::::new_helmholtz(3.0); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -181,7 +181,7 @@ fn test_helmholtz_adjoint_double_layer_dp0_dp0() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(c64, [ndofs, ndofs]); - let a = boundary::HelmholtzAdjointDoubleLayerAssembler::::new(3.0); + let a = boundary::AdjointDoubleLayerAssembler::::new_helmholtz(3.0); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -204,7 +204,7 @@ fn test_helmholtz_hypersingular_p1_p1() { let ndofs = space.global_size(); let mut matrix = rlst_dynamic_array2!(c64, [ndofs, ndofs]); - let a = boundary::HelmholtzHypersingularAssembler::::new(3.0); + let a = boundary::HypersingularAssembler::::new_helmholtz(3.0); a.assemble_into_dense(&mut matrix, &space, &space); // Compare to result from bempp-cl @@ -239,7 +239,7 @@ fn test_laplace_single_layer_potential_dp0() { *points.get_mut([1, 1]).unwrap() = 2.0; *points.get_mut([2, 2]).unwrap() = 2.0; - let a = potential::LaplaceSingleLayerPotentialAssembler::::default(); + let a = potential::SingleLayerPotentialAssembler::::new_laplace(); a.assemble_into_dense(&mut matrix, &space, &points); // Compare to result from bempp-cl @@ -267,7 +267,7 @@ fn test_helmholtz_single_layer_potential_dp0() { *points.get_mut([1, 1]).unwrap() = 2.0; *points.get_mut([2, 2]).unwrap() = 2.0; - let a = potential::HelmholtzSingleLayerPotentialAssembler::::new(3.0); + let a = potential::SingleLayerPotentialAssembler::::new_helmholtz(3.0); a.assemble_into_dense(&mut matrix, &space, &points); // Compare to result from bempp-cl @@ -295,7 +295,7 @@ fn test_laplace_double_layer_potential_dp0() { *points.get_mut([1, 1]).unwrap() = 2.0; *points.get_mut([2, 2]).unwrap() = 2.0; - let a = potential::LaplaceDoubleLayerPotentialAssembler::::default(); + let a = potential::DoubleLayerPotentialAssembler::::new_laplace(); a.assemble_into_dense(&mut matrix, &space, &points); // Compare to result from bempp-cl @@ -323,7 +323,7 @@ fn test_helmholtz_double_layer_potential_dp0() { *points.get_mut([1, 1]).unwrap() = 2.0; *points.get_mut([2, 2]).unwrap() = 2.0; - let a = potential::HelmholtzDoubleLayerPotentialAssembler::::new(3.0); + let a = potential::DoubleLayerPotentialAssembler::::new_helmholtz(3.0); a.assemble_into_dense(&mut matrix, &space, &points); // Compare to result from bempp-cl From c8e80ce63e66456b616ed577cc01369c1bdfae42 Mon Sep 17 00:00:00 2001 From: Matthew Scroggs Date: Mon, 19 Aug 2024 16:53:29 +0100 Subject: [PATCH 5/5] clippy --- src/assembly/boundary/hypersingular.rs | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/assembly/boundary/hypersingular.rs b/src/assembly/boundary/hypersingular.rs index 468cb9aa..fbf9ae2a 100644 --- a/src/assembly/boundary/hypersingular.rs +++ b/src/assembly/boundary/hypersingular.rs @@ -183,7 +183,7 @@ impl<'a, T: RlstScalar + MatrixInverse> BoundaryAssembler const TABLE_DERIVS: usize = 1; type T = T; fn options(&self) -> &BoundaryAssemblerOptions { - &self.options + self.options } fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { panic!("Cannot get mutable options") @@ -269,7 +269,7 @@ impl<'a, T: RlstScalar + MatrixInverse> BoundaryAssembler const TABLE_DERIVS: usize = 0; type T = T; fn options(&self) -> &BoundaryAssemblerOptions { - &self.options + self.options } fn options_mut(&mut self) -> &mut BoundaryAssemblerOptions { panic!("Cannot get mutable options")