module Algorithm.GoldenSectionSearch where

import Data.Coerce
import Data.Ord
import Data.Semigroup

-- | phi = (1+sqrt 5) / 2
phi :: Double
phi :: Double
phi = Double
1.618033988749895
{-# INLINE phi #-}

-- | resphi = 2-phi = 1/(1+phi)
resphi :: Double
resphi :: Double
resphi = Double
0.3819660112501051
{-# INLINE resphi #-}

-- | mid1 (mid1 low high) high = mid2 low high
mid1 :: Double -> Double -> Double
mid1 :: Double -> Double -> Double
mid1 Double
low Double
high = (Double
low Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
phi Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
high) Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
resphi
{-# INLINE mid1 #-}

-- | mid2 low (mid2 low high) = mid1 low high
mid2 :: Double -> Double -> Double
mid2 :: Double -> Double -> Double
mid2 Double
low Double
high = (Double
low Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
high Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
phi) Double -> Double -> Double
forall a. Num a => a -> a -> a
* Double
resphi
{-# INLINE mid2 #-}

epsGS :: Double
epsGS :: Double
epsGS = Double
1e-12
{-# INLINE epsGS #-}

{- |
@f@ should be unimodal

>>> goldenSectionSearchMin 0.0 10.0 (\x -> (x-1)^2)
Min {getMin = Arg 9.802089785436835e-28 0.9999999999999687}
-}
goldenSectionSearchMin ::
  (Ord a) => Double -> Double -> (Double -> a) -> ArgMin a Double
goldenSectionSearchMin :: forall a.
Ord a =>
Double -> Double -> (Double -> a) -> ArgMin a Double
goldenSectionSearchMin Double
low Double
high Double -> a
f =
  Int
-> Double
-> Double
-> Double
-> Double
-> a
-> a
-> Min (Arg a Double)
forall {t}.
(Num t, Eq t) =>
t
-> Double
-> Double
-> Double
-> Double
-> a
-> a
-> Min (Arg a Double)
go (Int
72 :: Int) Double
low Double
x01 Double
x02 Double
high (Double -> a
f Double
x01) (Double -> a
f Double
x02)
  where
    !x01 :: Double
x01 = Double -> Double -> Double
mid1 Double
low Double
high
    !x02 :: Double
x02 = Double -> Double -> Double
mid2 Double
low Double
high
    go :: t
-> Double
-> Double
-> Double
-> Double
-> a
-> a
-> Min (Arg a Double)
go !t
n !Double
x0 !Double
x1 !Double
x2 !Double
x3 !a
fx1 !a
fx2
      | t
n t -> t -> Bool
forall a. Eq a => a -> a -> Bool
== t
0 Bool -> Bool -> Bool
|| Double -> Double
forall a. Num a => a -> a
abs (Double
x3 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
x0) Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
< Double
epsGS = Arg a Double -> Min (Arg a Double)
forall a. a -> Min a
Min (a -> Double -> Arg a Double
forall a b. a -> b -> Arg a b
Arg a
fx1 Double
x1)
      | a
fx1 a -> a -> Bool
forall a. Ord a => a -> a -> Bool
< a
fx2 =
          let !x :: Double
x = Double -> Double -> Double
mid1 Double
x0 Double
x2 -- mid2 x0 x2 == x1
           in t
-> Double
-> Double
-> Double
-> Double
-> a
-> a
-> Min (Arg a Double)
go (t
n t -> t -> t
forall a. Num a => a -> a -> a
- t
1) Double
x0 Double
x Double
x1 Double
x2 (Double -> a
f Double
x) a
fx1
      | Bool
otherwise =
          let !x :: Double
x = Double -> Double -> Double
mid2 Double
x1 Double
x3 -- mid1 x1 x3 == x2
           in t
-> Double
-> Double
-> Double
-> Double
-> a
-> a
-> Min (Arg a Double)
go (t
n t -> t -> t
forall a. Num a => a -> a -> a
- t
1) Double
x1 Double
x2 Double
x Double
x3 a
fx2 (Double -> a
f Double
x)

{- |
>>> goldenSectionSearchMax 0.0 10.0 (\x -> log x/x)
Max {getMax = Arg 0.36787944117144233 2.7182818604248506}
-}
goldenSectionSearchMax ::
  (Ord a) => Double -> Double -> (Double -> a) -> ArgMax a Double
goldenSectionSearchMax :: forall a.
Ord a =>
Double -> Double -> (Double -> a) -> ArgMax a Double
goldenSectionSearchMax Double
low Double
high Double -> a
f =
  ArgMin (Down a) Double -> ArgMax a Double
forall a b. Coercible a b => a -> b
coerce (ArgMin (Down a) Double -> ArgMax a Double)
-> ArgMin (Down a) Double -> ArgMax a Double
forall a b. (a -> b) -> a -> b
$
    Double -> Double -> (Double -> Down a) -> ArgMin (Down a) Double
forall a.
Ord a =>
Double -> Double -> (Double -> a) -> ArgMin a Double
goldenSectionSearchMin Double
low Double
high (a -> Down a
forall a. a -> Down a
Down (a -> Down a) -> (Double -> a) -> Double -> Down a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Double -> a
f)

goldenSectionSearchMin2 ::
  (Ord a) =>
  Double ->
  Double ->
  (Double -> Double -> a) ->
  ArgMin a (Double, Double)
goldenSectionSearchMin2 :: forall a.
Ord a =>
Double
-> Double -> (Double -> Double -> a) -> ArgMin a (Double, Double)
goldenSectionSearchMin2 Double
low Double
high Double -> Double -> a
f = Arg a (Double, Double) -> Min (Arg a (Double, Double))
forall a. a -> Min a
Min (a -> (Double, Double) -> Arg a (Double, Double)
forall a b. a -> b -> Arg a b
Arg a
fxy (Double
x, Double
y))
  where
    Min (Arg (Min (Arg a
fxy Double
y)) Double
x) =
      Double
-> Double
-> (Double -> ArgMin a Double)
-> Min (Arg (ArgMin a Double) Double)
forall a.
Ord a =>
Double -> Double -> (Double -> a) -> ArgMin a Double
goldenSectionSearchMin Double
low Double
high ((Double -> ArgMin a Double) -> Min (Arg (ArgMin a Double) Double))
-> (Double -> ArgMin a Double)
-> Min (Arg (ArgMin a Double) Double)
forall a b. (a -> b) -> a -> b
$ \Double
cx ->
        Double -> Double -> (Double -> a) -> ArgMin a Double
forall a.
Ord a =>
Double -> Double -> (Double -> a) -> ArgMin a Double
goldenSectionSearchMin Double
low Double
high ((Double -> a) -> ArgMin a Double)
-> (Double -> a) -> ArgMin a Double
forall a b. (a -> b) -> a -> b
$ \Double
cy ->
          Double -> Double -> a
f Double
cx Double
cy

goldenSectionSearchMax2 ::
  (Ord a) =>
  Double ->
  Double ->
  (Double -> Double -> a) ->
  ArgMax a (Double, Double)
goldenSectionSearchMax2 :: forall a.
Ord a =>
Double
-> Double -> (Double -> Double -> a) -> ArgMax a (Double, Double)
goldenSectionSearchMax2 Double
low Double
high Double -> Double -> a
f =
  ArgMin (Down a) (Double, Double) -> ArgMax a (Double, Double)
forall a b. Coercible a b => a -> b
coerce (ArgMin (Down a) (Double, Double) -> ArgMax a (Double, Double))
-> ArgMin (Down a) (Double, Double) -> ArgMax a (Double, Double)
forall a b. (a -> b) -> a -> b
$
    Double
-> Double
-> (Double -> Double -> Down a)
-> ArgMin (Down a) (Double, Double)
forall a.
Ord a =>
Double
-> Double -> (Double -> Double -> a) -> ArgMin a (Double, Double)
goldenSectionSearchMin2 Double
low Double
high (\Double
x Double
y -> a -> Down a
forall a. a -> Down a
Down (Double -> Double -> a
f Double
x Double
y))