Mercurial > repos > pointsource_algs / file diff

--- a/src/kernels/mollifier.rs	Sun Apr 27 15:03:51 2025 -0500
+++ b/src/kernels/mollifier.rs	Thu Feb 26 11:38:43 2026 -0500
@@ -1,24 +1,17 @@
-
 //! Implementation of the standard mollifier
 
-use rgsl::hypergeometric::hyperg_U;
+use alg_tools::bisection_tree::{Bounds, Constant, GlobalAnalysis, LocalAnalysis, Support};
+use alg_tools::euclidean::Euclidean;
+use alg_tools::loc::Loc;
+use alg_tools::mapping::{Instance, Mapping};
+use alg_tools::maputil::array_init;
+use alg_tools::norms::*;
+use alg_tools::sets::Cube;
+use alg_tools::types::*;
 use float_extras::f64::tgamma as gamma;
 use numeric_literals::replace_float_literals;
+use rgsl::hypergeometric::hyperg_U;
 use serde::Serialize;
-use alg_tools::types::*;
-use alg_tools::euclidean::Euclidean;
-use alg_tools::norms::*;
-use alg_tools::loc::Loc;
-use alg_tools::sets::Cube;
-use alg_tools::bisection_tree::{
-    Support,
-    Constant,
-    Bounds,
-    LocalAnalysis,
-    GlobalAnalysis
-};
-use alg_tools::mapping::{Mapping, Instance};
-use alg_tools::maputil::array_init;
 
 /// Reresentation of the (unnormalised) standard mollifier.
 ///
@@ -29,20 +22,20 @@
 ///         0, & \text{otherwise}.
 ///     \end{cases}
 /// $$</div>
-#[derive(Copy,Clone,Serialize,Debug,Eq,PartialEq)]
-pub struct Mollifier<C : Constant, const N : usize> {
+#[derive(Copy, Clone, Serialize, Debug, Eq, PartialEq)]
+pub struct Mollifier<C: Constant, const N: usize> {
     /// The parameter $ε$ of the mollifier.
-    pub width : C,
+    pub width: C,
 }
 
 #[replace_float_literals(C::Type::cast_from(literal))]
-impl<C : Constant, const N : usize> Mapping<Loc<C::Type, N>> for Mollifier<C, N> {
+impl<C: Constant, const N: usize> Mapping<Loc<N, C::Type>> for Mollifier<C, N> {
     type Codomain = C::Type;
 
     #[inline]
-    fn apply<I : Instance<Loc<C::Type, N>>>(&self, x : I) -> Self::Codomain {
+    fn apply<I: Instance<Loc<N, C::Type>>>(&self, x: I) -> Self::Codomain {
         let ε = self.width.value();
-        let ε2 = ε*ε;
+        let ε2 = ε * ε;
         let n2 = x.eval(|x| x.norm2_squared());
         if n2 < ε2 {
             (n2 / (n2 - ε2)).exp()
@@ -52,27 +45,25 @@
     }
 }
 
-
-impl<'a, C : Constant, const N : usize> Support<C::Type, N> for Mollifier<C, N> {
+impl<'a, C: Constant, const N: usize> Support<N, C::Type> for Mollifier<C, N> {
     #[inline]
-    fn support_hint(&self) -> Cube<C::Type,N> {
+    fn support_hint(&self) -> Cube<N, C::Type> {
         let ε = self.width.value();
         array_init(|| [-ε, ε]).into()
     }
 
     #[inline]
-    fn in_support(&self, x : &Loc<C::Type,N>) -> bool {
+    fn in_support(&self, x: &Loc<N, C::Type>) -> bool {
         x.norm2() < self.width.value()
     }
-    
+
     /*fn fully_in_support(&self, _cube : &Cube<C::Type,N>) -> bool {
         todo!("Not implemented, but not used at the moment")
     }*/
 }
 
 #[replace_float_literals(C::Type::cast_from(literal))]
-impl<'a, C : Constant, const N : usize> GlobalAnalysis<C::Type, Bounds<C::Type>>
-for Mollifier<C, N> {
+impl<'a, C: Constant, const N: usize> GlobalAnalysis<C::Type, Bounds<C::Type>> for Mollifier<C, N> {
     #[inline]
     fn global_analysis(&self) -> Bounds<C::Type> {
         // The function is maximised/minimised where the 2-norm is minimised/maximised.
@@ -80,10 +71,11 @@
     }
 }
 
-impl<'a, C : Constant, const N : usize> LocalAnalysis<C::Type, Bounds<C::Type>, N>
-for Mollifier<C, N> {
+impl<'a, C: Constant, const N: usize> LocalAnalysis<C::Type, Bounds<C::Type>, N>
+    for Mollifier<C, N>
+{
     #[inline]
-    fn local_analysis(&self, cube : &Cube<C::Type, N>) -> Bounds<C::Type> {
+    fn local_analysis(&self, cube: &Cube<N, C::Type>) -> Bounds<C::Type> {
         // The function is maximised/minimised where the 2-norm is minimised/maximised.
         let lower = self.apply(cube.maxnorm_point());
         let upper = self.apply(cube.minnorm_point());
@@ -99,7 +91,7 @@
 /// If `rescaled` is `true`, return the integral of the scaled mollifier that has value one at the
 /// origin.
 #[inline]
-pub fn mollifier_norm1(n_ : usize, rescaled : bool) -> f64 {
+pub fn mollifier_norm1(n_: usize, rescaled: bool) -> f64 {
     assert!(n_ > 0);
     let n = n_ as f64;
     let q = 2.0;
@@ -108,23 +100,25 @@
                /*/ gamma(1.0 + n / p)
                * gamma(1.0 + n / q)*/
                * hyperg_U(1.0 + n / q, 2.0, 1.0);
-    if rescaled { base } else { base / f64::E }
+    if rescaled {
+        base
+    } else {
+        base / f64::E
+    }
 }
 
-impl<'a, C : Constant, const N : usize> Norm<C::Type, L1>
-for Mollifier<C, N> {
+impl<'a, C: Constant, const N: usize> Norm<L1, C::Type> for Mollifier<C, N> {
     #[inline]
-    fn norm(&self, _ : L1) -> C::Type {
+    fn norm(&self, _: L1) -> C::Type {
         let ε = self.width.value();
         C::Type::cast_from(mollifier_norm1(N, true)) * ε.powi(N as i32)
     }
 }
 
 #[replace_float_literals(C::Type::cast_from(literal))]
-impl<'a, C : Constant, const N : usize> Norm<C::Type, Linfinity>
-for Mollifier<C, N> {
+impl<'a, C: Constant, const N: usize> Norm<Linfinity, C::Type> for Mollifier<C, N> {
     #[inline]
-    fn norm(&self, _ : Linfinity) -> C::Type {
+    fn norm(&self, _: Linfinity) -> C::Type {
         1.0
     }
 }