{-# OPTIONS_GHC -Wall #-}
{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE ConstraintKinds #-}
{-# LANGUAGE CPP #-}
{-# LANGUAGE FlexibleContexts #-}
{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE ScopedTypeVariables #-}
{-# LANGUAGE TypeApplications #-}
module Numeric.Optimization
(
minimize
, IsProblem (..)
, HasGrad (..)
, HasHessian (..)
, Constraint (..)
, boundsUnconstrained
, isUnconstainedBounds
, WithGrad (..)
, WithHessian (..)
, WithBounds (..)
, WithConstraints (..)
, Method (..)
, isSupportedMethod
, Params (..)
, Result (..)
, Statistics (..)
, OptimizationException (..)
, Default (..)
, Optionally (..)
, hasOptionalDict
) where
import Control.Exception
import Control.Monad.Primitive
import Control.Monad.ST
import Data.Coerce
import Data.Constraint (Dict (..))
import Data.Default.Class
import Data.Functor.Contravariant
import Data.IORef
import Data.Maybe
import qualified Data.Vector as V
import Data.Vector.Storable (Vector)
import qualified Data.Vector.Generic as VG
import qualified Data.Vector.Generic.Mutable as VGM
import qualified Data.Vector.Storable.Mutable as VSM
import Foreign.C
import qualified Numeric.LBFGS.Vector as LBFGS
#ifdef WITH_CG_DESCENT
import qualified Numeric.Optimization.Algorithms.HagerZhang05 as CG
#endif
import Numeric.LinearAlgebra (Matrix)
import qualified Numeric.LinearAlgebra as LA
data Method
= CGDescent
| LBFGS
| Newton
deriving (Method -> Method -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: Method -> Method -> Bool
$c/= :: Method -> Method -> Bool
== :: Method -> Method -> Bool
$c== :: Method -> Method -> Bool
Eq, Eq Method
Method -> Method -> Bool
Method -> Method -> Ordering
Method -> Method -> Method
forall a.
Eq a
-> (a -> a -> Ordering)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> Bool)
-> (a -> a -> a)
-> (a -> a -> a)
-> Ord a
min :: Method -> Method -> Method
$cmin :: Method -> Method -> Method
max :: Method -> Method -> Method
$cmax :: Method -> Method -> Method
>= :: Method -> Method -> Bool
$c>= :: Method -> Method -> Bool
> :: Method -> Method -> Bool
$c> :: Method -> Method -> Bool
<= :: Method -> Method -> Bool
$c<= :: Method -> Method -> Bool
< :: Method -> Method -> Bool
$c< :: Method -> Method -> Bool
compare :: Method -> Method -> Ordering
$ccompare :: Method -> Method -> Ordering
Ord, Int -> Method
Method -> Int
Method -> [Method]
Method -> Method
Method -> Method -> [Method]
Method -> Method -> Method -> [Method]
forall a.
(a -> a)
-> (a -> a)
-> (Int -> a)
-> (a -> Int)
-> (a -> [a])
-> (a -> a -> [a])
-> (a -> a -> [a])
-> (a -> a -> a -> [a])
-> Enum a
enumFromThenTo :: Method -> Method -> Method -> [Method]
$cenumFromThenTo :: Method -> Method -> Method -> [Method]
enumFromTo :: Method -> Method -> [Method]
$cenumFromTo :: Method -> Method -> [Method]
enumFromThen :: Method -> Method -> [Method]
$cenumFromThen :: Method -> Method -> [Method]
enumFrom :: Method -> [Method]
$cenumFrom :: Method -> [Method]
fromEnum :: Method -> Int
$cfromEnum :: Method -> Int
toEnum :: Int -> Method
$ctoEnum :: Int -> Method
pred :: Method -> Method
$cpred :: Method -> Method
succ :: Method -> Method
$csucc :: Method -> Method
Enum, Int -> Method -> ShowS
[Method] -> ShowS
Method -> String
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [Method] -> ShowS
$cshowList :: [Method] -> ShowS
show :: Method -> String
$cshow :: Method -> String
showsPrec :: Int -> Method -> ShowS
$cshowsPrec :: Int -> Method -> ShowS
Show, Method
forall a. a -> a -> Bounded a
maxBound :: Method
$cmaxBound :: Method
minBound :: Method
$cminBound :: Method
Bounded)
isSupportedMethod :: Method -> Bool
isSupportedMethod :: Method -> Bool
isSupportedMethod Method
LBFGS = Bool
True
#ifdef WITH_CG_DESCENT
isSupportedMethod CGDescent = True
#else
isSupportedMethod Method
CGDescent = Bool
False
#endif
isSupportedMethod Method
Newton = Bool
True
data Params a
= Params
{ forall a. Params a -> Maybe (a -> IO Bool)
paramsCallback :: Maybe (a -> IO Bool)
, forall a. Params a -> Maybe Double
paramsTol :: Maybe Double
}
instance Default (Params a) where
def :: Params a
def =
Params
{ paramsCallback :: Maybe (a -> IO Bool)
paramsCallback = forall a. Maybe a
Nothing
, paramsTol :: Maybe Double
paramsTol = forall a. Maybe a
Nothing
}
instance Contravariant Params where
contramap :: forall a' a. (a' -> a) -> Params a -> Params a'
contramap a' -> a
f Params a
params =
Params a
params
{ paramsCallback :: Maybe (a' -> IO Bool)
paramsCallback = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap ((forall b c a. (b -> c) -> (a -> b) -> a -> c
. a' -> a
f)) (forall a. Params a -> Maybe (a -> IO Bool)
paramsCallback Params a
params)
}
data Result a
= Result
{ forall a. Result a -> Bool
resultSuccess :: Bool
, forall a. Result a -> String
resultMessage :: String
, forall a. Result a -> a
resultSolution :: a
, forall a. Result a -> Double
resultValue :: Double
, forall a. Result a -> Maybe a
resultGrad :: Maybe a
, forall a. Result a -> Maybe (Matrix Double)
resultHessian :: Maybe (Matrix Double)
, forall a. Result a -> Maybe (Matrix Double)
resultHessianInv :: Maybe (Matrix Double)
, forall a. Result a -> Statistics
resultStatistics :: Statistics
}
instance Functor Result where
fmap :: forall a b. (a -> b) -> Result a -> Result b
fmap a -> b
f Result a
result =
Result a
result
{ resultSolution :: b
resultSolution = a -> b
f (forall a. Result a -> a
resultSolution Result a
result)
, resultGrad :: Maybe b
resultGrad = forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap a -> b
f (forall a. Result a -> Maybe a
resultGrad Result a
result)
}
data Statistics
= Statistics
{ Statistics -> Int
totalIters :: Int
, Statistics -> Int
funcEvals :: Int
, Statistics -> Int
gradEvals :: Int
, Statistics -> Int
hessEvals :: Int
}
data OptimizationException
= UnsupportedProblem String
| UnsupportedMethod Method
| GradUnavailable
| HessianUnavailable
deriving (Int -> OptimizationException -> ShowS
[OptimizationException] -> ShowS
OptimizationException -> String
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [OptimizationException] -> ShowS
$cshowList :: [OptimizationException] -> ShowS
show :: OptimizationException -> String
$cshow :: OptimizationException -> String
showsPrec :: Int -> OptimizationException -> ShowS
$cshowsPrec :: Int -> OptimizationException -> ShowS
Show)
instance Exception OptimizationException
class IsProblem prob where
func :: prob -> Vector Double -> Double
bounds :: prob -> Maybe (V.Vector (Double, Double))
bounds prob
_ = forall a. Maybe a
Nothing
constraints :: prob -> [Constraint]
constraints prob
_ = []
{-# MINIMAL func #-}
class IsProblem prob => HasGrad prob where
grad :: prob -> Vector Double -> Vector Double
grad prob
prob = forall a b. (a, b) -> b
snd forall b c a. (b -> c) -> (a -> b) -> a -> c
. forall prob.
HasGrad prob =>
prob -> Vector Double -> (Double, Vector Double)
grad' prob
prob
grad' :: prob -> Vector Double -> (Double, Vector Double)
grad' prob
prob Vector Double
x = forall a. (forall s. ST s a) -> a
runST forall a b. (a -> b) -> a -> b
$ do
MVector s Double
gret <- forall (m :: * -> *) (v :: * -> * -> *) a.
(HasCallStack, PrimMonad m, MVector v a) =>
Int -> m (v (PrimState m) a)
VGM.new (forall (v :: * -> *) a. Vector v a => v a -> Int
VG.length Vector Double
x)
Double
y <- forall prob (m :: * -> *).
(HasGrad prob, PrimMonad m) =>
prob -> Vector Double -> MVector (PrimState m) Double -> m Double
grad'M prob
prob Vector Double
x MVector s Double
gret
Vector Double
g <- forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> m (v a)
VG.unsafeFreeze MVector s Double
gret
forall (m :: * -> *) a. Monad m => a -> m a
return (Double
y, Vector Double
g)
grad'M :: PrimMonad m => prob -> Vector Double -> VSM.MVector (PrimState m) Double -> m Double
grad'M prob
prob Vector Double
x MVector (PrimState m) Double
gvec = do
let y :: Double
y = forall prob. IsProblem prob => prob -> Vector Double -> Double
func prob
prob Vector Double
x
forall (m :: * -> *) (v :: * -> *) a b.
(Monad m, Vector v a) =>
(Int -> a -> m b) -> v a -> m ()
VG.imapM_ (forall (m :: * -> *) (v :: * -> * -> *) a.
(HasCallStack, PrimMonad m, MVector v a) =>
v (PrimState m) a -> Int -> a -> m ()
VGM.write MVector (PrimState m) Double
gvec) (forall prob. HasGrad prob => prob -> Vector Double -> Vector Double
grad prob
prob Vector Double
x)
forall (m :: * -> *) a. Monad m => a -> m a
return Double
y
{-# MINIMAL grad | grad' | grad'M #-}
class IsProblem prob => HasHessian prob where
hessian :: prob -> Vector Double -> Matrix Double
hessianProduct :: prob -> Vector Double -> Vector Double -> Vector Double
hessianProduct prob
prob Vector Double
x Vector Double
v = forall prob.
HasHessian prob =>
prob -> Vector Double -> Matrix Double
hessian prob
prob Vector Double
x forall t. Numeric t => Matrix t -> Vector t -> Vector t
LA.#> Vector Double
v
{-# MINIMAL hessian #-}
class Optionally c where
optionalDict :: Maybe (Dict c)
hasOptionalDict :: c => Maybe (Dict c)
hasOptionalDict :: forall (c :: Constraint). c => Maybe (Dict c)
hasOptionalDict = forall a. a -> Maybe a
Just forall (a :: Constraint). a => Dict a
Dict
data Constraint
boundsUnconstrained :: Int -> V.Vector (Double, Double)
boundsUnconstrained :: Int -> Vector (Double, Double)
boundsUnconstrained Int
n = forall a. Int -> a -> Vector a
V.replicate Int
n (-Double
1forall a. Fractional a => a -> a -> a
/Double
0, Double
1forall a. Fractional a => a -> a -> a
/Double
0)
isUnconstainedBounds :: V.Vector (Double, Double) -> Bool
isUnconstainedBounds :: Vector (Double, Double) -> Bool
isUnconstainedBounds = forall a. (a -> Bool) -> Vector a -> Bool
V.all forall {a} {a}. (RealFloat a, RealFloat a) => (a, a) -> Bool
p
where
p :: (a, a) -> Bool
p (a
lb, a
ub) = forall a. RealFloat a => a -> Bool
isInfinite a
lb Bool -> Bool -> Bool
&& a
lb forall a. Ord a => a -> a -> Bool
< a
0 Bool -> Bool -> Bool
&& forall a. RealFloat a => a -> Bool
isInfinite a
ub Bool -> Bool -> Bool
&& a
ub forall a. Ord a => a -> a -> Bool
> a
0
minimize
:: forall prob. (IsProblem prob, Optionally (HasGrad prob), Optionally (HasHessian prob))
=> Method
-> Params (Vector Double)
-> prob
-> Vector Double
-> IO (Result (Vector Double))
#ifdef WITH_CG_DESCENT
minimize CGDescent =
case optionalDict @(HasGrad prob) of
Just Dict -> minimize_CGDescent
Nothing -> \_ _ _ -> throwIO GradUnavailable
#endif
minimize :: forall prob.
(IsProblem prob, Optionally (HasGrad prob),
Optionally (HasHessian prob)) =>
Method
-> Params (Vector Double)
-> prob
-> Vector Double
-> IO (Result (Vector Double))
minimize Method
LBFGS =
case forall (c :: Constraint). Optionally c => Maybe (Dict c)
optionalDict @(HasGrad prob) of
Just Dict (HasGrad prob)
Dict -> forall prob.
HasGrad prob =>
Params (Vector Double)
-> prob -> Vector Double -> IO (Result (Vector Double))
minimize_LBFGS
Maybe (Dict (HasGrad prob))
Nothing -> \Params (Vector Double)
_ prob
_ Vector Double
_ -> forall e a. Exception e => e -> IO a
throwIO OptimizationException
GradUnavailable
minimize Method
Newton =
case forall (c :: Constraint). Optionally c => Maybe (Dict c)
optionalDict @(HasGrad prob) of
Maybe (Dict (HasGrad prob))
Nothing -> \Params (Vector Double)
_ prob
_ Vector Double
_ -> forall e a. Exception e => e -> IO a
throwIO OptimizationException
GradUnavailable
Just Dict (HasGrad prob)
Dict ->
case forall (c :: Constraint). Optionally c => Maybe (Dict c)
optionalDict @(HasHessian prob) of
Maybe (Dict (HasHessian prob))
Nothing -> \Params (Vector Double)
_ prob
_ Vector Double
_ -> forall e a. Exception e => e -> IO a
throwIO OptimizationException
HessianUnavailable
Just Dict (HasHessian prob)
Dict -> forall prob.
(HasGrad prob, HasHessian prob) =>
Params (Vector Double)
-> prob -> Vector Double -> IO (Result (Vector Double))
minimize_Newton
minimize Method
method = \Params (Vector Double)
_ prob
_ Vector Double
_ -> forall e a. Exception e => e -> IO a
throwIO (Method -> OptimizationException
UnsupportedMethod Method
method)
#ifdef WITH_CG_DESCENT
minimize_CGDescent :: HasGrad prob => Params (Vector Double) -> prob -> Vector Double -> IO (Result (Vector Double))
minimize_CGDescent _params prob _ | not (isNothing (bounds prob)) = throwIO (UnsupportedProblem "CGDescent does not support bounds")
minimize_CGDescent _params prob _ | not (null (constraints prob)) = throwIO (UnsupportedProblem "CGDescent does not support constraints")
minimize_CGDescent params prob x0 = do
let grad_tol = fromMaybe 1e-6 $ paramsTol params
cg_params =
CG.defaultParameters
{ CG.printFinal = False
}
mf :: forall m. PrimMonad m => CG.PointMVector m -> m Double
mf mx = do
x <- VG.unsafeFreeze mx
return $ func prob x
mg :: forall m. PrimMonad m => CG.PointMVector m -> CG.GradientMVector m -> m ()
mg mx mret = do
x <- VG.unsafeFreeze mx
_ <- grad'M prob x mret
return ()
mc :: forall m. PrimMonad m => CG.PointMVector m -> CG.GradientMVector m -> m Double
mc mx mret = do
x <- VG.unsafeFreeze mx
grad'M prob x mret
(x, result, stat) <-
CG.optimize
cg_params
grad_tol
x0
(CG.MFunction mf)
(CG.MGradient mg)
(Just (CG.MCombined mc))
let (success, msg) =
case result of
CG.ToleranceStatisfied -> (True, "convergence tolerance satisfied")
CG.FunctionChange -> (True, "change in func <= feps*|f|")
CG.MaxTotalIter -> (False, "total iterations exceeded maxit")
CG.NegativeSlope -> (False, "slope always negative in line search")
CG.MaxSecantIter -> (False, "number secant iterations exceed nsecant")
CG.NotDescent -> (False, "search direction not a descent direction")
CG.LineSearchFailsInitial -> (False, "line search fails in initial interval")
CG.LineSearchFailsBisection -> (False, "line search fails during bisection")
CG.LineSearchFailsUpdate -> (False, "line search fails during interval update")
CG.DebugTol -> (False, "debugger is on and the function value increases")
CG.FunctionValueNaN -> (False, "function value became nan")
CG.StartFunctionValueNaN -> (False, "starting function value is nan")
return $
Result
{ resultSuccess = success
, resultMessage = msg
, resultSolution = x
, resultValue = CG.finalValue stat
, resultGrad = Nothing
, resultHessian = Nothing
, resultHessianInv = Nothing
, resultStatistics =
Statistics
{ totalIters = fromIntegral $ CG.totalIters stat
, funcEvals = fromIntegral $ CG.funcEvals stat
, gradEvals = fromIntegral $ CG.gradEvals stat
, hessEvals = 0
}
}
#endif
minimize_LBFGS :: HasGrad prob => Params (Vector Double) -> prob -> Vector Double -> IO (Result (Vector Double))
minimize_LBFGS :: forall prob.
HasGrad prob =>
Params (Vector Double)
-> prob -> Vector Double -> IO (Result (Vector Double))
minimize_LBFGS Params (Vector Double)
_params prob
prob Vector Double
_ | Bool -> Bool
not (forall a. Maybe a -> Bool
isNothing (forall prob.
IsProblem prob =>
prob -> Maybe (Vector (Double, Double))
bounds prob
prob)) = forall e a. Exception e => e -> IO a
throwIO (String -> OptimizationException
UnsupportedProblem String
"LBFGS does not support bounds")
minimize_LBFGS Params (Vector Double)
_params prob
prob Vector Double
_ | Bool -> Bool
not (forall (t :: * -> *) a. Foldable t => t a -> Bool
null (forall prob. IsProblem prob => prob -> [Constraint]
constraints prob
prob)) = forall e a. Exception e => e -> IO a
throwIO (String -> OptimizationException
UnsupportedProblem String
"LBFGS does not support constraints")
minimize_LBFGS Params (Vector Double)
params prob
prob Vector Double
x0 = do
IORef Int
evalCounter <- forall a. a -> IO (IORef a)
newIORef (Int
0::Int)
IORef Int
iterRef <- forall a. a -> IO (IORef a)
newIORef (Int
0::Int)
let lbfgsParams :: LBFGSParameters
lbfgsParams =
LBFGS.LBFGSParameters
{ lbfgsPast :: Maybe Int
LBFGS.lbfgsPast = forall a. Maybe a
Nothing
, lbfgsDelta :: Double
LBFGS.lbfgsDelta = forall a. a -> Maybe a -> a
fromMaybe Double
0 forall a b. (a -> b) -> a -> b
$ forall a. Params a -> Maybe Double
paramsTol Params (Vector Double)
params
, lbfgsLineSearch :: LineSearchAlgorithm
LBFGS.lbfgsLineSearch = LineSearchAlgorithm
LBFGS.DefaultLineSearch
, lbfgsL1NormCoefficient :: Maybe Double
LBFGS.lbfgsL1NormCoefficient = forall a. Maybe a
Nothing
}
instanceData :: ()
instanceData :: ()
instanceData = ()
evalFun :: () -> VSM.IOVector CDouble -> VSM.IOVector CDouble -> CInt -> CDouble -> IO CDouble
evalFun :: ()
-> IOVector CDouble
-> IOVector CDouble
-> CInt
-> CDouble
-> IO CDouble
evalFun ()
_inst IOVector CDouble
xvec IOVector CDouble
gvec CInt
_n CDouble
_step = do
forall a. IORef a -> (a -> a) -> IO ()
modifyIORef' IORef Int
evalCounter (forall a. Num a => a -> a -> a
+Int
1)
#if MIN_VERSION_vector(0,13,0)
Vector Double
x <- forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> m (v a)
VG.unsafeFreeze (forall a b s. Coercible a b => MVector s a -> MVector s b
VSM.unsafeCoerceMVector IOVector CDouble
xvec :: VSM.IOVector Double)
Double
y <- forall prob (m :: * -> *).
(HasGrad prob, PrimMonad m) =>
prob -> Vector Double -> MVector (PrimState m) Double -> m Double
grad'M prob
prob Vector Double
x (forall a b s. Coercible a b => MVector s a -> MVector s b
VSM.unsafeCoerceMVector IOVector CDouble
gvec :: VSM.IOVector Double)
#else
x <- VG.unsafeFreeze (coerce xvec :: VSM.IOVector Double)
y <- grad'M prob x (coerce gvec :: VSM.IOVector Double)
#endif
forall (m :: * -> *) a. Monad m => a -> m a
return (coerce :: forall a b. Coercible a b => a -> b
coerce Double
y)
progressFun :: () -> VSM.IOVector CDouble -> VSM.IOVector CDouble -> CDouble -> CDouble -> CDouble -> CDouble -> CInt -> CInt -> CInt -> IO CInt
progressFun :: ()
-> IOVector CDouble
-> IOVector CDouble
-> CDouble
-> CDouble
-> CDouble
-> CDouble
-> CInt
-> CInt
-> CInt
-> IO CInt
progressFun ()
_inst IOVector CDouble
xvec IOVector CDouble
_gvec CDouble
_fx CDouble
_xnorm CDouble
_gnorm CDouble
_step CInt
_n CInt
iter CInt
_nev = do
forall a. IORef a -> a -> IO ()
writeIORef IORef Int
iterRef forall a b. (a -> b) -> a -> b
$! forall a b. (Integral a, Num b) => a -> b
fromIntegral CInt
iter
Bool
shouldStop <-
case forall a. Params a -> Maybe (a -> IO Bool)
paramsCallback Params (Vector Double)
params of
Maybe (Vector Double -> IO Bool)
Nothing -> forall (m :: * -> *) a. Monad m => a -> m a
return Bool
False
Just Vector Double -> IO Bool
callback -> do
#if MIN_VERSION_vector(0,13,0)
Vector Double
x <- forall (m :: * -> *) (v :: * -> *) a.
(PrimMonad m, Vector v a) =>
Mutable v (PrimState m) a -> m (v a)
VG.freeze (forall a b s. Coercible a b => MVector s a -> MVector s b
VSM.unsafeCoerceMVector IOVector CDouble
xvec :: VSM.IOVector Double)
#else
x <- VG.freeze (coerce xvec :: VSM.IOVector Double)
#endif
Vector Double -> IO Bool
callback Vector Double
x
forall (m :: * -> *) a. Monad m => a -> m a
return forall a b. (a -> b) -> a -> b
$ if Bool
shouldStop then CInt
1 else CInt
0
(LBFGSResult
result, [Double]
x_) <- forall a.
LBFGSParameters
-> EvaluateFun a
-> ProgressFun a
-> a
-> [Double]
-> IO (LBFGSResult, [Double])
LBFGS.lbfgs LBFGSParameters
lbfgsParams ()
-> IOVector CDouble
-> IOVector CDouble
-> CInt
-> CDouble
-> IO CDouble
evalFun ()
-> IOVector CDouble
-> IOVector CDouble
-> CDouble
-> CDouble
-> CDouble
-> CDouble
-> CInt
-> CInt
-> CInt
-> IO CInt
progressFun ()
instanceData (forall (v :: * -> *) a. Vector v a => v a -> [a]
VG.toList Vector Double
x0)
let x :: Vector Double
x = forall (v :: * -> *) a. Vector v a => [a] -> v a
VG.fromList [Double]
x_
(Bool
success, String
msg) =
case LBFGSResult
result of
LBFGSResult
LBFGS.Success -> (Bool
True, String
"Success")
LBFGSResult
LBFGS.Stop -> (Bool
True, String
"Stop")
LBFGSResult
LBFGS.AlreadyMinimized -> (Bool
True, String
"The initial variables already minimize the objective function.")
LBFGSResult
LBFGS.UnknownError -> (Bool
False, String
"Unknown error.")
LBFGSResult
LBFGS.LogicError -> (Bool
False, String
"Logic error.")
LBFGSResult
LBFGS.OutOfMemory -> (Bool
False, String
"Insufficient memory.")
LBFGSResult
LBFGS.Canceled -> (Bool
False, String
"The minimization process has been canceled.")
LBFGSResult
LBFGS.InvalidN -> (Bool
False, String
"Invalid number of variables specified.")
LBFGSResult
LBFGS.InvalidNSSE -> (Bool
False, String
"Invalid number of variables (for SSE) specified.")
LBFGSResult
LBFGS.InvalidXSSE -> (Bool
False, String
"The array x must be aligned to 16 (for SSE).")
LBFGSResult
LBFGS.InvalidEpsilon -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::epsilon specified.")
LBFGSResult
LBFGS.InvalidTestPeriod -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::past specified.")
LBFGSResult
LBFGS.InvalidDelta -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::delta specified.")
LBFGSResult
LBFGS.InvalidLineSearch -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::linesearch specified.")
LBFGSResult
LBFGS.InvalidMinStep -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::max_step specified.")
LBFGSResult
LBFGS.InvalidMaxStep -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::max_step specified.")
LBFGSResult
LBFGS.InvalidFtol -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::ftol specified.")
LBFGSResult
LBFGS.InvalidWolfe -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::wolfe specified.")
LBFGSResult
LBFGS.InvalidGtol -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::gtol specified.")
LBFGSResult
LBFGS.InvalidXtol -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::xtol specified.")
LBFGSResult
LBFGS.InvalidMaxLineSearch -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::max_linesearch specified.")
LBFGSResult
LBFGS.InvalidOrthantwise -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::orthantwise_c specified.")
LBFGSResult
LBFGS.InvalidOrthantwiseStart-> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::orthantwise_start specified.")
LBFGSResult
LBFGS.InvalidOrthantwiseEnd -> (Bool
False, String
"Invalid parameter lbfgs_parameter_t::orthantwise_end specified.")
LBFGSResult
LBFGS.OutOfInterval -> (Bool
False, String
"The line-search step went out of the interval of uncertainty.")
LBFGSResult
LBFGS.IncorrectTMinMax -> (Bool
False, String
"A logic error occurred; alternatively, the interval of uncertainty became too small.")
LBFGSResult
LBFGS.RoundingError -> (Bool
False, String
"A rounding error occurred; alternatively, no line-search step satisfies the sufficient decrease and curvature conditions.")
LBFGSResult
LBFGS.MinimumStep -> (Bool
False, String
"The line-search step became smaller than lbfgs_parameter_t::min_step.")
LBFGSResult
LBFGS.MaximumStep -> (Bool
False, String
"The line-search step became larger than lbfgs_parameter_t::max_step.")
LBFGSResult
LBFGS.MaximumLineSearch -> (Bool
False, String
"The line-search routine reaches the maximum number of evaluations.")
LBFGSResult
LBFGS.MaximumIteration -> (Bool
False, String
"The algorithm routine reaches the maximum number of iterations.")
LBFGSResult
LBFGS.WidthTooSmall -> (Bool
False, String
"Relative width of the interval of uncertainty is at most lbfgs_parameter_t::xtol.")
LBFGSResult
LBFGS.InvalidParameters -> (Bool
False, String
"A logic error (negative line-search step) occurred.")
LBFGSResult
LBFGS.IncreaseGradient -> (Bool
False, String
"The current search direction increases the objective function value.")
Int
nEvals <- forall a. IORef a -> IO a
readIORef IORef Int
evalCounter
forall (m :: * -> *) a. Monad m => a -> m a
return forall a b. (a -> b) -> a -> b
$
Result
{ resultSuccess :: Bool
resultSuccess = Bool
success
, resultMessage :: String
resultMessage = String
msg
, resultSolution :: Vector Double
resultSolution = Vector Double
x
, resultValue :: Double
resultValue = forall prob. IsProblem prob => prob -> Vector Double -> Double
func prob
prob Vector Double
x
, resultGrad :: Maybe (Vector Double)
resultGrad = forall a. Maybe a
Nothing
, resultHessian :: Maybe (Matrix Double)
resultHessian = forall a. Maybe a
Nothing
, resultHessianInv :: Maybe (Matrix Double)
resultHessianInv = forall a. Maybe a
Nothing
, resultStatistics :: Statistics
resultStatistics =
Statistics
{ totalIters :: Int
totalIters = forall a. HasCallStack => a
undefined
, funcEvals :: Int
funcEvals = Int
nEvals forall a. Num a => a -> a -> a
+ Int
1
, gradEvals :: Int
gradEvals = Int
nEvals forall a. Num a => a -> a -> a
+ Int
1
, hessEvals :: Int
hessEvals = Int
0
}
}
minimize_Newton :: (HasGrad prob, HasHessian prob) => Params (Vector Double) -> prob -> Vector Double -> IO (Result (Vector Double))
minimize_Newton :: forall prob.
(HasGrad prob, HasHessian prob) =>
Params (Vector Double)
-> prob -> Vector Double -> IO (Result (Vector Double))
minimize_Newton Params (Vector Double)
_params prob
prob Vector Double
_ | Bool -> Bool
not (forall a. Maybe a -> Bool
isNothing (forall prob.
IsProblem prob =>
prob -> Maybe (Vector (Double, Double))
bounds prob
prob)) = forall e a. Exception e => e -> IO a
throwIO (String -> OptimizationException
UnsupportedProblem String
"Newton does not support bounds")
minimize_Newton Params (Vector Double)
_params prob
prob Vector Double
_ | Bool -> Bool
not (forall (t :: * -> *) a. Foldable t => t a -> Bool
null (forall prob. IsProblem prob => prob -> [Constraint]
constraints prob
prob)) = forall e a. Exception e => e -> IO a
throwIO (String -> OptimizationException
UnsupportedProblem String
"Newton does not support constraints")
minimize_Newton Params (Vector Double)
params prob
prob Vector Double
x0 = do
let tol :: Double
tol = forall a. a -> Maybe a -> a
fromMaybe Double
1e-6 (forall a. Params a -> Maybe Double
paramsTol Params (Vector Double)
params)
loop :: Vector Double
-> Double
-> Vector Double
-> Matrix Double
-> Int
-> IO (Result (Vector Double))
loop !Vector Double
x !Double
y !Vector Double
g !Matrix Double
h !Int
n = do
Bool
shouldStop <-
case forall a. Params a -> Maybe (a -> IO Bool)
paramsCallback Params (Vector Double)
params of
Just Vector Double -> IO Bool
callback -> Vector Double -> IO Bool
callback Vector Double
x
Maybe (Vector Double -> IO Bool)
Nothing -> forall (m :: * -> *) a. Monad m => a -> m a
return Bool
False
if Bool
shouldStop then do
forall (m :: * -> *) a. Monad m => a -> m a
return forall a b. (a -> b) -> a -> b
$
Result
{ resultSuccess :: Bool
resultSuccess = Bool
False
, resultMessage :: String
resultMessage = String
"The minimization process has been canceled."
, resultSolution :: Vector Double
resultSolution = Vector Double
x
, resultValue :: Double
resultValue = Double
y
, resultGrad :: Maybe (Vector Double)
resultGrad = forall a. a -> Maybe a
Just Vector Double
g
, resultHessian :: Maybe (Matrix Double)
resultHessian = forall a. a -> Maybe a
Just Matrix Double
h
, resultHessianInv :: Maybe (Matrix Double)
resultHessianInv = forall a. Maybe a
Nothing
, resultStatistics :: Statistics
resultStatistics =
Statistics
{ totalIters :: Int
totalIters = Int
n
, funcEvals :: Int
funcEvals = Int
n
, gradEvals :: Int
gradEvals = Int
n
, hessEvals :: Int
hessEvals = Int
n
}
}
else do
let p :: Vector Double
p = Matrix Double
h forall (c :: * -> *) t.
(LSDiv c, Field t) =>
Matrix t -> c t -> c t
LA.<\> Vector Double
g
x' :: Vector Double
x' = forall (v :: * -> *) a b c.
(Vector v a, Vector v b, Vector v c) =>
(a -> b -> c) -> v a -> v b -> v c
VG.zipWith (-) Vector Double
x Vector Double
p
if forall a. Normed a => a -> Double
LA.norm_Inf (forall (v :: * -> *) a b c.
(Vector v a, Vector v b, Vector v c) =>
(a -> b -> c) -> v a -> v b -> v c
VG.zipWith (-) Vector Double
x' Vector Double
x) forall a. Ord a => a -> a -> Bool
> Double
tol then do
let (Double
y', Vector Double
g') = forall prob.
HasGrad prob =>
prob -> Vector Double -> (Double, Vector Double)
grad' prob
prob Vector Double
x'
h' :: Matrix Double
h' = forall prob.
HasHessian prob =>
prob -> Vector Double -> Matrix Double
hessian prob
prob Vector Double
x'
Vector Double
-> Double
-> Vector Double
-> Matrix Double
-> Int
-> IO (Result (Vector Double))
loop Vector Double
x' Double
y' Vector Double
g' Matrix Double
h' (Int
nforall a. Num a => a -> a -> a
+Int
1)
else do
forall (m :: * -> *) a. Monad m => a -> m a
return forall a b. (a -> b) -> a -> b
$
Result
{ resultSuccess :: Bool
resultSuccess = Bool
True
, resultMessage :: String
resultMessage = String
"success"
, resultSolution :: Vector Double
resultSolution = Vector Double
x
, resultValue :: Double
resultValue = Double
y
, resultGrad :: Maybe (Vector Double)
resultGrad = forall a. a -> Maybe a
Just Vector Double
g
, resultHessian :: Maybe (Matrix Double)
resultHessian = forall a. a -> Maybe a
Just Matrix Double
h
, resultHessianInv :: Maybe (Matrix Double)
resultHessianInv = forall a. Maybe a
Nothing
, resultStatistics :: Statistics
resultStatistics =
Statistics
{ totalIters :: Int
totalIters = Int
n
, funcEvals :: Int
funcEvals = Int
n
, gradEvals :: Int
gradEvals = Int
n
, hessEvals :: Int
hessEvals = Int
n
}
}
let (Double
y0, Vector Double
g0) = forall prob.
HasGrad prob =>
prob -> Vector Double -> (Double, Vector Double)
grad' prob
prob Vector Double
x0
h0 :: Matrix Double
h0 = forall prob.
HasHessian prob =>
prob -> Vector Double -> Matrix Double
hessian prob
prob Vector Double
x0
Vector Double
-> Double
-> Vector Double
-> Matrix Double
-> Int
-> IO (Result (Vector Double))
loop Vector Double
x0 Double
y0 Vector Double
g0 Matrix Double
h0 Int
1
instance IsProblem (Vector Double -> Double) where
func :: (Vector Double -> Double) -> Vector Double -> Double
func Vector Double -> Double
f = Vector Double -> Double
f
instance Optionally (HasGrad (Vector Double -> Double)) where
optionalDict :: Maybe (Dict (HasGrad (Vector Double -> Double)))
optionalDict = forall a. Maybe a
Nothing
instance Optionally (HasHessian (Vector Double -> Double)) where
optionalDict :: Maybe (Dict (HasHessian (Vector Double -> Double)))
optionalDict = forall a. Maybe a
Nothing
data WithGrad prob = WithGrad prob (Vector Double -> Vector Double)
instance IsProblem prob => IsProblem (WithGrad prob) where
func :: WithGrad prob -> Vector Double -> Double
func (WithGrad prob
prob Vector Double -> Vector Double
_g) = forall prob. IsProblem prob => prob -> Vector Double -> Double
func prob
prob
bounds :: WithGrad prob -> Maybe (Vector (Double, Double))
bounds (WithGrad prob
prob Vector Double -> Vector Double
_g) = forall prob.
IsProblem prob =>
prob -> Maybe (Vector (Double, Double))
bounds prob
prob
constraints :: WithGrad prob -> [Constraint]
constraints (WithGrad prob
prob Vector Double -> Vector Double
_g) = forall prob. IsProblem prob => prob -> [Constraint]
constraints prob
prob
instance IsProblem prob => HasGrad (WithGrad prob) where
grad :: WithGrad prob -> Vector Double -> Vector Double
grad (WithGrad prob
_prob Vector Double -> Vector Double
g) = Vector Double -> Vector Double
g
instance HasHessian prob => HasHessian (WithGrad prob) where
hessian :: WithGrad prob -> Vector Double -> Matrix Double
hessian (WithGrad prob
prob Vector Double -> Vector Double
_g) = forall prob.
HasHessian prob =>
prob -> Vector Double -> Matrix Double
hessian prob
prob
hessianProduct :: WithGrad prob -> Vector Double -> Vector Double -> Vector Double
hessianProduct (WithGrad prob
prob Vector Double -> Vector Double
_g) = forall prob.
HasHessian prob =>
prob -> Vector Double -> Vector Double -> Vector Double
hessianProduct prob
prob
instance IsProblem prob => Optionally (HasGrad (WithGrad prob)) where
optionalDict :: Maybe (Dict (HasGrad (WithGrad prob)))
optionalDict = forall (c :: Constraint). c => Maybe (Dict c)
hasOptionalDict
instance Optionally (HasHessian prob) => Optionally (HasHessian (WithGrad prob)) where
optionalDict :: Maybe (Dict (HasHessian (WithGrad prob)))
optionalDict =
case forall (c :: Constraint). Optionally c => Maybe (Dict c)
optionalDict @(HasHessian prob) of
Just Dict (HasHessian prob)
Dict -> forall (c :: Constraint). c => Maybe (Dict c)
hasOptionalDict
Maybe (Dict (HasHessian prob))
Nothing -> forall a. Maybe a
Nothing
data WithHessian prob = WithHessian prob (Vector Double -> Matrix Double)
instance IsProblem prob => IsProblem (WithHessian prob) where
func :: WithHessian prob -> Vector Double -> Double
func (WithHessian prob
prob Vector Double -> Matrix Double
_hess) = forall prob. IsProblem prob => prob -> Vector Double -> Double
func prob
prob
bounds :: WithHessian prob -> Maybe (Vector (Double, Double))
bounds (WithHessian prob
prob Vector Double -> Matrix Double
_hess) = forall prob.
IsProblem prob =>
prob -> Maybe (Vector (Double, Double))
bounds prob
prob
constraints :: WithHessian prob -> [Constraint]
constraints (WithHessian prob
prob Vector Double -> Matrix Double
_hess) = forall prob. IsProblem prob => prob -> [Constraint]
constraints prob
prob
instance HasGrad prob => HasGrad (WithHessian prob) where
grad :: WithHessian prob -> Vector Double -> Vector Double
grad (WithHessian prob
prob Vector Double -> Matrix Double
_) = forall prob. HasGrad prob => prob -> Vector Double -> Vector Double
grad prob
prob
instance IsProblem prob => HasHessian (WithHessian prob) where
hessian :: WithHessian prob -> Vector Double -> Matrix Double
hessian (WithHessian prob
_prob Vector Double -> Matrix Double
hess) = Vector Double -> Matrix Double
hess
instance Optionally (HasGrad prob) => Optionally (HasGrad (WithHessian prob)) where
optionalDict :: Maybe (Dict (HasGrad (WithHessian prob)))
optionalDict =
case forall (c :: Constraint). Optionally c => Maybe (Dict c)
optionalDict @(HasGrad prob) of
Just Dict (HasGrad prob)
Dict -> forall (c :: Constraint). c => Maybe (Dict c)
hasOptionalDict
Maybe (Dict (HasGrad prob))
Nothing -> forall a. Maybe a
Nothing
instance IsProblem prob => Optionally (HasHessian (WithHessian prob)) where
optionalDict :: Maybe (Dict (HasHessian (WithHessian prob)))
optionalDict = forall (c :: Constraint). c => Maybe (Dict c)
hasOptionalDict
data WithBounds prob = WithBounds prob (V.Vector (Double, Double))
instance IsProblem prob => IsProblem (WithBounds prob) where
func :: WithBounds prob -> Vector Double -> Double
func (WithBounds prob
prob Vector (Double, Double)
_bounds) = forall prob. IsProblem prob => prob -> Vector Double -> Double
func prob
prob
bounds :: WithBounds prob -> Maybe (Vector (Double, Double))
bounds (WithBounds prob
_prob Vector (Double, Double)
bounds) = forall a. a -> Maybe a
Just Vector (Double, Double)
bounds
constraints :: WithBounds prob -> [Constraint]
constraints (WithBounds prob
prob Vector (Double, Double)
_bounds) = forall prob. IsProblem prob => prob -> [Constraint]
constraints prob
prob
instance HasGrad prob => HasGrad (WithBounds prob) where
grad :: WithBounds prob -> Vector Double -> Vector Double
grad (WithBounds prob
prob Vector (Double, Double)
_bounds) = forall prob. HasGrad prob => prob -> Vector Double -> Vector Double
grad prob
prob
grad' :: WithBounds prob -> Vector Double -> (Double, Vector Double)
grad' (WithBounds prob
prob Vector (Double, Double)
_bounds) = forall prob.
HasGrad prob =>
prob -> Vector Double -> (Double, Vector Double)
grad' prob
prob
grad'M :: forall (m :: * -> *).
PrimMonad m =>
WithBounds prob
-> Vector Double -> MVector (PrimState m) Double -> m Double
grad'M (WithBounds prob
prob Vector (Double, Double)
_bounds) = forall prob (m :: * -> *).
(HasGrad prob, PrimMonad m) =>
prob -> Vector Double -> MVector (PrimState m) Double -> m Double
grad'M prob
prob
instance HasHessian prob => HasHessian (WithBounds prob) where
hessian :: WithBounds prob -> Vector Double -> Matrix Double
hessian (WithBounds prob
prob Vector (Double, Double)
_bounds) = forall prob.
HasHessian prob =>
prob -> Vector Double -> Matrix Double
hessian prob
prob
hessianProduct :: WithBounds prob -> Vector Double -> Vector Double -> Vector Double
hessianProduct (WithBounds prob
prob Vector (Double, Double)
_bounds) = forall prob.
HasHessian prob =>
prob -> Vector Double -> Vector Double -> Vector Double
hessianProduct prob
prob
instance Optionally (HasGrad prob) => Optionally (HasGrad (WithBounds prob)) where
optionalDict :: Maybe (Dict (HasGrad (WithBounds prob)))
optionalDict =
case forall (c :: Constraint). Optionally c => Maybe (Dict c)
optionalDict @(HasGrad prob) of
Just Dict (HasGrad prob)
Dict -> forall (c :: Constraint). c => Maybe (Dict c)
hasOptionalDict
Maybe (Dict (HasGrad prob))
Nothing -> forall a. Maybe a
Nothing
instance Optionally (HasHessian prob) => Optionally (HasHessian (WithBounds prob)) where
optionalDict :: Maybe (Dict (HasHessian (WithBounds prob)))
optionalDict =
case forall (c :: Constraint). Optionally c => Maybe (Dict c)
optionalDict @(HasHessian prob) of
Just Dict (HasHessian prob)
Dict -> forall (c :: Constraint). c => Maybe (Dict c)
hasOptionalDict
Maybe (Dict (HasHessian prob))
Nothing -> forall a. Maybe a
Nothing
data WithConstraints prob = WithConstraints prob [Constraint]
instance IsProblem prob => IsProblem (WithConstraints prob) where
func :: WithConstraints prob -> Vector Double -> Double
func (WithConstraints prob
prob [Constraint]
_constraints) = forall prob. IsProblem prob => prob -> Vector Double -> Double
func prob
prob
bounds :: WithConstraints prob -> Maybe (Vector (Double, Double))
bounds (WithConstraints prob
prob [Constraint]
_constraints) = forall prob.
IsProblem prob =>
prob -> Maybe (Vector (Double, Double))
bounds prob
prob
constraints :: WithConstraints prob -> [Constraint]
constraints (WithConstraints prob
_prob [Constraint]
constraints) = [Constraint]
constraints
instance HasGrad prob => HasGrad (WithConstraints prob) where
grad :: WithConstraints prob -> Vector Double -> Vector Double
grad (WithConstraints prob
prob [Constraint]
_constraints) = forall prob. HasGrad prob => prob -> Vector Double -> Vector Double
grad prob
prob
grad' :: WithConstraints prob -> Vector Double -> (Double, Vector Double)
grad' (WithConstraints prob
prob [Constraint]
_constraints) = forall prob.
HasGrad prob =>
prob -> Vector Double -> (Double, Vector Double)
grad' prob
prob
grad'M :: forall (m :: * -> *).
PrimMonad m =>
WithConstraints prob
-> Vector Double -> MVector (PrimState m) Double -> m Double
grad'M (WithConstraints prob
prob [Constraint]
_constraints) = forall prob (m :: * -> *).
(HasGrad prob, PrimMonad m) =>
prob -> Vector Double -> MVector (PrimState m) Double -> m Double
grad'M prob
prob
instance HasHessian prob => HasHessian (WithConstraints prob) where
hessian :: WithConstraints prob -> Vector Double -> Matrix Double
hessian (WithConstraints prob
prob [Constraint]
_constraints) = forall prob.
HasHessian prob =>
prob -> Vector Double -> Matrix Double
hessian prob
prob
hessianProduct :: WithConstraints prob
-> Vector Double -> Vector Double -> Vector Double
hessianProduct (WithConstraints prob
prob [Constraint]
_constraints) = forall prob.
HasHessian prob =>
prob -> Vector Double -> Vector Double -> Vector Double
hessianProduct prob
prob
instance Optionally (HasGrad prob) => Optionally (HasGrad (WithConstraints prob)) where
optionalDict :: Maybe (Dict (HasGrad (WithConstraints prob)))
optionalDict =
case forall (c :: Constraint). Optionally c => Maybe (Dict c)
optionalDict @(HasGrad prob) of
Just Dict (HasGrad prob)
Dict -> forall (c :: Constraint). c => Maybe (Dict c)
hasOptionalDict
Maybe (Dict (HasGrad prob))
Nothing -> forall a. Maybe a
Nothing
instance Optionally (HasHessian prob) => Optionally (HasHessian (WithConstraints prob)) where
optionalDict :: Maybe (Dict (HasHessian (WithConstraints prob)))
optionalDict =
case forall (c :: Constraint). Optionally c => Maybe (Dict c)
optionalDict @(HasHessian prob) of
Just Dict (HasHessian prob)
Dict -> forall (c :: Constraint). c => Maybe (Dict c)
hasOptionalDict
Maybe (Dict (HasHessian prob))
Nothing -> forall a. Maybe a
Nothing