~ubuntu-branches/ubuntu/trusty/haskell-hmatrix/trusty

The vector Broyden-Fletcher-Goldfarb-Shanno (BFGS) algorithm is a quasi-Newton method which builds up an approximation to the second derivatives of the function f using the difference between successive gradient vectors. By combining the first and second derivatives the algorithm is able to take Newton-type steps towards the function minimum, assuming quadratic behavior in that region.

The bfgs2 version of this minimizer is the most efficient version available, and is a faithful implementation of the line minimization scheme described in Fletcher's Practical Methods of Optimization, Algorithms 2.6.2 and 2.6.4. It supercedes the original bfgs routine and requires substantially fewer function and gradient evaluations. The user-supplied tolerance tol corresponds to the parameter \sigma used by Fletcher. A value of 0.1 is recommended for typical use (larger values correspond to less accurate line searches).

The nmsimplex2 version is a new O(N) implementation of the earlier O(N^2) nmsimplex minimiser. It calculates the size of simplex as the rms distance of each vertex from the center rather than the mean distance, which has the advantage of allowing a linear update.

-----------------------------------------------------------------------------

module Numeric.GSL.Minimization (

minimize, minimizeV, MinimizeMethod(..),

minimizeD, minimizeVD, MinimizeMethodD(..),

minimizeNMSimplex,

minimizeConjugateGradient,

minimizeVectorBFGS2

) where

import Data.Packed.Internal

import Data.Packed.Matrix

import Numeric.GSL.Internal

import Foreign.Ptr(Ptr, FunPtr, freeHaskellFunPtr)

import Foreign.C.Types

import System.IO.Unsafe(unsafePerformIO)

------------------------------------------------------------------------

{-# DEPRECATED minimizeNMSimplex "use minimize NMSimplex2 eps maxit sizes f xi" #-}

minimizeNMSimplex f xi szs eps maxit = minimize NMSimplex eps maxit szs f xi

{-# DEPRECATED minimizeConjugateGradient "use minimizeD ConjugateFR eps maxit step tol f g xi" #-}

minimizeConjugateGradient step tol eps maxit f g xi = minimizeD ConjugateFR eps maxit step tol f g xi

{-# DEPRECATED minimizeVectorBFGS2 "use minimizeD VectorBFGS2 eps maxit step tol f g xi" #-}

minimizeVectorBFGS2 step tol eps maxit f g xi = minimizeD VectorBFGS2 eps maxit step tol f g xi

-------------------------------------------------------------------------

data MinimizeMethod = NMSimplex

| NMSimplex2

deriving (Enum,Eq,Show,Bounded)

-- | Minimization without derivatives

minimize :: MinimizeMethod

-> Double -- ^ desired precision of the solution (size test)

-> Int -- ^ maximum number of iterations allowed

-> [Double] -- ^ sizes of the initial search box

-> ([Double] -> Double) -- ^ function to minimize

-> [Double] -- ^ starting point

-> ([Double], Matrix Double) -- ^ solution vector and optimization path

-- | Minimization without derivatives (vector version)

minimizeV :: MinimizeMethod

-> Double -- ^ desired precision of the solution (size test)

100

-> Int -- ^ maximum number of iterations allowed

101

-> Vector Double -- ^ sizes of the initial search box

102

-> (Vector Double -> Double) -- ^ function to minimize

103

-> Vector Double -- ^ starting point

104

-> (Vector Double, Matrix Double) -- ^ solution vector and optimization path

105

106

minimize method eps maxit sz f xi = v2l $ minimizeV method eps maxit (fromList sz) (f.toList) (fromList xi)

107

where v2l (v,m) = (toList v, m)

108

109

ww2 w1 o1 w2 o2 f = w1 o1 $ \a1 -> w2 o2 $ \a2 -> f a1 a2

110

111

minimizeV method eps maxit szv f xiv = unsafePerformIO $ do

112

let n = dim xiv

113

fp <- mkVecfun (iv f)

114

rawpath <- ww2 vec xiv vec szv $ \xiv' szv' ->

115

createMIO maxit (n+3)

116

(c_minimize (fi (fromEnum method)) fp eps (fi maxit) // xiv' // szv')

117

"minimize"

118

let it = round (rawpath @@> (maxit-1,0))

119

path = takeRows it rawpath

120

sol = cdat $ dropColumns 3 $ dropRows (it-1) path

121

freeHaskellFunPtr fp

122

return (sol, path)

123

124

125

foreign import ccall safe "gsl-aux.h minimize"

126

c_minimize:: CInt -> FunPtr (CInt -> Ptr Double -> Double) -> Double -> CInt -> TVVM

127

128

----------------------------------------------------------------------------------

129

130

131

data MinimizeMethodD = ConjugateFR

132

| ConjugatePR

133

| VectorBFGS

134

| VectorBFGS2

135

| SteepestDescent

136

deriving (Enum,Eq,Show,Bounded)

137

138

-- | Minimization with derivatives.

139

minimizeD :: MinimizeMethodD

140

-> Double -- ^ desired precision of the solution (gradient test)

141

-> Int -- ^ maximum number of iterations allowed

142

-> Double -- ^ size of the first trial step

143

-> Double -- ^ tol (precise meaning depends on method)

144

-> ([Double] -> Double) -- ^ function to minimize

145

-> ([Double] -> [Double]) -- ^ gradient

146

-> [Double] -- ^ starting point

147

-> ([Double], Matrix Double) -- ^ solution vector and optimization path

148

149

-- | Minimization with derivatives (vector version)

150

minimizeVD :: MinimizeMethodD

151

-> Double -- ^ desired precision of the solution (gradient test)

152

-> Int -- ^ maximum number of iterations allowed

153

-> Double -- ^ size of the first trial step

154

-> Double -- ^ tol (precise meaning depends on method)

155

-> (Vector Double -> Double) -- ^ function to minimize

156

-> (Vector Double -> Vector Double) -- ^ gradient

157

-> Vector Double -- ^ starting point

158

-> (Vector Double, Matrix Double) -- ^ solution vector and optimization path

159

160

minimizeD method eps maxit istep tol f df xi = v2l $ minimizeVD

161

method eps maxit istep tol (f.toList) (fromList.df.toList) (fromList xi)

162

where v2l (v,m) = (toList v, m)

163

164

165

minimizeVD method eps maxit istep tol f df xiv = unsafePerformIO $ do

166

let n = dim xiv

167

f' = f

168

df' = (checkdim1 n . df)

169

fp <- mkVecfun (iv f')

170

dfp <- mkVecVecfun (aux_vTov df')

171

rawpath <- vec xiv $ \xiv' ->

172

createMIO maxit (n+2)

173

(c_minimizeD (fi (fromEnum method)) fp dfp istep tol eps (fi maxit) // xiv')

174

"minimizeD"

175

let it = round (rawpath @@> (maxit-1,0))

176

path = takeRows it rawpath

177

sol = cdat $ dropColumns 2 $ dropRows (it-1) path

178

freeHaskellFunPtr fp

179

freeHaskellFunPtr dfp

180

return (sol,path)

181

182

foreign import ccall safe "gsl-aux.h minimizeD"

183

c_minimizeD :: CInt

184

-> FunPtr (CInt -> Ptr Double -> Double)

185

-> FunPtr TVV

186

-> Double -> Double -> Double -> CInt

187

-> TVM

188

189

---------------------------------------------------------------------

190

191

checkdim1 n v

192

| dim v == n = v

193

| otherwise = error $ "Error: "++ show n

194

++ " components expected in the result of the gradient supplied to minimizeD"

Older »