Skip to content

Commit ca99065

Browse files
committed
Added helmholtz stubs and tracing
1 parent 6da0b19 commit ca99065

File tree

9 files changed

+126
-109
lines changed

9 files changed

+126
-109
lines changed

Cargo.toml

+3
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,7 @@
11
[features]
22
strict = []
33
default = []
4+
enable_tracing = ["rlst/enable_tracing", "dep:log"]
45

56
[package]
67
name = "bempp"
@@ -44,12 +45,14 @@ kifmm = { git = "https://github.com/bempp/kifmm.git", features = ["mpi"] }
4445
bempp-distributed-tools = { git = "https://github.com/bempp/distributed_tools.git"}
4546
rand = "0.8"
4647
rand_chacha = "0.3"
48+
log = { version = "0.4", optional = true}
4749

4850
[dev-dependencies]
4951
approx = "0.5"
5052
cauchy = "0.4.*"
5153
criterion = { version = "0.5.*", features = ["html_reports"] }
5254
# kifmm = { version = "1.0" }
55+
env_logger = "0.11"
5356

5457
[build-dependencies]
5558
cbindgen = "0.27.0"

examples/laplace_evaluator.rs

+4
Original file line numberDiff line numberDiff line change
@@ -62,6 +62,10 @@ fn main() {
6262
let world = universe.world();
6363
let rank = world.rank();
6464

65+
if rank == 0 {
66+
env_logger::init();
67+
}
68+
6569
let refinement_level = 5;
6670

6771
let grid = bempp::shapes::regular_sphere::<f64, _>(refinement_level, 1, &world);

src/boundary_assemblers.rs

+4-2
Original file line numberDiff line numberDiff line change
@@ -27,8 +27,9 @@ use ndgrid::traits::{Entity, Grid, Topology};
2727
use ndgrid::types::Ownership;
2828
use rayon::prelude::*;
2929
use rlst::{
30-
rlst_dynamic_array2, rlst_dynamic_array4, DefaultIterator, DistributedCsrMatrix, DynamicArray,
31-
MatrixInverse, RandomAccessMut, RawAccess, RawAccessMut, RlstScalar, Shape,
30+
measure_duration, rlst_dynamic_array2, rlst_dynamic_array4, DefaultIterator,
31+
DistributedCsrMatrix, DynamicArray, MatrixInverse, RandomAccessMut, RawAccess, RawAccessMut,
32+
RlstScalar, Shape,
3233
};
3334
use std::collections::HashMap;
3435

@@ -142,6 +143,7 @@ impl<'o, T: RlstScalar + MatrixInverse, Integrand: BoundaryIntegrand<T = T>, K:
142143
BoundaryAssembler<'o, T, Integrand, K>
143144
{
144145
/// Assemble the singular part into a CSR matrix.
146+
#[measure_duration(id = "assemble_singular")]
145147
pub fn assemble_singular<'a, C: Communicator, Space: FunctionSpaceTrait<T = T, C = C>>(
146148
&self,
147149
trial_space: &'a Space,

src/evaluator_tools.rs

+4
Original file line numberDiff line numberDiff line change
@@ -14,6 +14,7 @@ use ndgrid::traits::{Entity, GeometryMap, Grid, ParallelGrid, Topology};
1414

1515
use num::Zero;
1616
use rlst::{
17+
measure_duration,
1718
operator::{
1819
interface::{
1920
distributed_sparse_operator::DistributedCsrMatrixOperatorImpl,
@@ -33,6 +34,7 @@ use crate::function::{FunctionSpaceTrait, LocalFunctionSpaceTrait};
3334
/// Create a linear operator from the map of a basis to points. The points are sorted by global
3435
/// index of the corresponding cells in the support. The output space is the space obtained through
3536
/// the owned support cells on each process.
37+
#[measure_duration(id = "space_to_point_map")]
3638
pub fn space_to_point_map<
3739
'a,
3840
T: RlstScalar + Equivalence,
@@ -229,6 +231,7 @@ pub struct NeighbourEvaluator<'a, Space: FunctionSpaceTrait, K: Kernel<T = Space
229231

230232
impl<'a, K: Kernel<T = Space::T>, Space: FunctionSpaceTrait> NeighbourEvaluator<'a, Space, K> {
231233
/// Create a new neighbour evaluator.
234+
#[measure_duration(id = "neighbour_evaluator_from_spaces_and_kernel")]
232235
pub fn from_spaces_and_kernel(
233236
source_space: &'a Space,
234237
target_space: &'a Space,
@@ -387,6 +390,7 @@ where
387390
Space::LocalGrid: Sync,
388391
for<'a> <Space::LocalGrid as Grid>::GeometryMap<'a>: Sync,
389392
{
393+
#[measure_duration(id = "neighbour_evaluator_apply")]
390394
fn apply_extended<
391395
ContainerIn: rlst::ElementContainer<E = <Self::Domain as rlst::LinearSpace>::E>,
392396
ContainerOut: rlst::ElementContainerMut<E = <Self::Range as rlst::LinearSpace>::E>,

src/greens_function_evaluators/kifmm_evaluator.rs

+3
Original file line numberDiff line numberDiff line change
@@ -24,6 +24,7 @@ use kifmm::{
2424
use mpi::traits::{Communicator, Equivalence};
2525
use num::Float;
2626
use rlst::{
27+
measure_duration,
2728
operator::{interface::DistributedArrayVectorSpace, Operator},
2829
rlst_dynamic_array1, AsApply, Element, IndexLayout, MatrixSvd, OperatorBase, RawAccess,
2930
RawAccessMut, RlstScalar,
@@ -73,6 +74,7 @@ where
7374
KiFmmMulti<T, Laplace3dKernel<T>, FftFieldTranslation<T>>: SourceToTargetTranslationMetadata,
7475
{
7576
/// Instantiate a new KiFmm Evaluator.
77+
#[measure_duration(id = "kifmm_new")]
7678
pub fn new(
7779
sources: &[T::Real],
7880
targets: &[T::Real],
@@ -217,6 +219,7 @@ where
217219
KiFmmMulti<T, Laplace3dKernel<T>, FftFieldTranslation<T>>: DataAccessMulti<Scalar = T>,
218220
KiFmmMulti<T, Laplace3dKernel<T>, FftFieldTranslation<T>>: GhostExchange,
219221
{
222+
#[measure_duration(id = "kifmm_apply_extended")]
220223
fn apply_extended<
221224
ContainerIn: rlst::ElementContainer<E = <Self::Domain as rlst::LinearSpace>::E>,
222225
ContainerOut: rlst::ElementContainerMut<E = <Self::Range as rlst::LinearSpace>::E>,

src/helmholtz.rs

+2-107
Original file line numberDiff line numberDiff line change
@@ -1,108 +1,3 @@
11
//! Helmholtz operators
2-
3-
/// Assemblers for Helmholtz problems
4-
pub mod assembler {
5-
use green_kernels::{helmholtz_3d::Helmholtz3dKernel, types::GreenKernelEvalType};
6-
use rlst::{MatrixInverse, RlstScalar};
7-
8-
use crate::boundary_assemblers::{
9-
helpers::KernelEvaluator,
10-
integrands::{
11-
AdjointDoubleLayerBoundaryIntegrand, BoundaryIntegrandSum,
12-
BoundaryIntegrandTimesScalar, DoubleLayerBoundaryIntegrand,
13-
HypersingularCurlCurlBoundaryIntegrand, HypersingularNormalNormalBoundaryIntegrand,
14-
SingleLayerBoundaryIntegrand,
15-
},
16-
BoundaryAssembler, BoundaryAssemblerOptions,
17-
};
18-
19-
/// Helmholtz single layer assembler type.
20-
pub type SingleLayer3dAssembler<'o, T> =
21-
BoundaryAssembler<'o, T, SingleLayerBoundaryIntegrand<T>, Helmholtz3dKernel<T>>;
22-
23-
/// Helmholtz double layer assembler type.
24-
pub type DoubleLayer3dAssembler<'o, T> =
25-
BoundaryAssembler<'o, T, DoubleLayerBoundaryIntegrand<T>, Helmholtz3dKernel<T>>;
26-
27-
/// Helmholtz adjoint double layer assembler type.
28-
pub type AdjointDoubleLayer3dAssembler<'o, T> =
29-
BoundaryAssembler<'o, T, AdjointDoubleLayerBoundaryIntegrand<T>, Helmholtz3dKernel<T>>;
30-
31-
/// Helmholtz hypersingular double layer assembler type.
32-
pub type Hypersingular3dAssembler<'o, T> = BoundaryAssembler<
33-
'o,
34-
T,
35-
BoundaryIntegrandSum<
36-
T,
37-
HypersingularCurlCurlBoundaryIntegrand<T>,
38-
BoundaryIntegrandTimesScalar<T, HypersingularNormalNormalBoundaryIntegrand<T>>,
39-
>,
40-
Helmholtz3dKernel<T>,
41-
>;
42-
43-
/// Assembler for the Helmholtz single layer operator.
44-
pub fn single_layer<T: RlstScalar<Complex = T> + MatrixInverse>(
45-
wavenumber: T::Real,
46-
options: &BoundaryAssemblerOptions,
47-
) -> SingleLayer3dAssembler<T> {
48-
let kernel = KernelEvaluator::new(
49-
Helmholtz3dKernel::new(wavenumber),
50-
GreenKernelEvalType::Value,
51-
);
52-
53-
BoundaryAssembler::new(SingleLayerBoundaryIntegrand::new(), kernel, options, 1, 0)
54-
}
55-
56-
/// Assembler for the Helmholtz double layer operator.
57-
pub fn double_layer<T: RlstScalar<Complex = T> + MatrixInverse>(
58-
wavenumber: T::Real,
59-
options: &BoundaryAssemblerOptions,
60-
) -> DoubleLayer3dAssembler<T> {
61-
let kernel = KernelEvaluator::new(
62-
Helmholtz3dKernel::new(wavenumber),
63-
GreenKernelEvalType::ValueDeriv,
64-
);
65-
66-
BoundaryAssembler::new(DoubleLayerBoundaryIntegrand::new(), kernel, options, 4, 0)
67-
}
68-
69-
/// Assembler for the Helmholtz adjoint double layer operator.
70-
pub fn adjoint_double_layer<T: RlstScalar<Complex = T> + MatrixInverse>(
71-
wavenumber: T::Real,
72-
options: &BoundaryAssemblerOptions,
73-
) -> AdjointDoubleLayer3dAssembler<T> {
74-
let kernel = KernelEvaluator::new(
75-
Helmholtz3dKernel::new(wavenumber),
76-
GreenKernelEvalType::ValueDeriv,
77-
);
78-
79-
BoundaryAssembler::new(
80-
AdjointDoubleLayerBoundaryIntegrand::new(),
81-
kernel,
82-
options,
83-
4,
84-
0,
85-
)
86-
}
87-
88-
/// Assembler for the Helmholtz hypersingular operator.
89-
pub fn hypersingular<T: RlstScalar<Complex = T> + MatrixInverse>(
90-
wavenumber: T::Real,
91-
options: &BoundaryAssemblerOptions,
92-
) -> Hypersingular3dAssembler<T> {
93-
let kernel = KernelEvaluator::new(
94-
Helmholtz3dKernel::new(wavenumber),
95-
GreenKernelEvalType::ValueDeriv,
96-
);
97-
98-
let integrand = BoundaryIntegrandSum::new(
99-
HypersingularCurlCurlBoundaryIntegrand::new(),
100-
BoundaryIntegrandTimesScalar::new(
101-
num::cast::<T::Real, T>(-wavenumber.powi(2)).unwrap(),
102-
HypersingularNormalNormalBoundaryIntegrand::new(),
103-
),
104-
);
105-
106-
BoundaryAssembler::new(integrand, kernel, options, 4, 1)
107-
}
108-
}
2+
pub mod assembler;
3+
pub mod evaluator;

src/helmholtz/assembler.rs

+104
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,104 @@
1+
//! Implementation of Helmholtz assemblers
2+
3+
use green_kernels::{helmholtz_3d::Helmholtz3dKernel, types::GreenKernelEvalType};
4+
use rlst::{MatrixInverse, RlstScalar};
5+
6+
use crate::boundary_assemblers::{
7+
helpers::KernelEvaluator,
8+
integrands::{
9+
AdjointDoubleLayerBoundaryIntegrand, BoundaryIntegrandSum, BoundaryIntegrandTimesScalar,
10+
DoubleLayerBoundaryIntegrand, HypersingularCurlCurlBoundaryIntegrand,
11+
HypersingularNormalNormalBoundaryIntegrand, SingleLayerBoundaryIntegrand,
12+
},
13+
BoundaryAssembler, BoundaryAssemblerOptions,
14+
};
15+
16+
/// Helmholtz single layer assembler type.
17+
pub type SingleLayer3dAssembler<'o, T> =
18+
BoundaryAssembler<'o, T, SingleLayerBoundaryIntegrand<T>, Helmholtz3dKernel<T>>;
19+
20+
/// Helmholtz double layer assembler type.
21+
pub type DoubleLayer3dAssembler<'o, T> =
22+
BoundaryAssembler<'o, T, DoubleLayerBoundaryIntegrand<T>, Helmholtz3dKernel<T>>;
23+
24+
/// Helmholtz adjoint double layer assembler type.
25+
pub type AdjointDoubleLayer3dAssembler<'o, T> =
26+
BoundaryAssembler<'o, T, AdjointDoubleLayerBoundaryIntegrand<T>, Helmholtz3dKernel<T>>;
27+
28+
/// Helmholtz hypersingular double layer assembler type.
29+
pub type Hypersingular3dAssembler<'o, T> = BoundaryAssembler<
30+
'o,
31+
T,
32+
BoundaryIntegrandSum<
33+
T,
34+
HypersingularCurlCurlBoundaryIntegrand<T>,
35+
BoundaryIntegrandTimesScalar<T, HypersingularNormalNormalBoundaryIntegrand<T>>,
36+
>,
37+
Helmholtz3dKernel<T>,
38+
>;
39+
40+
/// Assembler for the Helmholtz single layer operator.
41+
pub fn single_layer<T: RlstScalar<Complex = T> + MatrixInverse>(
42+
wavenumber: T::Real,
43+
options: &BoundaryAssemblerOptions,
44+
) -> SingleLayer3dAssembler<T> {
45+
let kernel = KernelEvaluator::new(
46+
Helmholtz3dKernel::new(wavenumber),
47+
GreenKernelEvalType::Value,
48+
);
49+
50+
BoundaryAssembler::new(SingleLayerBoundaryIntegrand::new(), kernel, options, 1, 0)
51+
}
52+
53+
/// Assembler for the Helmholtz double layer operator.
54+
pub fn double_layer<T: RlstScalar<Complex = T> + MatrixInverse>(
55+
wavenumber: T::Real,
56+
options: &BoundaryAssemblerOptions,
57+
) -> DoubleLayer3dAssembler<T> {
58+
let kernel = KernelEvaluator::new(
59+
Helmholtz3dKernel::new(wavenumber),
60+
GreenKernelEvalType::ValueDeriv,
61+
);
62+
63+
BoundaryAssembler::new(DoubleLayerBoundaryIntegrand::new(), kernel, options, 4, 0)
64+
}
65+
66+
/// Assembler for the Helmholtz adjoint double layer operator.
67+
pub fn adjoint_double_layer<T: RlstScalar<Complex = T> + MatrixInverse>(
68+
wavenumber: T::Real,
69+
options: &BoundaryAssemblerOptions,
70+
) -> AdjointDoubleLayer3dAssembler<T> {
71+
let kernel = KernelEvaluator::new(
72+
Helmholtz3dKernel::new(wavenumber),
73+
GreenKernelEvalType::ValueDeriv,
74+
);
75+
76+
BoundaryAssembler::new(
77+
AdjointDoubleLayerBoundaryIntegrand::new(),
78+
kernel,
79+
options,
80+
4,
81+
0,
82+
)
83+
}
84+
85+
/// Assembler for the Helmholtz hypersingular operator.
86+
pub fn hypersingular<T: RlstScalar<Complex = T> + MatrixInverse>(
87+
wavenumber: T::Real,
88+
options: &BoundaryAssemblerOptions,
89+
) -> Hypersingular3dAssembler<T> {
90+
let kernel = KernelEvaluator::new(
91+
Helmholtz3dKernel::new(wavenumber),
92+
GreenKernelEvalType::ValueDeriv,
93+
);
94+
95+
let integrand = BoundaryIntegrandSum::new(
96+
HypersingularCurlCurlBoundaryIntegrand::new(),
97+
BoundaryIntegrandTimesScalar::new(
98+
num::cast::<T::Real, T>(-wavenumber.powi(2)).unwrap(),
99+
HypersingularNormalNormalBoundaryIntegrand::new(),
100+
),
101+
);
102+
103+
BoundaryAssembler::new(integrand, kernel, options, 4, 1)
104+
}

src/helmholtz/evaluator.rs

+1
Original file line numberDiff line numberDiff line change
@@ -0,0 +1 @@
1+
//! Helmholtz evaluators

src/laplace/evaluator.rs

+1
Original file line numberDiff line numberDiff line change
@@ -12,6 +12,7 @@ use crate::{
1212
};
1313

1414
/// Implement a single layer operator
15+
#[measure_duration(id = "laplace_single_layer_evaluate")]
1516
#[allow(clippy::type_complexity)]
1617
pub fn single_layer<
1718
'a,

0 commit comments

Comments
 (0)