diff options
author | Alberto Ruiz <aruiz@um.es> | 2007-10-01 15:04:16 +0000 |
---|---|---|
committer | Alberto Ruiz <aruiz@um.es> | 2007-10-01 15:04:16 +0000 |
commit | c99b8fd6e3f8a2fb365ec12baf838f864b118ece (patch) | |
tree | 11b5b8515861fe88d547253ae10c2182d5fadaf2 /lib/Numeric/GSL/Integration.hs | |
parent | 768f08d4134a066d773d56a9c03ae688e3850352 (diff) |
LinearAlgebra and GSL moved to Numeric
Diffstat (limited to 'lib/Numeric/GSL/Integration.hs')
-rw-r--r-- | lib/Numeric/GSL/Integration.hs | 90 |
1 files changed, 90 insertions, 0 deletions
diff --git a/lib/Numeric/GSL/Integration.hs b/lib/Numeric/GSL/Integration.hs new file mode 100644 index 0000000..d756417 --- /dev/null +++ b/lib/Numeric/GSL/Integration.hs | |||
@@ -0,0 +1,90 @@ | |||
1 | {-# OPTIONS #-} | ||
2 | ----------------------------------------------------------------------------- | ||
3 | {- | | ||
4 | Module : Numeric.GSL.Integration | ||
5 | Copyright : (c) Alberto Ruiz 2006 | ||
6 | License : GPL-style | ||
7 | |||
8 | Maintainer : Alberto Ruiz (aruiz at um dot es) | ||
9 | Stability : provisional | ||
10 | Portability : uses ffi | ||
11 | |||
12 | Numerical integration routines. | ||
13 | |||
14 | <http://www.gnu.org/software/gsl/manual/html_node/Numerical-Integration.html#Numerical-Integration> | ||
15 | -} | ||
16 | ----------------------------------------------------------------------------- | ||
17 | |||
18 | module Numeric.GSL.Integration ( | ||
19 | integrateQNG, | ||
20 | integrateQAGS | ||
21 | ) where | ||
22 | |||
23 | import Foreign | ||
24 | import Data.Packed.Internal(mkfun,check,(//)) | ||
25 | |||
26 | |||
27 | -------------------------------------------------------------------- | ||
28 | {- | Numerical integration using /gsl_integration_qags/ (adaptive integration with singularities). For example: | ||
29 | |||
30 | @\> let quad = integrateQAGS 1E-9 1000 | ||
31 | \> let f a x = x**(-0.5) * log (a*x) | ||
32 | \> quad (f 1) 0 1 | ||
33 | (-3.999999999999974,4.871658632055187e-13)@ | ||
34 | |||
35 | -} | ||
36 | |||
37 | integrateQAGS :: Double -- ^ precision (e.g. 1E-9) | ||
38 | -> Int -- ^ size of auxiliary workspace (e.g. 1000) | ||
39 | -> (Double -> Double) -- ^ function to be integrated on the interval (a,b) | ||
40 | -> Double -- ^ a | ||
41 | -> Double -- ^ b | ||
42 | -> (Double, Double) -- ^ result of the integration and error | ||
43 | integrateQAGS prec n f a b = unsafePerformIO $ do | ||
44 | r <- malloc | ||
45 | e <- malloc | ||
46 | fp <- mkfun (\x _ -> f x) | ||
47 | c_integrate_qags fp a b prec n r e // check "integrate_qags" [] | ||
48 | vr <- peek r | ||
49 | ve <- peek e | ||
50 | let result = (vr,ve) | ||
51 | free r | ||
52 | free e | ||
53 | freeHaskellFunPtr fp | ||
54 | return result | ||
55 | |||
56 | foreign import ccall "gsl-aux.h integrate_qags" | ||
57 | c_integrate_qags :: FunPtr (Double-> Ptr() -> Double) -> Double -> Double -> Double -> Int | ||
58 | -> Ptr Double -> Ptr Double -> IO Int | ||
59 | |||
60 | ----------------------------------------------------------------- | ||
61 | {- | Numerical integration using /gsl_integration_qng/ (useful for fast integration of smooth functions). For example: | ||
62 | |||
63 | @\> let quad = integrateQNG 1E-6 | ||
64 | \> quad (\\x -> 4\/(1+x*x)) 0 1 | ||
65 | (3.141592653589793,3.487868498008632e-14)@ | ||
66 | |||
67 | -} | ||
68 | |||
69 | integrateQNG :: Double -- ^ precision (e.g. 1E-9) | ||
70 | -> (Double -> Double) -- ^ function to be integrated on the interval (a,b) | ||
71 | -> Double -- ^ a | ||
72 | -> Double -- ^ b | ||
73 | -> (Double, Double) -- ^ result of the integration and error | ||
74 | integrateQNG prec f a b = unsafePerformIO $ do | ||
75 | r <- malloc | ||
76 | e <- malloc | ||
77 | fp <- mkfun (\x _ -> f x) | ||
78 | c_integrate_qng fp a b prec r e // check "integrate_qng" [] | ||
79 | vr <- peek r | ||
80 | ve <- peek e | ||
81 | let result = (vr,ve) | ||
82 | free r | ||
83 | free e | ||
84 | freeHaskellFunPtr fp | ||
85 | return result | ||
86 | |||
87 | foreign import ccall "gsl-aux.h integrate_qng" | ||
88 | c_integrate_qng :: FunPtr (Double-> Ptr() -> Double) -> Double -> Double -> Double | ||
89 | -> Ptr Double -> Ptr Double -> IO Int | ||
90 | |||