numerical-integration 0.1.2.2 → 0.1.2.3

raw patch · 5 files changed

+19/−15 lines, 5 filesPVP: major bump suggested

API removals or changes: PVP suggests a major version bump

API changes (from Hackage documentation)

Files

• CHANGELOG.md +5/−0
• README.md +5/−8
• cpp/integration.cpp +1/−2
• numerical-integration.cabal +1/−1
• src/Numerical/Integration.hs +7/−4

@@ -27,3 +27,8 @@ ## 0.1.2.2 - 2023-09-19
 
 - Updated the cabal file so that the Github actions for Mac work with a recent version of GHC.
+
+
+## 0.1.2.3 - 2023-09-20
+
+- Added the choice of the desired absolute accuracy (it was fixed to 0 before).

@@ -5,12 +5,12 @@ 
 ___
 
-***Example.*** Integrate x² between 0 and 1 with desired relative error 1e-5 and 
-using 200 subdivisions. Exact value: 1/3.
+***Example.*** Integrate x² between 0 and 1 with desired absolute error 0, 
+desired relative error 1e-5 and using 200 subdivisions. Exact value: 1/3.
 
 ```haskell
 example :: IO IntegralResult -- value, error estimate, error code
-example = integration (\x -> x*x) 0 1 1e-5 200
+example = integration (\x -> x*x) 0 1 0.0 1e-5 200
 -- IntegralResult {
 --   _value = 0.3333333333333334, 
 --   _error = 3.7007434154171895e-15, 
@@ -43,7 +43,7 @@ 
 ```haskell
 intgr :: IO IntegralResult 
-intgr = integration (\t -> 5 * cos(2*t/(1-t)) / (25*(1-t)**2 + t**2)) 0 1 1e-4 100000
+intgr = integration (\t -> 5 * cos(2*t/(1-t)) / (25*(1-t)**2 + t**2)) 0 1 0.0 1e-4 100000
 -- IntegralResult {
 --    _value = 7.131328051415349e-5, 
 --    _error = 4.991435083852171e-7, 
@@ -56,7 +56,7 @@ 
 ```haskell
 intgr :: IO IntegralResult 
-intgr = integration (\t -> 5 * cos(2*t/(1-t)) / (25*(1-t)**2 + t**2)) 0 1 1e-4 250000
+intgr = integration (\t -> 5 * cos(2*t/(1-t)) / (25*(1-t)**2 + t**2)) 0 1 0.0 1e-4 250000
 -- IntegralResult {
 --   _value = 7.131328051415349e-5, 
 --   _error = 4.991435083852171e-7, 
@@ -65,9 +65,6 @@ ```
 
 The result is the same!
-
-Note that in the C++ code of 'numerical-integration', I fixed the desired 
-absolute error to 0. I'm not sure it is a good idea.
 
 ___
 

@@ -22,6 +22,7 @@ double integration(double f(double),
                    double lower,
                    double upper,
+                   double absError,
                    double relError,
                    int subdiv,
                    double* errorEstimate,
@@ -34,8 +35,6 @@   // Define a quadrature rule.
   Eigen::Integrator<double>::QuadratureRule rule =
       Eigen::Integrator<double>::GaussKronrod201;
-  // Define the desired absolute error.
-  double absError = 0.0;
   // Integrate.
   double result = integrator.quadratureAdaptive(integrand, lower, upper,
                                                 absError, relError, rule);

@@ -1,6 +1,6 @@ cabal-version:       2.2
 name:                numerical-integration
-version:             0.1.2.2
+version:             0.1.2.3
 synopsis:            Numerical integration.
 description:         One-dimensional numerical integration using the 'NumericalIntegration' C++ library.
 homepage:            https://github.com/stla/numerical-integration#readme

@@ -36,26 +36,29 @@     :: (CDouble -> CDouble) -> IO (FunPtr (CDouble -> CDouble))
 
 foreign import ccall safe "integration" c_integration
-    :: FunPtr (CDouble -> CDouble) -> CDouble -> CDouble -> CDouble -> CInt
-    -> Ptr CDouble -> Ptr CInt -> IO CDouble
+    :: FunPtr (CDouble -> CDouble) -> CDouble -> CDouble -> CDouble
+    -> CDouble -> CInt -> Ptr CDouble -> Ptr CInt -> IO CDouble
 
 -- | Numerical integration.
 integration :: (CDouble -> CDouble)   -- ^ integrand
             -> Double                 -- ^ lower bound
             -> Double                 -- ^ upper bound
+            -> Double                 -- ^ desired absolute error
             -> Double                 -- ^ desired relative error
             -> Int                    -- ^ number of subdivisions
             -> IO IntegralResult      -- ^ value, error estimate, error code
-integration f lower upper relError subdiv = do
+integration f lower upper absError relError subdiv = do
   errorEstimatePtr <- mallocBytes (sizeOf (0 :: CDouble))
   errorCodePtr <- mallocBytes (sizeOf (0 :: CInt))
   fPtr <- funPtr f
   let lower' = double2Cdouble lower
       upper' = double2Cdouble upper 
+      absError' = double2Cdouble absError
       relError' = double2Cdouble relError
       subdiv' = fromIntegral subdiv
   result <-
-    c_integration fPtr lower' upper' relError' subdiv' errorEstimatePtr errorCodePtr
+    c_integration fPtr lower' upper' absError' relError' 
+    subdiv' errorEstimatePtr errorCodePtr
   let result' = if isNaN result 
       then nanDouble 
       else realToFrac result