module Ipopt.Sparsity (
nanPropagate1,
nanPropagateG,
nanPropagateJ,
nanPropagateH,
nanPropagateG1,
nanPropagateHF,
jacobianSS,
) where
import Control.Applicative
import Control.Monad
import qualified Data.Vector as V
import Data.Vector ((!))
import Numeric.AD
import qualified Data.IntMap as M
import qualified Data.IntSet as S
import Data.Monoid
import Ipopt.AnyRF
nanPropagateG1 :: Int
-> (forall a. AnyRFCxt a => V.Vector a -> a)
-> [(Int, V.Vector Int)]
nanPropagateG1 nx f | x0:_ <- dropWhile (isNaN . f) (trialV nx [0,0.5,1])
= [ (i,j) | i <- [0 .. nx1],
let g = grad f (x0 V.// [(i,(0/0))])
j = V.findIndices isNaN g,
not (V.null j) ]
nanPropagateG :: Int
-> (forall a. AnyRFCxt a => V.Vector a -> a)
-> V.Vector Int
nanPropagateG nx f = V.findIndices (\x -> x /= 0)
(grad f (V.replicate nx (0/0)))
nanPropagate1 :: Int
-> (forall a. AnyRFCxt a => V.Vector a -> a)
-> [Int]
nanPropagate1 nx f | x0:_ <- dropWhile (isNaN . f) (trialV nx [0,0.5,1])
= [ i | i <- [0 .. nx1],
let g = f (x0 V.// [(i,(0/0))]),
isNaN g ]
nanPropagateH :: Int
-> (forall a. AnyRFCxt a => V.Vector a -> a)
-> V.Vector (Int, Int)
nanPropagateH nx f = nonzeroIxs $ hessian f (V.replicate nx (0/0::Double))
nanPropagateJ :: Int
-> (forall a. AnyRFCxt a => V.Vector a -> V.Vector a)
-> V.Vector (Int,Int)
nanPropagateJ nx f = nonzeroIxs $ jacobian f (V.replicate nx (0/0::Double))
nanPropagateHF :: Int
-> (forall a. AnyRFCxt a => V.Vector a -> V.Vector a)
-> V.Vector (V.Vector (Int,Int))
nanPropagateHF nx f =
V.map nonzeroIxs $ hessianF f (V.replicate nx (0/0::Double))
collapse :: V.Vector (V.Vector (Int,Int)) -> V.Vector (Int,Int)
collapse = V.fromList . concatMap (\(a,bs) -> map (a,) (S.toList bs)) . M.toList .
M.mapWithKey (\k s -> fst (S.split (k+1) s)) .
M.fromListWith (<>) . map (\(a,b) -> (a, S.singleton b)) . concatMap V.toList . V.toList
nonzeroIxs = V.concatMap id . V.imap (\i -> V.map (i,) . V.findIndices (/= 0))
trialV :: Int -> [Double] -> [V.Vector Double]
trialV nx cand = V.fromList <$> replicateM nx cand
jacobianSS :: AnyRFCxt a => V.Vector (Int,Int)
-> (forall a. AnyRFCxt a => V.Vector a -> V.Vector a)
-> V.Vector a
-> V.Vector a
jacobianSS ijs f x = let v = jacobian f x
in V.map (\(i,j) -> v V.! i V.! j) ijs