← Back to branch summary

~ubuntu-branches/ubuntu/vivid/haskell-hmatrix/vivid

« back to all changes in this revision

Viewing changes to lib/Numeric/GSL/Root.hs

  • Committer: Package Import Robot
  • Author(s): Denis Laxalde
  • Date: 2013-07-06 15:37:50 UTC
  • mfrom: (1.1.1)
  • Revision ID: package-import@ubuntu.com-20130706153750-wxxplc788jedqvv5
Tags: 0.15.0.0-1
* New upstream release.
* Make it clear in copyright that the license is GPL-3 (as stated by the
  author in <https://github.com/albertoruiz/hmatrix/issues/45>) although
  there is no proper license file yet.

Show diffs side-by-side

added added

removed removed

Lines of Context:
lib/Numeric/GSL/Root.hs
45
45
-----------------------------------------------------------------------------
46
46
 
47
47
module Numeric.GSL.Root (
 
48
    uniRoot, UniRootMethod(..),
 
49
    uniRootJ, UniRootMethodJ(..),
48
50
    root, RootMethod(..),
49
51
    rootJ, RootMethodJ(..),
50
52
) where
58
60
 
59
61
-------------------------------------------------------------------------
60
62
 
 
63
data UniRootMethod = Bisection
 
64
                   | FalsePos
 
65
                   | Brent
 
66
                   deriving (Enum, Eq, Show, Bounded)
 
67
 
 
68
uniRoot :: UniRootMethod
 
69
        -> Double
 
70
        -> Int
 
71
        -> (Double -> Double)
 
72
        -> Double
 
73
        -> Double
 
74
        -> (Double, Matrix Double)
 
75
uniRoot method epsrel maxit fun xl xu = uniRootGen (fi (fromEnum method)) fun xl xu epsrel maxit
 
76
 
 
77
uniRootGen m f xl xu epsrel maxit = unsafePerformIO $ do
 
78
    fp <- mkDoublefun f
 
79
    rawpath <- createMIO maxit 4
 
80
                         (c_root m fp epsrel (fi maxit) xl xu)
 
81
                         "root"
 
82
    let it = round (rawpath @@> (maxit-1,0))
 
83
        path = takeRows it rawpath
 
84
        [sol] = toLists $ dropRows (it-1) path
 
85
    freeHaskellFunPtr fp
 
86
    return (sol !! 1, path)
 
87
 
 
88
foreign import ccall safe "root"
 
89
    c_root:: CInt -> FunPtr (Double -> Double) -> Double -> CInt -> Double -> Double -> TM
 
90
 
 
91
-------------------------------------------------------------------------
 
92
data UniRootMethodJ = UNewton
 
93
                    | Secant
 
94
                    | Steffenson
 
95
                    deriving (Enum, Eq, Show, Bounded)
 
96
 
 
97
uniRootJ :: UniRootMethodJ
 
98
        -> Double
 
99
        -> Int
 
100
        -> (Double -> Double)
 
101
        -> (Double -> Double)
 
102
        -> Double
 
103
        -> (Double, Matrix Double)
 
104
uniRootJ method epsrel maxit fun dfun x = uniRootJGen (fi (fromEnum method)) fun
 
105
    dfun x epsrel maxit
 
106
 
 
107
uniRootJGen m f df x epsrel maxit = unsafePerformIO $ do
 
108
    fp <- mkDoublefun f
 
109
    dfp <- mkDoublefun df
 
110
    rawpath <- createMIO maxit 2
 
111
                         (c_rootj m fp dfp epsrel (fi maxit) x)
 
112
                         "rootj"
 
113
    let it = round (rawpath @@> (maxit-1,0))
 
114
        path = takeRows it rawpath
 
115
        [sol] = toLists $ dropRows (it-1) path
 
116
    freeHaskellFunPtr fp
 
117
    return (sol !! 1, path)
 
118
 
 
119
foreign import ccall safe "rootj"
 
120
    c_rootj :: CInt -> FunPtr (Double -> Double) -> FunPtr (Double -> Double)
 
121
            -> Double -> CInt -> Double -> TM
 
122
 
 
123
-------------------------------------------------------------------------
 
124
 
61
125
data RootMethod = Hybrids
62
126
                | Hybrid
63
127
                | DNewton
82
146
    fp <- mkVecVecfun (aux_vTov (checkdim1 n . fromList . f . toList))
83
147
    rawpath <- vec xiv $ \xiv' ->
84
148
                   createMIO maxit (2*n+1)
85
 
                         (c_root m fp epsabs (fi maxit) // xiv')
86
 
                         "root"
 
149
                         (c_multiroot m fp epsabs (fi maxit) // xiv')
 
150
                         "multiroot"
87
151
    let it = round (rawpath @@> (maxit-1,0))
88
152
        path = takeRows it rawpath
89
153
        [sol] = toLists $ dropRows (it-1) path
91
155
    return (take n $ drop 1 sol, path)
92
156
 
93
157
 
94
 
foreign import ccall safe "root"
95
 
    c_root:: CInt -> FunPtr TVV -> Double -> CInt -> TVM
 
158
foreign import ccall safe "multiroot"
 
159
    c_multiroot:: CInt -> FunPtr TVV -> Double -> CInt -> TVM
96
160
 
97
161
-------------------------------------------------------------------------
98
162
 
120
184
    jp <- mkVecMatfun (aux_vTom (checkdim2 n . fromLists . jac . toList))
121
185
    rawpath <- vec xiv $ \xiv' ->
122
186
                   createMIO maxit (2*n+1)
123
 
                         (c_rootj m fp jp epsabs (fi maxit) // xiv')
124
 
                         "root"
 
187
                         (c_multirootj m fp jp epsabs (fi maxit) // xiv')
 
188
                         "multiroot"
125
189
    let it = round (rawpath @@> (maxit-1,0))
126
190
        path = takeRows it rawpath
127
191
        [sol] = toLists $ dropRows (it-1) path
129
193
    freeHaskellFunPtr jp
130
194
    return (take n $ drop 1 sol, path)
131
195
 
132
 
 
133
 
foreign import ccall safe "rootj"
134
 
    c_rootj:: CInt -> FunPtr TVV -> FunPtr TVM -> Double -> CInt -> TVM
 
196
foreign import ccall safe "multirootj"
 
197
    c_multirootj:: CInt -> FunPtr TVV -> FunPtr TVM -> Double -> CInt -> TVM
135
198
 
136
199
-------------------------------------------------------
137
200