diff --git a/.gitignore b/.gitignore new file mode 100644 --- /dev/null +++ b/.gitignore @@ -0,0 +1,2 @@ +*~ +dist diff --git a/LICENSE b/LICENSE new file mode 100644 --- /dev/null +++ b/LICENSE @@ -0,0 +1,19 @@ +Copyright (c) 2010 Andrew Miller + +Permission is hereby granted, free of charge, to any person obtaining a copy +of this software and associated documentation files (the "Software"), to deal +in the Software without restriction, including without limitation the rights +to use, copy, modify, merge, publish, distribute, sublicense, and/or sell +copies of the Software, and to permit persons to whom the Software is +furnished to do so, subject to the following conditions: + +The above copyright notice and this permission notice shall be included in +all copies or substantial portions of the Software. + +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, +OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN +THE SOFTWARE. diff --git a/Setup.hs b/Setup.hs new file mode 100644 --- /dev/null +++ b/Setup.hs @@ -0,0 +1,3 @@ +#!/usr/bin/env runhaskell +import Distribution.Simple +main = defaultMain diff --git a/bindings-eskit.cabal b/bindings-eskit.cabal new file mode 100644 --- /dev/null +++ b/bindings-eskit.cabal @@ -0,0 +1,37 @@ +name: bindings-eskit +version: 0.0.1 +synopsis: Bindings to ESKit. +description: + Bindings to ESKit . +license: BSD3 +license-file: LICENSE +author: Andrew Miller +maintainer: Andrew Miller +homepage: http://github.com/a1kmm/bindings-eskit +bug-reports: http://github.com/a1kmm/bindings-eskit/issues +category: Numerical, FFI + +build-type: Simple +cabal-version: >= 1.2.3 + +library + hs-source-dirs: + src + exposed-modules: + Bindings.ESKitLowlevel, + Bindings.ESKit + build-depends: + base >= 3 && < 5, + bindings-DSL >= 1.0 && < 1.1, + array >= 0.3 && < 0.4 + ghc-options: + -Wall + extensions: + ForeignFunctionInterface + includes: + eskit.h + extra-libraries: + eskit + pkgconfig-depends: + eskit >= 1.1.1 + c-sources: src/Bindings/eskit_support.c diff --git a/src/Bindings/ESKit.hs b/src/Bindings/ESKit.hs new file mode 100644 --- /dev/null +++ b/src/Bindings/ESKit.hs @@ -0,0 +1,117 @@ +module Bindings.ESKit (OptimisationProblem(OptimisationProblem), + DistributionHandler(CMAHandler, SepCMAHandler, CSAHandler), + n, mu, lambda, fitnessFunction, distributionHandlerSettings, + optimize) +where + +import qualified Bindings.ESKitLowlevel as LL +import Control.Monad +import Data.Maybe +import Foreign +import Foreign.C.Types + +data OptimisationProblem = OptimisationProblem { + n :: Int, + mu :: Maybe Int, + lambda :: Maybe Int, + fitnessFunction :: [Double] -> Double, + distributionHandlerSettings :: DistributionHandler + }; + +data DistributionHandler = + CMAHandler (Double, Double) | + SepCMAHandler (Double, Double) | + CSAHandler (Double, Double) + +setupDistribHandler :: CSize -> Ptr LL.C'ekOptimizer -> DistributionHandler -> IO (Ptr CChar) +setupDistribHandler np eo (CMAHandler (sigmalow, sigmahigh)) = do + p <- LL.c'mkekCMA np + LL.c'ekCMA_setSigma p (realToFrac sigmahigh) (realToFrac sigmalow) + LL.c'ekCMA_setOptimizer p eo + return $ castPtr p + +setupDistribHandler np eo (SepCMAHandler (sigmalow, sigmahigh)) = do + p <- LL.c'mkekSepCMA np + LL.c'ekSepCMA_setSigma p (realToFrac sigmahigh) (realToFrac sigmalow) + LL.c'ekSepCMA_setOptimizer p eo + return $ castPtr p + +setupDistribHandler np eo (CSAHandler (sigmalow, sigmahigh)) = do + p <- LL.c'mkekCSA np + LL.c'ekCSA_setSigma p (realToFrac sigmahigh) (realToFrac sigmalow) + LL.c'ekCSA_setOptimizer p eo + return $ castPtr p + +destroyDistribHandler :: Ptr CChar -> DistributionHandler -> IO () +destroyDistribHandler dh (CMAHandler _) = + LL.c'delekCMA (castPtr dh) +destroyDistribHandler dh (SepCMAHandler _) = + LL.c'delekSepCMA (castPtr dh) +destroyDistribHandler dh (CSAHandler _) = + LL.c'delekCSA (castPtr dh) + +optimize :: OptimisationProblem -> Either String [Double] +optimize (OptimisationProblem { n = np, mu = mmu, lambda = mlambda, + fitnessFunction = fitness, + distributionHandlerSettings = dhsettings }) = + unsafePerformIO $ do + eo <- LL.c'mkekOptimizer (fromIntegral np) + dh <- setupDistribHandler (fromIntegral np) eo dhsettings + eod <- peek eo + let lambdav = (flip fromMaybe (fmap fromIntegral mlambda) (LL.c'ekOptimizer'lambda eod)) + (LL.c'ekOptimizer_setMuLambda eo) + (flip fromMaybe (fmap fromIntegral mmu) (LL.c'ekOptimizer'mu eod)) + lambdav + LL.c'ekOptimizer_start eo + let continueOptimization = + do + LL.c'ekOptimizer_sampleCloud eo + allPts <- liftM LL.c'ekOptimizer'points (peek eo) + forM_ [0..((fromIntegral lambdav)-1)] $ \i -> + do + ptp <- peekElemOff allPts i + pt <- peek ptp + let xcoordp = LL.c'ekPoint'x pt + xcoord <- mapM (liftM (realToFrac :: CDouble -> Double) . peekElemOff xcoordp) [0..(np-1)] + poke ptp (pt {LL.c'ekPoint'fitness = realToFrac $ fitness xcoord }) + LL.c'ekOptimizer_update eo + stop <- LL.c'ekOptimizer_stop eo + if stop == 0 + then + continueOptimization + else + return stop + whyStop <- continueOptimization + destroyDistribHandler dh dhsettings + ret <- + do + bestPoint <- liftM LL.c'ekOptimizer'bestPoint (peek eo) + let xcoordp = LL.c'ekPoint'x bestPoint + xcoord <- mapM (liftM (realToFrac :: CDouble -> Double) . peekElemOff xcoordp) [0..(np-1)] + return $ Right xcoord + {- + case () + of + () | whyStop == LL.c'ekStopCriterionId_LowSigma -> + () | whyStop == LL.c'ekStopCriterionId_DistributionNotSet -> + return $ Left "No point distribution handler has been associated with the opimiser" + () | whyStop == LL.c'ekStopCriterionId_NoEffectAxis -> + return $ Left "Axes of the Gaussian distribution are beyond what numerical precision can handle" + () | whyStop == LL.c'ekStopCriterionId_NoEffectCoord -> + return $ Left "Eigen vector basis of the Gaussian distribution are beyond what numerical precision can handle" + () | whyStop == LL.c'ekStopCriterionId_ConditionCov -> + return $ Left "Covariance matrix conditionning is beyond what numerical precision can handle." + () -> + return $ Left "The eigenvalue solver failed." + -} + LL.c'delekOptimizer eo + return ret + +{- +To do: tests: + optimize (OptimisationProblem { n = 2, mu = Nothing, lambda = Nothing, fitnessFunction = \(x:y:[]) -> (1-x)^2 + 100 * (y - x^2)^2, distributionHandlerSettings = CMAHandler (0.0001,1000)}) + => Should be close to [1,1] + + optimize (OptimisationProblem { n = 1, mu = Nothing, lambda = Nothing, fitnessFunction = \(x:[]) -> abs((x-100) * (x - 200) * (x - 300) * (x + 400)), distributionHandlerSettings = CMAHandler (0.0001,1000)}) + => Several solutions (obviously). +-} diff --git a/src/Bindings/ESKitLowlevel.hsc b/src/Bindings/ESKitLowlevel.hsc new file mode 100644 --- /dev/null +++ b/src/Bindings/ESKitLowlevel.hsc @@ -0,0 +1,63 @@ +#include +#include + +module Bindings.ESKitLowlevel where +#strict_import + +#opaque_t ekCMA +#opaque_t ekSepCMA +#opaque_t ekCSA + +#starttype ekPoint +#field x , Ptr CDouble +#field fitness, CDouble +#stoptype + +#starttype ekOptimizer +#field N, CSize +#field mu, CSize +#field lambda, CSize +#field nbUpdates, CSize +#field xMean, Ptr CDouble +#field points, Ptr (Ptr ) +#field bestPoint, +#stoptype + +#define hsc_sizeof(type) \ + bc_varid(#type) printf("Size = %u\n", sizeof(type)); \ + bc_varid(#type) printf("Size :: (Num t) => t\n"); + +#sizeof ekPoint +#ccall mkekOptimizer, CSize -> IO (Ptr ) +#ccall delekOptimizer, Ptr -> IO () +#ccall mkekCMA, CSize -> IO (Ptr ) +#ccall delekCMA, Ptr -> IO () +#ccall mkekSepCMA, CSize -> IO (Ptr ) +#ccall delekSepCMA, Ptr -> IO () +#ccall mkekCSA, CSize -> IO (Ptr ) +#ccall delekCSA, Ptr -> IO () + +#ccall ekOptimizer_setMuLambda, Ptr -> CSize -> CSize -> IO () +#ccall ekOptimizer_start, Ptr -> IO () +#ccall ekOptimizer_sampleCloud, Ptr -> IO () +#ccall ekOptimizer_samplePoint, Ptr -> CSize -> IO () +#ccall ekOptimizer_update, Ptr -> IO () +#integral_t enum ekStopCriterionId +#num ekStopCriterionId_DistributionNotSet +#num ekStopCriterionId_LowSigma +#num ekStopCriterionId_NoEffectAxis +#num ekStopCriterionId_NoEffectCoord +#num ekStopCriterionId_ConditionCov +#num ekStopCriterionId_EigenSolverFailure +#ccall ekOptimizer_stop, Ptr -> IO +#callback ekEvalFunction, Ptr CDouble -> CSize -> IO CDouble +#ccall ekOptimizer_evaluateFunction, Ptr -> -> IO () +#ccall ekCMA_setSigma, Ptr -> CDouble -> CDouble -> IO () +#ccall ekCMA_setOptimizer, Ptr -> Ptr -> IO () +#ccall ekCMA_sigma, Ptr -> CDouble +#ccall ekSepCMA_setSigma, Ptr -> CDouble -> CDouble -> IO () +#ccall ekSepCMA_setOptimizer, Ptr -> Ptr -> IO () +#ccall ekSepCMA_sigma, Ptr -> CDouble +#ccall ekCSA_setSigma, Ptr -> CDouble -> CDouble -> IO () +#ccall ekCSA_setOptimizer, Ptr -> Ptr -> IO () +#ccall ekCSA_sigma, Ptr -> CDouble diff --git a/src/Bindings/eskit_support.c b/src/Bindings/eskit_support.c new file mode 100644 --- /dev/null +++ b/src/Bindings/eskit_support.c @@ -0,0 +1,22 @@ +#include +#include + +#define MKDEL_DATA(type, paramTypes, params) \ + type*\ + mk##type(paramTypes)\ + {\ + type* eo = malloc(sizeof(type)); \ + type##_init(eo, params); \ + return eo; \ + } \ + void \ + del##type(type* eo) \ + { \ + type##_destroy(eo); \ + free(eo); \ + } + +MKDEL_DATA(ekOptimizer, size_t N, N) +MKDEL_DATA(ekCMA, size_t N, N) +MKDEL_DATA(ekSepCMA, size_t N, N) +MKDEL_DATA(ekCSA, size_t N, N)