text:
{-# LANGUAGE DeriveGeneric,
             MultiParamTypeClasses,
             FlexibleInstances,
             FlexibleContexts,
             GeneralizedNewtypeDeriving,
             FunctionalDependencies,
             UndecidableInstances,
             TypeSynonymInstances,
             TypeFamilies
#-}
{-# OPTIONS_GHC -Wno-unused-top-binds #-}

module Main where

import Control.Exception.Safe (assert, MonadThrow)
import Control.Monad
import Control.Monad.Loops
import Control.Monad.IO.Class
import Control.Monad.Primitive
import Control.Monad.RWS.Class
import Control.Monad.Trans.Class
import Control.Monad.Trans.Reader (ReaderT(..))
import Control.Monad.Trans.RWS (RWST, evalRWST)
import Data.List (maximumBy, sort)
import Data.Ord (comparing)
import Data.Vector ((!), Vector)
import Data.IORef

import Numeric.Ketamine.Agent
import Numeric.Ketamine.Environment
import Numeric.Ketamine.Episode
import qualified Numeric.Ketamine.Outcome as O

import qualified Control.Monad.Trans.Reader as TR
import qualified Data.DList as D
import qualified Data.Map.Strict as Map
import qualified Data.Vector as V
import qualified Statistics.Distribution as Dist
import qualified Statistics.Distribution.Normal as N
import qualified System.Random.MWC as R

-- | The slot machine index whose arm will be pulled
type Action = Int
type Reward = Double
type Outcome = O.Outcome Action Reward

-- TODO: separate transformer stacks and log in Episode
type SharedState = RWST (IORef BanditState) (D.DList Outcome) EnvState IO
type Environment = EnT Action Outcome SharedState
type Agent = AgT Action Outcome SharedState

-------------------------------------------------------------------------------
-- | Episode

episodeLength :: Int
episodeLength = 100

narms :: Int
narms = 10

main :: IO ()
main = do
    casino <- envState narms
    bandit <- defaultBanditState narms
    let ep :: SharedState ()
        ep = runEpisode $ stepEnv >~> stepBandit 0.1 $ 9
    res <- snd <$> evalRWST ep bandit casino
    mapM_ print $ D.toList res

-------------------------------------------------------------------------------
-- | Environment
--
data EnvState =
    EnvState { arms :: Vector N.NormalDistribution , gen :: R.GenIO, pulls :: Int }

instance Show EnvState where
    show c = "EnvState" ++
      "{ means = " ++
        show (fmap Dist.mean . V.toList $ arms c) ++
          ", pulls = " ++ show (pulls c) ++ " }"

mkEnvState :: Int -> [Double] -> IO EnvState
mkEnvState n vars = do
    gen <- R.createSystemRandom
    means <- replicateM n $ R.uniform gen
    let arms = V.fromList $ zipWith N.normalDistr (reverse $ sort means) vars
    return $ EnvState arms gen 0

-- | Default config of a n-armed bandit
envState :: Int -> IO EnvState
envState n = mkEnvState n $ replicate n 1.0

genContVar :: Dist.ContGen d => d -> Environment Double
genContVar d = do
    g <- gets gen
    liftIO $ Dist.genContVar d g

stepEnv :: Action -> Environment ()
stepEnv = loop where
    loop act = do
        rwd <- genContVar =<< (! act) . arms <$> get
        modify $ \(EnvState a g p) -> EnvState a g (p+1)
        r <- ask
        liftIO $ modifyIORef r $ addStats act (Stats act rwd)
        n <- gets pulls
        if n >= episodeLength
            then return ()
            else outcome (O.Outcome act rwd) >>= loop

------------------------------------------------------------------------------
-- | Bandit

-- Statistics observed for a particular candidate.
data Stats = Stats
    { armCount :: !Int  -- ^ Number of times this candidate was observed
    , armTotal :: !Double -- ^ Total reward over all observations
    } deriving (Show, Eq)

instance Semigroup Stats where
    s <> s' = Stats { armCount = armCount s + armCount s'
                    , armTotal = armTotal s + armTotal s' }

instance Monoid Stats where
    mempty = Stats 0 0

-- | A record of statistics for all possibilities
type BanditState = Map.Map Action Stats

-- | Average reward over all observations for that arm.
mean :: Stats -> Reward
mean (Stats 0 _) = 0
mean ss = armTotal ss / toEnum (armCount ss)

addStats :: Action -> Stats -> BanditState -> BanditState
addStats = Map.insertWith (<>)

epsilonGreedy :: (R.Variate e, Ord e, Ord r) => [(a, r)] -> e -> SharedState a
epsilonGreedy acts =
    epsilonChoice (fst $ maximumBy (comparing snd) acts) acts

epsilonChoice :: (R.Variate e, Ord e) => a -> [(a, r)] -> e -> SharedState a
epsilonChoice a acts eps = do
    g <- gets gen
    compare eps <$> R.uniform g >>= \case
        LT -> pure a
        _  -> do
            i <- R.uniformR (0, length acts) g
            pure . fst . head $ drop (i-1) acts

banditState :: Int -> Stats -> IO (IORef BanditState)
banditState n s = do
    r <- newIORef mempty
    let init = Map.fromList $ take n $ zip [0..] (repeat s)
    writeIORef r init
    return r

defaultBanditState :: Int -> IO (IORef BanditState)
defaultBanditState n = banditState n mempty

stepBandit :: Float -> Outcome -> Agent ()
stepBandit eps = loop where
    loop o@(O.Outcome act rwd) = do
        tell . pure $ o
        r <- ask
        liftIO $ modifyIORef r $ addStats act (Stats act rwd)
        s <- liftIO $ readIORef r
        let rwds = Map.map mean s
        a <- lift $ epsilonGreedy (Map.toList rwds) eps
        action a >>= loop

max_stars_repo_name: cmk/dopamine

max_issues_repo_name: cmk/dopamine

max_forks_repo_name: cmk/dopamine


text:
-- |
-- Module      : Cartesian.Internal.Types
-- Description :
-- Copyright   : (c) Jonatan H Sundqvist, 2015
-- License     : MIT
-- Maintainer  : Jonatan H Sundqvist
-- Stability   : experimental|stable
-- Portability : POSIX (not sure)
--

-- Created October 31 2015

-- TODO | - Use TemplateHaskell (?)
--        - Strictness
--        - Performance, inlining

-- SPEC | -
--        -

------------------------------------------------------------------------------------------------------------------------------------------------------
-- GHC Pragmas
------------------------------------------------------------------------------------------------------------------------------------------------------
{-# LANGUAGE TemplateHaskell        #-}
{-# LANGUAGE MultiParamTypeClasses  #-}
{-# LANGUAGE FunctionalDependencies #-}
{-# LANGUAGE RankNTypes             #-}
{-# LANGUAGE FlexibleInstances      #-}

------------------------------------------------------------------------------------------------------------------------------------------------------
-- API
------------------------------------------------------------------------------------------------------------------------------------------------------
module Cartesian.Internal.Types where

------------------------------------------------------------------------------------------------------------------------------------------------------
-- We'll need these
------------------------------------------------------------------------------------------------------------------------------------------------------
import Control.Lens (Simple, Lens, lens)

import Data.Complex (Complex(..))

import Linear.V1
import Linear.V2
import Linear.V3
import Linear.V4

------------------------------------------------------------------------------------------------------------------------------------------------------
-- Types
------------------------------------------------------------------------------------------------------------------------------------------------------

-- Synonyms ------------------------------------------------------------------------------------------------------------------------------------------

-- TODO: Add some aliased lenses for these aliased types (?)

-- | A lens focusing on a single [vector-]component in a BoundingBox
type BoxLens v v' f f' = Lens (BoundingBox (v f)) (BoundingBox (v' f')) f f'

-- | An axis represented as (begin, length)
type Axis a = (a, a)

-- | A vector where each component represents a single axis (cf. 'Axis')
type Axes v a = v (Axis a)

-- |
-- type Domain

-- |
-- TODO: Rename (eg. 'Shape') (?)
type Polygon m v f = m (v f)

-- | Coordinate system wrappers
newtype Normalised v = Normalised { absolute   :: v } --
newtype Absolute   v = Absoloute  { normalised :: v } --

-- Types ---------------------------------------------------------------------------------------------------------------------------------------------

-- |
-- TODO: Anchors (eg. C, N, S, E W and combinations thereof, perhaps represented as relative Vectors)
-- TODO: Define some standard instances (eg. Functor, Applicative)
data BoundingBox v = BoundingBox { cornerOf :: v, sizeOf :: v } deriving (Show, Eq)

-- |
-- TODO: Use record (eg. from, to) (?)
data Line v = Line v v deriving (Show, Eq)

-- |
data Linear f = Linear { interceptOf :: f, slopeOf :: f } deriving (Show, Eq)

-- |
data Inclusivity r = Inclusive r | Exclusive r -- TODO: Rename (?)
data Interval r    = Interval (Inclusivity r) (Inclusivity r)

-- |
-- TODO: Use existing type instead (?)
-- data Side = SideLeft | SideRight | SideTop | SideBottom

-- Classes -------------------------------------------------------------------------------------------------------------------------------------------

-- TODO: How do you generate lenses for non-record types (?)
class HasX a f | a -> f where { x :: Simple Lens a f }
class HasY a f | a -> f where { y :: Simple Lens a f }
class HasZ a f | a -> f where { z :: Simple Lens a f }

-- Instances -----------------------------------------------------------------------------------------------------------------------------------------

-- TODO: Generate these instances with TemplateHaskell (?)

instance HasX (V1 f) f where
  x = lens (\(V1 x') -> x') (\_ x' -> V1 x')

instance HasX (V2 f) f where
  x = lens (\(V2 x' _) -> x') (\(V2 _ y') x' -> V2 x' y')

instance HasY (V2 f) f where
  y = lens (\(V2 _ y') -> y') (\(V2 x' _) y' -> V2 x' y')

instance HasX (V3 f) f where
  x = lens (\(V3 x' _ _) -> x') (\(V3 _ y' z') x' -> V3 x' y' z')

instance HasY (V3 f) f where
  y = lens (\(V3 _ y' _) -> y') (\(V3 x' _ z') y' -> V3 x' y' z')

instance HasZ (V3 f) f where
  z = lens (\(V3 _ _ z') -> z') (\(V3 x' y' _) z' -> V3 x' y' z')

instance HasX (V4 f) f where
  x = lens (\(V4 x' _ _ _) -> x') (\(V4 _ y' z' w') x' -> V4 x' y' z' w')

instance HasY (V4 f) f where
  y = lens (\(V4 _ y' _ _) -> y') (\(V4 x' _ z' w') y' -> V4 x' y' z' w')

instance HasZ (V4 f) f where
  z = lens (\(V4 _ _ z' _) -> z') (\(V4 x' y' _ w') z' -> V4 x' y' z' w')

instance HasX (Complex f) f where
  x = lens (\(x':+_) -> x') (\(_:+y') x' -> x':+y')

instance HasY (Complex f) f where
  y = lens (\(_':+y') -> y') (\(x':+_) y' -> x':+y')

ext: hs

max_stars_repo_name: SwiftsNamesake/Cartesian

max_issues_repo_name: SwiftsNamesake/Cartesian

max_forks_repo_name: SwiftsNamesake/Cartesian


text:
module AugLag where

import Numeric.LinearAlgebra
import Numeric.NLOPT

{-
See https://www.mat.univie.ac.at/~neum/glopt/coconut/Benchmark/Library1_new_v1.html

Alkyl problem in AMPL format:
var x2 := 1.745, >= 0, <= 2;
var x3 := 1.2, >= 0, <= 1.6;
var x4 := 1.1, >= 0, <= 1.2;
var x5 := 3.048, >= 0, <= 5;
var x6 := 1.974, >= 0, <= 2;
var x7 := 0.893, >= 0.85, <= 0.93;
var x8 := 0.928, >= 0.9, <= 0.95;
var x9 := 8, >= 3, <= 12;
var x10 := 3.6, >= 1.2, <= 4;
var x11 := 1.45, >= 1.45, <= 1.62;
var x12 := 1, >= 0.99, <= 1.01010101010101;
var x13 := 1, >= 0.99, <= 1.01010101010101;
var x14 := 1, >= 0.9, <= 1.11111111111111;
var x15 := 1, >= 0.99, <= 1.01010101010101;

minimize obj:  - 6.3*x5*x8 + 5.04*x2 + 0.35*x3 + x4 + 3.36*x6;

subject to
e2:  - 0.819672131147541*x2 + x5 - 0.819672131147541*x6 = 0;
e3: 0.98*x4 - x7*(0.01*x5*x10 + x4) = 0;
e4:  - x2*x9 + 10*x3 + x6 = 0;
e5: x5*x12 - x2*(1.12 + 0.13167*x9 - 0.0067*x9*x9) = 0;
e6: x8*x13 - 0.01*(1.098*x9 - 0.038*x9*x9) - 0.325*x7 - 0.57425 = 0;
e7: x10*x14 + 22.2*x11 - 35.82 = 0;
e8: x11*x15 - 3*x8 + 1.33 = 0;
-}
alkyl :: [Double] -> Either Result Solution
alkyl x0 =
  let lb = fromList [0, 0, 0, 0, 0, 0.85, 0.9, 3, 1.2, 1.45, 0.99, 0.99, 0.9, 0.99]
      ub = fromList [2, 1.6, 1.2, 5, 2, 0.93, 0.95, 12, 4, 1.62, 1.01010101010101, 1.01010101010101, 1.11111111111111, 1.01010101010101]
      objf x = - 6.3 * x ! 3 * x ! 6 + 5.04 * x ! 0 + 0.35 * x ! 1 + x ! 2 + 3.36 * x ! 4
      stop = ObjectiveRelativeTolerance 1e-12 :| []
      eqCons =
        [ EqualityConstraint
            { eqConstraintFunctions = Scalar (\x -> -0.819672131147541 * x ! 0 + x ! 3 - 0.819672131147541 * x ! 4),
              eqConstraintTolerance = 1e-10
            },
          EqualityConstraint
            { eqConstraintFunctions = Scalar (\x -> 0.98 * x ! 2 - x ! 5 * (0.01 * x ! 3 * x ! 8 + x ! 2)),
              eqConstraintTolerance = 1e-10
            },
          EqualityConstraint
            { eqConstraintFunctions = Scalar (\x -> - x ! 0 * x ! 7 + 10 * x ! 1 + x ! 4),
              eqConstraintTolerance = 1e-10
            },
          EqualityConstraint
            { eqConstraintFunctions = Scalar (\x -> x ! 3 * x ! 10 - x ! 0 * (1.12 + 0.13167 * x ! 7 - 0.0067 * x ! 7 * x ! 7)),
              eqConstraintTolerance = 1e-10
            },
          EqualityConstraint
            { eqConstraintFunctions = Scalar (\x -> x ! 6 * x ! 11 - 0.01 * (1.098 * x ! 7 - 0.038 * x ! 7 * x ! 7) - 0.325 * x ! 5 - 0.57425),
              eqConstraintTolerance = 1e-10
            },
          EqualityConstraint
            { eqConstraintFunctions = Scalar (\x -> x ! 8 * x ! 12 + 22.2 * x ! 9 - 35.82),
              eqConstraintTolerance = 1e-10
            },
          EqualityConstraint
            { eqConstraintFunctions = Scalar (\x -> x ! 9 * x ! 13 - 3 * x ! 6 + 1.33),
              eqConstraintTolerance = 1e-10
            }
        ]
      problem = AugLagProblem eqCons [] (AUGLAG_EQ_LOCAL subproblem)
      subproblem = LocalProblem (fromIntegral $ size lb) stop algorithm
      algorithm = SBPLX objf [LowerBounds lb, UpperBounds ub] Nothing
   in minimizeAugLag problem (fromList x0)


text:
{-# LANGUAGE TemplateHaskell #-}
{-# LANGUAGE TypeSynonymInstances #-}
{-# LANGUAGE FlexibleInstances #-}

module FEECa.FiniteElementTest where

import Control.Monad    ( liftM )

import FEECa.Internal.Form hiding ( inner )
import FEECa.Internal.Spaces
import qualified FEECa.Internal.MultiIndex  as MI
import qualified FEECa.Internal.Vector      as V
import qualified FEECa.Internal.Simplex     as S

import FEECa.Utility.Combinatorics
import FEECa.Utility.Print
import FEECa.Utility.Utility

import FEECa.FiniteElementSpace
import qualified FEECa.Polynomial                 as P
import qualified FEECa.PolynomialDifferentialForm as DF
import qualified FEECa.Bernstein                  as B

import FEECa.Utility.Test
import FEECa.Internal.SimplexTest

import qualified Test.QuickCheck                 as Q
import           Test.QuickCheck  ( (==>) )

import qualified Numeric.LinearAlgebra.HMatrix    as M

--------------------------------------------------------------------------------
-- Random Finite Element Space
--------------------------------------------------------------------------------

max_n = 3
max_r = 3
max_k = 3

arbitraryPL :: (Int -> Int -> Simplex -> FiniteElementSpace)
            -> Q.Gen FiniteElementSpace
arbitraryPL cpl = do
  [n,r,k] <- mapM Q.choose [(1,max_n), (0,max_r), (0,max_k)]
  liftM (cpl r k) (arbitrarySimplex n)

arbitraryPrLk :: Q.Gen FiniteElementSpace
arbitraryPrLk = arbitraryPL PrLk

arbitraryPrmLk :: Q.Gen FiniteElementSpace
arbitraryPrmLk = arbitraryPL PrmLk

instance Q.Arbitrary FiniteElementSpace where
  arbitrary = arbitraryPrLk

n = 4

--------------------------------------------------------------------------------
-- Whitney Forms
--------------------------------------------------------------------------------

data WhitneyTest = WhitneyTest Simplex
  deriving Show

instance Q.Arbitrary WhitneyTest where
  arbitrary = do k <- Q.choose (1,n)
                 liftM WhitneyTest (arbitrarySubsimplex k n)

prop_whitney_integral :: WhitneyTest -> Bool
prop_whitney_integral (WhitneyTest t) =
    abs (DF.integrate t (whitneyForm t [0..k]))
      `eqNum` (1 / fromInt (factorial k))
  where k = S.topologicalDimension t

--------------------------------------------------------------------------------
-- Psi Forms
--------------------------------------------------------------------------------

data PsiTest = PsiTest Simplex Simplex MI.MultiIndex Vector
  deriving Show

instance Q.Arbitrary PsiTest where
  arbitrary = do
    k <- Q.choose (1,n)
    t <- arbitrarySimplex n
    f <- Q.elements $ S.subsimplices t k
    r <- Q.choose (0,10)
    mi <- liftM (MI.extend n (S.sigma f)) (arbitraryMI (k+1) r)
    liftM (PsiTest t f mi) (arbitraryVector n)

prop_psi :: PsiTest -> Q.Property
prop_psi (PsiTest t f mi v) =
    MI.degree mi > 0 ==> all (eqNum 0.0 . evaluate v) ps
  where ps  = [DF.apply (psi' t f mi i) [l] | i <- is, l <- axs]
        axs = projectionAxes t f mi
        is  = S.sigma f
        k   = S.topologicalDimension f

convexCombination :: Simplex -> MI.MultiIndex -> Vector
convexCombination t mi = sumV (zipWith sclV mi' vs)
  where vs    = S.vertices t
        zero  = zeroV (head vs)
        (l,r) = (MI.toList mi, fromInt $ MI.degree mi)
        mi'   = [ let i = fromInt (l !! s)
                  in i / r  -- XXX: mulInv???...
                | s <- S.sigma t]

projectionAxes :: Simplex -> Simplex -> MI.MultiIndex -> [Vector]
projectionAxes t f mi = map (subV xmi) (S.complement t f)
  where xmi = convexCombination f mi

prop_basis :: FiniteElementSpace -> Q.Property
prop_basis s =
    length bs > 0 ==> length bs == dim s && linearIndependent DF.inner bs
  where bs  = basis s
        -- XXX: why the normalisation?? ... removed for now, seems unnecessary
        --      ==> we use bs instead
        --bs' = map (\x -> fmap (sclV (mulInv (sqrt (DF.inner x x)))) x) bs
        --n   = vspaceDim s

--------------------------------------------------------------------------------
-- PrmLk Form associated to a face.
--------------------------------------------------------------------------------

linearIndependent :: Module v => (v -> v -> Double) -> [v] -> Bool
linearIndependent f bs = M.rank mat == n
  where es  = M.eigenvaluesSH' mat
        -- TODO: changed eigenvaluesSH to eigenvaluesSH' for loading; check!
        mat = M.matrix n [ f omega eta | omega <- bs, eta <- bs ]
        n   = length bs

--------------------------------------------------------------------------------
-- DoFs of PrLk space
--------------------------------------------------------------------------------

prop_dof_basis :: FiniteElementSpace -> Q.Property
prop_dof_basis s = dim s > 0 && degree s > 0 ==> M.rank mat == n
  where mat = M.matrix (n) [d b | d <- dofs s, b <- basis s]
        n   = dim s

prop_ndofs :: FiniteElementSpace -> Q.Property
prop_ndofs s = degree s > 0 ==> length (dofs s) == sum [nDofs s k | k <- [0..n]]
  where n = S.topologicalDimension $ simplex s

return []
testFiniteElement = $quickCheckWithAll

space = PrmLk 3 0 (S.referenceSimplex 3)

max_stars_repo_name: Airini/FEECa

max_issues_repo_name: Airini/FEECa

max_forks_repo_name: Airini/FEECa


text:
---------------------------------------------------------------------------------
-- |
-- Utility functions for matrix operations
---------------------------------------------------------------------------------
module Utils where

import           Numeric.LinearAlgebra (Container, Element, Matrix, R, cmap,
                                        sumElements)
import           Prelude               (Num, Ord, (**), (.), (*), (-), abs, flip,
                                        max, signum)

frobeniusNorm :: Matrix R -> R
frobeniusNorm = flip (**) 0.5 . sumElements . cmap ((flip (**) 2) . abs)

absL1Norm :: Matrix R -> R
absL1Norm = sumElements . cmap abs

softThreshold
  :: (Ord e, Num e, Element e, Container c e)
  => c e
  -> e
  -> c e
softThreshold xs mu' = cmap f xs
  where f x = (signum x) * max ((abs x) - mu') 0

max_stars_repo_name: jamesthompson/robust-pca

max_issues_repo_name: jamesthompson/robust-pca

max_forks_repo_name: jamesthompson/robust-pca


text:
------------------------------------------------------------
--
-- ShowExample.hs
-- Code sample accompanying topic 1.3.4 "Type classes"
-- See README.md for details
--
-- Fundamentals of Practical Haskell Programming
-- By Richard Cook
--
------------------------------------------------------------

import Data.Complex
import Data.Ratio

main :: IO ()
main = do
    print 5
    print 5.0
    print "five"
    print (5 :+ 6)
    print (5 % 6)

max_stars_repo_path: haskell/original/Fundamentals of Practical Haskell Programming [Video]/Volume_01_Section_03_Topic_04/ShowExample.hs


text:
{-# LANGUAGE DataKinds #-}
{-# LANGUAGE FlexibleContexts #-}
{-# LANGUAGE GADTs #-}
{-# LANGUAGE RankNTypes #-}
{-# LANGUAGE ScopedTypeVariables #-}
{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE TypeOperators #-}
{-# LANGUAGE ViewPatterns #-}

{-# OPTIONS_GHC -fplugin GHC.TypeLits.KnownNat.Solver #-}
{-# OPTIONS_GHC -fplugin GHC.TypeLits.Normalise #-}

{-|
Module      : Numeric.Proximal.Static
Description : Solve proximal problems using FISTA
Copyright   : (c) Ryan Orendorff, 2020
License     : BSD3
Stability   : experimental
-}
module Numeric.Proximal.Static
    ( fista
    , fistaCost
    , indicatorBox
    , l1Cost
    , l1Prox
    , lsqLipschitz
    , proxIndicatorBox
    , quadraticLipschitz
    , residuals
    )
where

import           GHC.TypeLits
import           Numeric.LinearAlgebra.Static
import           Numeric.LinearAlgebra.Static.Util
import           Data.Maybe
import           Data.Tuple.HT                  ( fst3 )

------------------------------------------------------------------
--  Lipschitz calculations on actual matrixes by max eigenvalue --
------------------------------------------------------------------

-- TODO: Implement solution for extracting k largest eigenvalues directly.

-- Note that these are not particularly efficient as each eigenvalue must be extracted and then searched through to find the maximum.
-- | Find the Lipschitz constant for a least squares problem.
lsqLipschitz :: (KnownNat m, KnownNat n) => L m n -> Double
lsqLipschitz a = foldR1 max . eigenvalues $ mTm a

-- | Find the Lipschitz constant for a quadratic program.
quadraticLipschitz :: (KnownNat n) => Sym n -> Double
quadraticLipschitz = foldR1 max . eigenvalues

---------------------------------------------------------
--                  Proximal Functions                 --
---------------------------------------------------------

-- | Returns the cost function for the box indicator function.
-- I(l, u) = {0 if l < xᵢ < u ∀ xᵢ ∈ x
--            ∞ otherwise}
indicatorBox :: Double -> Double -> (forall n . KnownNat n => R n -> Double)
indicatorBox l u =
    foldR (\a b -> a + (if (b > u) || (b < l) then 1 / 0 else 0)) 0

-- | Proximal operator for the box indicator function
proxIndicatorBox :: Double -> Double -> (forall n . KnownNat n => R n -> R n)
proxIndicatorBox l u = dvmap bound
  where
    bound x = case (x < l, x > u) of
        (True, _   ) -> l
        (_   , True) -> u
        _            -> x

-- | Fast Iterative Soft Thresholding Algorithm
-- Used to solve problems of the form F(x) = f(x) + g(x), where
-- f is a smooth convex function with grad F being L-Lipschitz continuous and
-- g is a (potentially non-smooth) convex function. Common problems solved
-- with FISTA include min ||Ax - b|| st. x \in C or min ||Ax - b|| + ||x||_1
fista
    :: (KnownNat n)
    => (R n -> R n)    -- ^ gradient of f
    -> (R n -> R n)    -- ^ proximal operator of g
    -> Double          -- ^ Lipschitz constant
    -> Maybe (R n)     -- ^ Initial guess at a solution
    -> [R n]           -- ^ Output cost function and optimal solution over time
fista grad_f prox_g lipschitz (fromMaybe (konst 0) -> x0) = map fst3
    $ iterate go (x0, x0, 1)
  where
        -- Sometimes update is called p_L
    gradient_step = \x -> x - konst (1 / lipschitz) * grad_f x
    prox_step     = prox_g
    update        = prox_step . gradient_step

    -- Algorithm here
    go (x, y, t) =
        let
            x_update = update y
            t_update = (1 + sqrt (1 + 4 * t ^^ (2 :: Int))) / 2
            y_update =
                x_update + konst ((t - 1) / t_update) * (x_update - x)
        in
            (x_update, y_update, t_update)

-- | Residuals for a list of vectors.
residuals :: forall n . KnownNat n => [R n] -> [Double]
residuals x = zipWith (\y yn -> norm_2 (y - yn)) x (tail x)

-----------------------------------------------------------------
--                        Cost Functions                       --
-----------------------------------------------------------------

-- | Takes the most smooth convex function \(f\) and the convex (potentially
-- non-smooth) function \(g\) and calculates the cost of a given input.
fistaCost :: (R n -> Double) -> (R n -> Double) -> R n -> Double
fistaCost f g x = f x + g x

-- | Cost function for an ℓ₁ term
l1Cost :: (KnownNat n) => Double -> R n -> Double
l1Cost lambda = (lambda *) . norm_1

-- | Proximal equation a the ℓ₁ term
l1Prox :: (KnownNat n) => Double -> R n -> R n
l1Prox lambda = dvmap f where f xi = signum xi * max (abs xi - lambda) 0

max_stars_repo_name: ryanorendorff/convex

max_issues_repo_name: ryanorendorff/convex

max_forks_repo_name: ryanorendorff/convex


text:
{-# LANGUAGE FlexibleInstances #-}

-- | 'Graphable' Instances for lists of real or complex floating point numbers.
--   Users FlexibleInstances.
module System.Console.Ansigraph.Internal.FlexInstances where

import System.Console.Ansigraph.Core

import Data.Complex

-- | 1-dimensional real vector graph.
instance Graphable [Double] where
  graphWith = displayRV
  graphHeight _ = 2

-- | 1-dimensional complex vector graph.
instance Graphable [Complex Double] where
  graphWith = displayCV
  graphHeight _ = 4

-- | 2-dimensional real matrix graph.
instance Graphable [[Double]] where
  graphWith = displayMat
  graphHeight = length

-- | 2-dimensional complex matrix graph.
instance Graphable [[Complex Double]] where
  graphWith = displayCMat
  graphHeight = length

max_stars_repo_name: fieldstrength/ansigraph

max_issues_repo_name: fieldstrength/ansigraph

max_forks_repo_name: fieldstrength/ansigraph


text:
{-# LANGUAGE ExistentialQuantification #-}

module Main where

import           System.Environment
import           System.IO
import           Control.Monad
import           Control.Monad.Except
import           Numeric
import           Data.Ratio
import           Data.Complex
import           Data.Array
import           Data.IORef
import           Text.ParserCombinators.Parsec
                                         hiding ( spaces )

data LispVal =
    Atom String
    | List [LispVal]
    | DottedList [LispVal] LispVal
    | Vector (Array Int LispVal)
    | Number Integer
    | Float Double
    | Ratio Rational
    | Complex (Complex Double)
    | String String
    | Character Char
    | Bool Bool
    | PrimitiveFunc ([LispVal] -> ThrowsError LispVal)
    | Func {
        params :: [String],
        vararg :: Maybe String,
        body :: [LispVal],
        closure :: Env
    }
    | IOFunc ([LispVal] -> IOThrowsError LispVal)
    | Port Handle
    -- deriving Eq

unwordList :: [LispVal] -> String
unwordList = unwords . map showVal

showVal :: LispVal -> String
showVal (String contents) = "\"" ++ contents ++ "\""
showVal (Atom   name    ) = name
showVal (Number num     ) = show num
showVal (Bool   True    ) = "#t"
showVal (Bool   False   ) = "#f"
showVal (List   contents) = "(" ++ unwordList contents ++ ")"
showVal (DottedList head tail) =
    "(" ++ unwordList head ++ " . " ++ showVal tail ++ ")"
showVal (Character     char) = show char
showVal (PrimitiveFunc _   ) = "<primitive>"
showVal (Func { params = args, vararg = varargs, body = body, closure = env })
    = "(lambda ("
        ++ unwords (map show args)
        ++ (case varargs of
               Nothing  -> ""
               Just arg -> " . " ++ arg
           )
        ++ ") ...)"
showVal (Port   _) = "<IO port>"
showVal (IOFunc _) = "<IO primitive>"

instance Show LispVal where
    show = showVal

data LispError =
    NumArgs Integer [LispVal]
    | TypeMismatch String LispVal
    | Parser ParseError
    | BadSpecialForm String LispVal
    | NotFunction String String
    | UnboundVar String String
    | Default String

showError :: LispError -> String
showError (UnboundVar     message varname) = message ++ ": " ++ varname
showError (BadSpecialForm message form   ) = message ++ ": " ++ show form
showError (NotFunction    message func   ) = message ++ ": " ++ show func
showError (NumArgs expected found) =
    "Expected " ++ show expected ++ " args; found values " ++ unwordList found
showError (TypeMismatch expected found) =
    "Invalid type: expected " ++ expected ++ ", found " ++ show found
showError (Parser  parseErr) = "Parse error at " ++ show parseErr
showError (Default message ) = message

instance Show LispError where
    show = showError

type ThrowsError = Either LispError

trapError action = catchError action (return . show)

extractValue :: ThrowsError a -> a
extractValue (Right val) = val

escapedChars :: Parser Char
escapedChars = do
    char '\\'
    x <- oneOf "\\\"nrt"
    return $ case x of
        '\\' -> x
        '"'  -> x
        'n'  -> '\n'
        'r'  -> '\r'
        't'  -> '\t'

parseString :: Parser LispVal
parseString = do
    char '"'
    x <- many $ escapedChars <|> noneOf "\"\\"
    char '"'
    return $ String x

parseBool :: Parser LispVal
parseBool = do
    _ <- char '#'
    (char 't' >> return (Bool True)) <|> (char 'f' >> return (Bool False))

parseAtom :: Parser LispVal
parseAtom = do
    first <- letter <|> symbol
    rest  <- many (letter <|> digit <|> symbol)
    return $ Atom (first : rest)

parseDecimal1 :: Parser LispVal
parseDecimal1 = many1 digit >>= (return . Number . read)

parseDecimal2 :: Parser LispVal
parseDecimal2 = do
    try $ string "#d"
    x <- many1 digit
    (return . Number . read) x

hex2digit :: String -> Integer
hex2digit = fst . head . readHex

parseHex :: Parser LispVal
parseHex = do
    try $ string "#x"
    x <- many1 hexDigit
    (return . Number . hex2digit) x

oct2digit :: String -> Integer
oct2digit x = fst $ readOct x !! 0

parseOct :: Parser LispVal
parseOct = do
    try $ string "#o"
    x <- many1 octDigit
    (return . Number . oct2digit) x

bin2digit :: String -> Integer
bin2digit = bin2digit' 0
  where
    bin2digit' digitInt "" = digitInt
    bin2digit' digitInt (x : xs) =
        let old = 2 * digitInt + (if x == '0' then 0 else 1)
        in  bin2digit' old xs

parseBin :: Parser LispVal
parseBin = do
    try $ string "#b"
    x <- many1 (oneOf "01")
    (return . Number . bin2digit) x

parseNumber :: Parser LispVal
parseNumber =
    parseDecimal1 <|> parseDecimal2 <|> parseHex <|> parseOct <|> parseBin

parseCharacter :: Parser LispVal
parseCharacter = do
    try $ string "#\\"
    value <- try (string "newline" <|> string "space") <|> do
        x <- anyChar
        notFollowedBy alphaNum
        return [x]
    return $ Character $ case value of
        "space"   -> ' '
        "newline" -> '\n'
        otherwise -> value !! 0

parseFloat :: Parser LispVal
parseFloat = do
    x <- many1 digit
    char '.'
    y <- many1 digit
    (return . Float . fst . head . readFloat) (x ++ "." ++ y)

parseRatio :: Parser LispVal
parseRatio = do
    x <- many1 digit
    char '/'
    y <- many1 digit
    (return . Ratio) ((read x) % (read y))

toDouble :: LispVal -> Double
toDouble (Float  f) = f
toDouble (Number n) = fromIntegral n

parseComplex :: Parser LispVal
parseComplex = do
    x <- (try parseFloat <|> parseNumber)
    char '+'
    y <- (try parseFloat <|> parseNumber)
    char 'i'
    (return . Complex) (toDouble x :+ toDouble y)

parseList :: Parser LispVal
parseList = liftM List $ sepBy parseExpr spaces

parseDottedList :: Parser LispVal
parseDottedList = do
    head <- endBy parseExpr spaces
    tail <- char '.' >> spaces >> parseExpr
    return $ DottedList head tail

parseVector :: Parser LispVal
parseVector = do
    arrayValues <- sepBy parseExpr spaces
    return $ Vector $ listArray (0, length arrayValues - 1) arrayValues

parseQuoted :: Parser LispVal
parseQuoted = do
    char '\''
    x <- parseExpr
    return $ List [Atom "quote", x]

parseQuasiQuoted :: Parser LispVal
parseQuasiQuoted = do
    char '`'
    x <- parseExpr
    return $ List [Atom "quasiquote", x]

parseUnQuote :: Parser LispVal
parseUnQuote = do
    char ','
    x <- parseExpr
    return $ List [Atom "unquote", x]

parseExpr :: Parser LispVal
parseExpr =
    parseAtom
        <|> parseString
        <|> try parseFloat
        <|> try parseRatio
        <|> try parseComplex
        <|> try parseNumber
        <|> try parseBool
        <|> try parseCharacter
        <|> parseQuoted
        <|> parseQuasiQuoted
        <|> parseUnQuote
        <|> try
                (do
                    string "#("
                    x <- parseVector
                    char ')'
                    return x
                )
        <|> do
                char '('
                x <- try parseList <|> parseDottedList
                char ')'
                return x

symbol :: Parser Char
symbol = oneOf "!$%&|*+-/:<=>?@^_~"

spaces :: Parser ()
spaces = skipMany1 space

numericBinop
    :: (Integer -> Integer -> Integer) -> [LispVal] -> ThrowsError LispVal
numericBinop op []            = throwError $ NumArgs 2 []
numericBinop op singleVal@[_] = throwError $ NumArgs 2 singleVal
numericBinop op params        = Number . foldl1 op <$> mapM unpackNum params

unpackNum :: LispVal -> ThrowsError Integer
unpackNum (Number n) = return n
unpackNum (String n) =
    let parsed = reads n
    in  if null parsed
            then throwError $ TypeMismatch "number" $ String n
            else return $ fst $ parsed !! 0
unpackNum (List [n]) = unpackNum n
unpackNum notNum     = throwError $ TypeMismatch "number" notNum

unpackStr :: LispVal -> ThrowsError String
unpackStr (String s) = return s
unpackStr (Number n) = return $ show n
unpackStr (Bool   b) = return $ show b
unpackStr notString  = throwError $ TypeMismatch "string" notString

unpackBool :: LispVal -> ThrowsError Bool
unpackBool (Bool b) = return b
unpackBool notBool  = throwError $ TypeMismatch "boolean" notBool

unaryOp :: (LispVal -> LispVal) -> [LispVal] -> ThrowsError LispVal
unaryOp f [v]    = return $ f v
unaryOp f []     = throwError $ NumArgs 1 []
unaryOp f params = throwError $ NumArgs 1 params

symbolp :: LispVal -> LispVal
symbolp (Atom _) = Bool True
symbolp _        = Bool False

stringp :: LispVal -> LispVal
stringp (String _) = Bool True
stringp _          = Bool False

numberp :: LispVal -> LispVal
numberp (Number _) = Bool True
numberp _          = Bool False

boolp :: LispVal -> LispVal
boolp (Bool _) = Bool True
boolp _        = Bool False

listp :: LispVal -> LispVal
listp (List _        ) = Bool True
listp (DottedList _ _) = Bool True
listp _                = Bool False

symbol2string :: LispVal -> LispVal
symbol2string (Atom s) = String s
symbol2string _        = String ""

string2symbol :: LispVal -> LispVal
string2symbol (String s) = Atom s
string2symbol _          = Atom ""

boolBinop
    :: (LispVal -> ThrowsError a)
    -> (a -> a -> Bool)
    -> [LispVal]
    -> ThrowsError LispVal
boolBinop unpacker op args = if length args /= 2
    then throwError $ NumArgs 2 args
    else do
        left  <- unpacker $ args !! 0
        right <- unpacker $ args !! 1
        return $ Bool $ left `op` right

numBoolBinop :: (Integer -> Integer -> Bool) -> [LispVal] -> ThrowsError LispVal
numBoolBinop = boolBinop unpackNum

strBoolBinop :: (String -> String -> Bool) -> [LispVal] -> ThrowsError LispVal
strBoolBinop = boolBinop unpackStr

boolBoolBinop :: (Bool -> Bool -> Bool) -> [LispVal] -> ThrowsError LispVal
boolBoolBinop = boolBinop unpackBool

car :: [LispVal] -> ThrowsError LispVal
car [List (x : xs)        ] = return x
car [DottedList (x : xs) _] = return x
car [badArg               ] = throwError $ TypeMismatch "pair" badArg
car badArgList              = throwError $ NumArgs 1 badArgList

cdr :: [LispVal] -> ThrowsError LispVal
cdr [List (x : xs)        ] = return $ List xs
cdr [DottedList [_     ] x] = return x
cdr [DottedList (_ : xs) x] = return $ DottedList xs x
cdr [badArg               ] = throwError $ TypeMismatch "pair" badArg
cdr badArgList              = throwError $ NumArgs 1 badArgList

cons :: [LispVal] -> ThrowsError LispVal
cons [x1, List []            ] = return $ List [x1]
cons [x , List xs            ] = return $ List $ x : xs
cons [x , DottedList xs xlast] = return $ DottedList (x : xs) xlast
cons [x1, x2                 ] = return $ DottedList [x1] x2
cons badArgList                = throwError $ NumArgs 2 badArgList

eqv :: [LispVal] -> ThrowsError LispVal
eqv [Bool   arg1, Bool arg2  ] = return $ Bool $ arg1 == arg2
eqv [Number arg1, Number arg2] = return $ Bool $ arg1 == arg2
eqv [String arg1, String arg2] = return $ Bool $ arg1 == arg2
eqv [Atom   arg1, Atom arg2  ] = return $ Bool $ arg1 == arg2
eqv [DottedList xs x, DottedList ys y] =
    eqv [List $ xs ++ [x], List $ ys ++ [y]]
eqv [List arg1, List arg2] =
    return $ Bool $ (length arg1 == length arg2) && all eqvPair (zip arg1 arg2)
  where
    eqvPair (x1, x2) = case eqv [x1, x2] of
        Left  err        -> False
        Right (Bool val) -> val
eqv [_, _]     = return $ Bool False
eqv badArgList = throwError $ NumArgs 2 badArgList

data Unpacker = forall a. Eq a => AnyUnpacker (LispVal -> ThrowsError a)

unpackEquals :: LispVal -> LispVal -> Unpacker -> ThrowsError Bool
unpackEquals arg1 arg2 (AnyUnpacker unpacker) =
    do
            unpacked1 <- unpacker arg1
            unpacked2 <- unpacker arg2
            return $ unpacked1 == unpacked2
        `catchError` (const $ return False)

equal :: [LispVal] -> ThrowsError LispVal
equal [(List arg1), (List arg2)] =
    return $ Bool $ (length arg1 == length arg2) && all equalPair
                                                        (zip arg1 arg2)
  where
    equalPair (x1, x2) = case equal [x1, x2] of
        Left  err        -> False
        Right (Bool val) -> val
equal [(DottedList xs x), (DottedList ys y)] =
    equal [List $ xs ++ [x], List $ ys ++ [y]]
equal [arg1, arg2] = do
    primitiveEquals <- fmap or $ mapM
        (unpackEquals arg1 arg2)
        [AnyUnpacker unpackNum, AnyUnpacker unpackStr, AnyUnpacker unpackBool]
    eqvEquals <- eqv [arg1, arg2]
    return $ Bool $ (primitiveEquals || let (Bool x) = eqvEquals in x)
equal badArgList = throwError $ NumArgs 2 badArgList

stringLen :: [LispVal] -> ThrowsError LispVal
stringLen [(String s)] = return $ Number $ fromIntegral $ length s
stringLen [notString ] = throwError $ TypeMismatch "string" notString
stringLen badArgList   = throwError $ NumArgs 1 badArgList

stringRef :: [LispVal] -> ThrowsError LispVal
stringRef [(String s), (Number k)]
    | length s < k' + 1 = throwError $ Default "Out of bound error"
    | otherwise         = return $ Character $ s !! k'
    where k' = fromIntegral k
stringRef [(String s), notNum] = throwError $ TypeMismatch "number" notNum
stringRef [notString , _     ] = throwError $ TypeMismatch "string" notString
stringRef badArgList           = throwError $ NumArgs 2 badArgList

primitives :: [(String, [LispVal] -> ThrowsError LispVal)]
primitives =
    [ ("+"             , numericBinop (+))
    , ("-"             , numericBinop (-))
    , ("*"             , numericBinop (*))
    , ("/"             , numericBinop div)
    , ("mod"           , numericBinop mod)
    , ("quotient"      , numericBinop quot)
    , ("remainder"     , numericBinop rem)
    , ("symbol?"       , unaryOp symbolp)
    , ("string?"       , unaryOp stringp)
    , ("number?"       , unaryOp numberp)
    , ("bool?"         , unaryOp boolp)
    , ("list?"         , unaryOp listp)
    , ("symbol->string", unaryOp symbol2string)
    , ("string->symbol", unaryOp string2symbol)
    , ("="             , numBoolBinop (==))
    , ("<"             , numBoolBinop (<))
    , (">"             , numBoolBinop (>))
    , ("/="            , numBoolBinop (/=))
    , (">="            , numBoolBinop (>=))
    , ("<="            , numBoolBinop (<=))
    , ("&&"            , boolBoolBinop (&&))
    , ("||"            , boolBoolBinop (||))
    , ("string=?"      , strBoolBinop (==))
    , ("string<?"      , strBoolBinop (<))
    , ("string>?"      , strBoolBinop (>))
    , ("string<=?"     , strBoolBinop (<=))
    , ("string>=?"     , strBoolBinop (>=))
    , ("car"           , car)
    , ("cdr"           , cdr)
    , ("cons"          , cons)
    , ("eq?"           , eqv)
    , ("eqv?"          , eqv)
    , ("equal?"        , equal)
    , ("string-length" , stringLen)
    , ("string-ref"    , stringRef)
    ]

applyProc :: [LispVal] -> IOThrowsError LispVal
applyProc [func, List args] = apply func args
applyProc (func : args)     = apply func args

makePort :: IOMode -> [LispVal] -> IOThrowsError LispVal
makePort mode [String filename] = liftM Port $ liftIO $ openFile filename mode

closePort :: [LispVal] -> IOThrowsError LispVal
closePort [Port port] = liftIO $ hClose port >> (return $ Bool True)
closePort _           = return $ Bool False

readProc :: [LispVal] -> IOThrowsError LispVal
readProc []          = readProc [Port stdin]
readProc [Port port] = liftIO (hGetLine port) >>= liftThrows . readExpr

writeProc :: [LispVal] -> IOThrowsError LispVal
writeProc [obj]            = writeProc [obj, Port stdout]
writeProc [obj, Port port] = liftIO $ hPrint port obj >> (return $ Bool True)

readContents :: [LispVal] -> IOThrowsError LispVal
readContents [String filename] = fmap String $ liftIO $ readFile filename

load :: String -> IOThrowsError [LispVal]
load filename = (liftIO $ readFile filename) >>= liftThrows . readExprList

readAll :: [LispVal] -> IOThrowsError LispVal
readAll [String filename] = liftM List $ load filename

ioPrimitives :: [(String, [LispVal] -> IOThrowsError LispVal)]
ioPrimitives =
    [ ("apply"            , applyProc)
    , ("open-input-file"  , makePort ReadMode)
    , ("open-output-file" , makePort WriteMode)
    , ("close-input-port" , closePort)
    , ("close-output-port", closePort)
    , ("read"             , readProc)
    , ("write"            , writeProc)
    , ("read-contents"    , readContents)
    , ("read-all"         , readAll)
    ]

primitiveBindings :: IO Env
primitiveBindings =
    nullEnv
        >>= (  flip bindVars
            $  map (makeFunc IOFunc)        ioPrimitives
            ++ map (makeFunc PrimitiveFunc) primitives
            )
    where makeFunc constructor (var, func) = (var, constructor func)

apply :: LispVal -> [LispVal] -> IOThrowsError LispVal
apply (PrimitiveFunc func) args = liftThrows $ func args
apply (Func params varargs body closure) args =
    if num params /= num args && varargs == Nothing
        then throwError $ NumArgs (num params) args
        else
            (liftIO $ bindVars closure $ zip params args)
            >>= bindVarArgs varargs
            >>= evalBody
  where
    remainingArgs = drop (length params) args
    num           = toInteger . length
    evalBody env = liftM last $ mapM (eval env) body
    bindVarArgs arg env = case arg of
        Just argName -> liftIO $ bindVars env [(argName, List $ remainingArgs)]
        Nothing      -> return env
apply (IOFunc func) args = func args

type Env = IORef [(String, IORef LispVal)]

nullEnv :: IO Env
nullEnv = newIORef []

type IOThrowsError = ExceptT LispError IO

liftThrows :: ThrowsError a -> IOThrowsError a
liftThrows (Left  err) = throwError err
liftThrows (Right val) = return val

runIOThrows :: IOThrowsError String -> IO String
runIOThrows action = runExceptT (trapError action) >>= return . extractValue

isBound :: Env -> String -> IO Bool
isBound envRef var =
    readIORef envRef >>= return . maybe False (const True) . lookup var

getVar :: Env -> String -> IOThrowsError LispVal
getVar envRef var = do
    env <- liftIO $ readIORef envRef
    maybe (throwError $ UnboundVar "Getting an unbound variable" var)
          (liftIO . readIORef)
          (lookup var env)

setVar :: Env -> String -> LispVal -> IOThrowsError LispVal
setVar envRef var value = do
    env <- liftIO $ readIORef envRef
    maybe (throwError $ UnboundVar "Setting an unbound variable" var)
          (liftIO . (flip writeIORef value))
          (lookup var env)
    return value

defineVar :: Env -> String -> LispVal -> IOThrowsError LispVal
defineVar envRef var value = do
    alreadyDefined <- liftIO $ isBound envRef var
    if alreadyDefined
        then setVar envRef var value >> return value
        else liftIO $ do
            valueRef <- newIORef value
            env      <- readIORef envRef
            writeIORef envRef ((var, valueRef) : env)
            return value

bindVars :: Env -> [(String, LispVal)] -> IO Env
bindVars envRef bindings = readIORef envRef >>= extendEnv bindings >>= newIORef
  where
    extendEnv bindings env = liftM (++ env) (mapM addBinding bindings)
    addBinding (var, value) = do
        ref <- newIORef value
        return (var, ref)

makeFunc varargs env params body =
    return $ Func (map showVal params) varargs body env
makeNormalFunc = makeFunc Nothing
makeVarArgs = makeFunc . Just . showVal

eval :: Env -> LispVal -> IOThrowsError LispVal
eval env val@(String _                             ) = return val
eval env val@(Number _                             ) = return val
eval env val@(Bool   _                             ) = return val
eval env (    Atom   id                            ) = getVar env id
eval env (    List   [Atom "quote", val]           ) = return val
eval env (    List   [Atom "if", pred, conseq, alt]) = do
    result <- eval env pred
    case result of
        Bool True  -> eval env conseq
        Bool False -> eval env alt
        otherwise  -> throwError $ TypeMismatch "boolean" otherwise
eval env form@(List (Atom "cond" : clauses)) = if null clauses
    then throwError $ BadSpecialForm "no true clause in cond expression: " form
    else case head clauses of
        List [Atom "else", expr] -> eval env expr
        List [test, expr] ->
            eval env $ List
                [Atom "if", test, expr, List (Atom "cond" : tail clauses)]
        _ -> throwError $ BadSpecialForm "ill-formed cond expression" form
-- eval env form@(List (Atom "case" : key : clauses)) =
--     if null clauses
--         then throwError $ BadSpecialForm "no true clause in case expression: " form
--         else case head clauses of
--             List (Atom "else" : exprs) -> mapM (eval env) exprs >>= return . last
--             List ((List datums) : exprs) ->
--                 do
--                     result <- eval env key
--                     equality <- liftThrows $ mapM (\x -> eqv [result, x]) datums
--                     if Bool True `elem` equality
--                         then mapM (eval env) exprs >>= return . last
--                         else eval env $ List (Atom "case" : key : tail clauses)
--             _ -> throwError $ BadSpecialForm "ill-formed case expression" form
eval env (List [Atom "set!", Atom var, form]) =
    eval env form >>= setVar env var
eval env (List [Atom "define", Atom var, form]) =
    eval env form >>= defineVar env var
eval env (List (Atom "define" : List (Atom var : params) : body)) =
    makeNormalFunc env params body >>= defineVar env var
eval env (List (Atom "define" : DottedList (Atom var : params) varargs : body))
    = makeVarArgs varargs env params body >>= defineVar env var
eval env (List (Atom "lambda" : List params : body)) =
    makeNormalFunc env params body
eval env (List (Atom "lambda" : DottedList params varargs : body)) =
    makeVarArgs varargs env params body
eval env (List (Atom "lambda" : varargs@(Atom _) : body)) =
    makeVarArgs varargs env [] body
eval env (List [Atom "load", String filename]) =
    load filename >>= liftM last . mapM (eval env)
eval env (List (func : args)) = do
    func    <- eval env func
    argVals <- mapM (eval env) args
    apply func argVals
eval env badForm =
    throwError $ BadSpecialForm "Unrecognized special form" badForm

readOrThrow :: Parser a -> String -> ThrowsError a
readOrThrow parser input = case parse parser "lisp" input of
    Left  err -> throwError $ Parser err
    Right val -> return val

readExpr = readOrThrow parseExpr
readExprList = readOrThrow (endBy parseExpr spaces)

flushStr :: String -> IO ()
flushStr str = putStr str >> hFlush stdout

readPrompt :: String -> IO String
readPrompt prompt = flushStr prompt >> getLine

evalString :: Env -> String -> IO String
evalString env expr =
    runIOThrows $ liftM show $ (liftThrows $ readExpr expr) >>= eval env

evalAndPrint :: Env -> String -> IO ()
evalAndPrint env expr = evalString env expr >>= putStrLn

runOne :: [String] -> IO ()
runOne args = do
    env <- primitiveBindings
        >>= flip bindVars [("args", List $ map String $ drop 1 args)]
    (runIOThrows $ liftM show $ eval
            env
            (List [Atom "load", String (args !! 0)])
        )
        >>= hPutStrLn stderr

until_ :: Monad m => (a -> Bool) -> m a -> (a -> m ()) -> m ()
until_ pred prompt action = do
    result <- prompt
    if pred result
        then return ()
        else action result >> until_ pred prompt action

runRepl :: IO ()
runRepl =
    primitiveBindings
        >>= until_ (== "quit") (readPrompt "Lisp>>> ")
        .   evalAndPrint

main :: IO ()
main = do
    args <- getArgs
    if null args then runRepl else runOne $ args


text:
module ComplexFunctions where

import Data.Complex
import RealFunctions

fibonacci :: (Fractional a, Floating a, Num a) => a -> a
fibonacci z = ((phi ** z) - (cos (pi*z)) * (phi ** (-z))) / (sqrt 5)
    where phi = (1 + (sqrt 5)) / 2

weierstrass :: (Floating a) => Int -> a -> a -> a -> a
weierstrass limit a b x = sum $ map (\n -> (a^n) * cos (b^n * pi * x)) [0..limit]

foldPhase :: (RealFloat a) => a -> Complex a -> Complex a
foldPhase p z = mkPolar (magnitude z) (triangle (phase z) (p) (p))

foldMag :: (RealFloat a) => a -> Complex a -> Complex a
foldMag p z = mkPolar (triangle (magnitude z) p p) (phase z)

foldX :: (RealFloat a) => a -> Complex a -> Complex a
foldX p z = (x:+y)
    where x = triangle (realPart z) p p
          y = imagPart z

foldY :: (RealFloat a) => a -> Complex a -> Complex a
foldY p z = (x:+y)
    where x = realPart z
          y = triangle (imagPart z) p p

-- remaps x and y in [-1,1]
remap :: (RealFloat a) => Complex a -> a -> a -> Complex a
remap z w h = ((realPart z') / a) :+ ((imagPart z') / a)
    where a = 0.5 * (min w h)
          z'= z - (w/2 :+ h/2) -- centering

-- STOP CONDITIONS
isNotInfinite :: (RealFloat a, Fractional a ) => Complex a -> Bool
isNotInfinite z = not $ (isNaN . magnitude) z || (isInfinite . magnitude) z

magnitudeLimit :: (Ord a, RealFloat a, Fractional a) => a -> Complex a -> Bool
magnitudeLimit limit z = magnitude z < limit

diamond :: (RealFloat a, Fractional a) => Complex a -> a
diamond z = y + x
    where x = abs (realPart z)
          y = abs (imagPart z)

compareBy :: (Ord a, RealFloat a) => (Complex a -> a) ->  Complex a -> Complex a -> Ordering
compareBy func z1 z2 = compare (func z1) (func z2)
{--
hofstadterQ :: (Integral a) => a -> a
hofstadterQ 1 = 1
hofstadterQ 2 = 1
hofstadterQ n = (hofstadterQ (n - hofstadterQ (n - 1))) +
    (hofstadterQ (n - hofstadterQ (n - 2)))
--}

max_stars_repo_name: Vetii/Fractal

max_issues_repo_name: Vetii/Fractal

max_forks_repo_name: Vetii/Fractal


text:
{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE NoMonomorphismRestriction #-}
{-# LANGUAGE ScopedTypeVariables #-}
{-# LANGUAGE ConstraintKinds #-}
{-# LANGUAGE FlexibleContexts #-}

{-# OPTIONS -Wall #-}
module Statistics.Iteratee.Sample (
  minMaxNBy
, range
, mean
, harmonicMean
, variance
, stdDev
) where

import Statistics.Iteratee.Compat
import Control.Arrow
import Control.Monad
import Data.Iteratee as I

import Data.Heap as Heap

--  | /O(n)/ minMaxNBy. Calculate the 'n' highest and lowest elements of the
--  stream, according to the given priority function.
--  Returns /([minimum],[maximum]/ with the /(minimum,maximum)/
--  elements listed first.
minMaxNBy
  :: forall m prio s el. (Monad m, Ord prio, ListLikey s el)
  => Int
  -> (el -> prio)
  -> Iteratee s m ([(prio,el)],[(prio,el)])
minMaxNBy ns prio = finalize `liftM` I.foldl' step (Heap.empty,Heap.empty)
  where
    finalize :: (MaxPrioHeap prio el, MinPrioHeap prio el)
             -> ([(prio,el)],[(prio,el)])
    finalize = Heap.toDescList *** Heap.toDescList
    addHeap val = Heap.insert (prio val, val)
    step :: (MaxPrioHeap prio el, MinPrioHeap prio el) -> el
            -> (MaxPrioHeap prio el, MinPrioHeap prio el)
    step (!mins,!maxes) val = let sz = Heap.size mins
                                  adj hp = if sz >= ns
                                             then Heap.drop 1 hp
                                             else hp
                              in (adj $ addHeap val mins
                                 ,adj $ addHeap val maxes)
{-# INLINE minMaxNBy #-}

-- | /O(n)/ Range. The difference between the largest and smallest elements of
-- a stream.
range :: (Monad m, ListLikey s el, Num el, Ord el)
      => Iteratee s m el
range = finalize `liftM` minMaxNBy 1 id
  where
    finalize ([mins],[maxes]) = snd maxes - snd mins
    finalize _ = 0
{-# INLINE range #-}

-- | /O(n)/ Arithmetic mean.  Uses Welford's algorithm.
mean :: forall s m el. (Fractional el, Monad m, ListLikey s el)
     => Iteratee s m el
mean = fst `liftM` I.foldl' step (0,0)
  where
    step :: (el,Integer) -> el -> (el,Integer)
    step (!m,!n) x = let m' = m + (x-m) / fromIntegral n'
                         n' = n + 1
                     in  (m',n')
{-# INLINE mean #-}

-- | /O(n)/ Harmonic mean.
harmonicMean :: (Fractional el, Monad m, ListLikey s el) => Iteratee s m el
harmonicMean = finalize `liftM` I.foldl' step (0,0 :: Integer)
  where
    finalize (m,n) = fromIntegral n / m
    step (!m,!n) val = (m+(1/val),n+1)
{-# INLINE harmonicMean #-}

-- | /O(n)/ variance, using Knuth's algorithm.
var :: (Fractional el, Integral t, Monad m, ListLikey s el)
    => Iteratee s m (t, el, el)
var = I.foldl' step (0,0,0)
  where
    step (!n,!m,!s) x = let n' = n+1
                            m' = m+d/fromIntegral n'
                            s' = s+d* (x-m')
                            d  = x-m
                        in (n',m',s')
{-# INLINE var #-}

-- | /O(n)/ Maximum likelihood estimate of a sample's variance, using Knuth's
--   algorithm.
variance :: (Fractional b, Monad m, ListLikey s b) => Iteratee s m b
variance = finalize `liftM` var
  where
    finalize (n,_,s)
      | n > 1 = s / fromInteger n
      | otherwise = 0
{-# INLINE variance #-}

-- | /O(n) Standard deviation, using Knuth's algorithm.
stdDev :: (Floating b, Monad m, Functor m, ListLikey s b) => Iteratee s m b
stdDev = sqrt `liftM` variance
{-# INLINE stdDev #-}

max_stars_repo_name: JohnLato/iter-stats

max_issues_repo_name: JohnLato/iter-stats

max_forks_repo_name: JohnLato/iter-stats


text:
{-# LANGUAGE DataKinds #-}
{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE FlexibleContexts #-}
module Numeric.Layer where

import Dependent.Size
import Numeric.Vector.Sized
import Data.Kind(Type)

class Layer layer where
  type Inputs layer :: Size
  type Outputs layer :: Size
  type Tape layer :: Type
  type Gradient layer :: Type

  forward :: layer
          -> SizedArray (Inputs layer)
          -> (SizedArray (Outputs layer), Tape layer)

  backward :: layer
           -> Tape layer
           -> SizedArray (Outputs layer)
           -> (Gradient layer, SizedArray (Inputs layer))

  applyGradient :: layer -> Gradient layer -> layer

max_stars_repo_name: mixed-signals/mixed-signals

max_issues_repo_name: mixed-signals/mixed-signals

max_forks_repo_name: mixed-signals/mixed-signals


text:
{-# LANGUAGE GADTs #-}
module Math.IRT.MLE.Internal.Generic
  ( MLEResult (..)
  , generic_mleEst
  , logLike
  ) where

import Numeric.AD (Mode, Scalar)
import Numeric.AD.Halley

import Statistics.Distribution (Distribution, ContDistr)

import Math.IRT.Fisher
import Math.IRT.Internal.Distribution
import Math.IRT.Internal.LogLikelihood
import Math.IRT.MLE.Internal.MLEResult

generic_mleEst :: (Distribution d, ContDistr d, DensityDeriv d, LogLikelihood d) => [Bool] -> [d] -> Int -> Double -> MLEResult
generic_mleEst rs params steps vTheta =
    let est    = (!! steps) $ extremumNoEq (logLike rs params) vTheta
        fisher = fisherInfoObserved est rs params
    in case fisher of
         (FisherInfo _ t s) -> MLEResult est t s

logLike :: (Mode a, Floating a, Scalar a ~ Double, Distribution d, LogLikelihood d) => [Bool] -> [d] -> a -> a
logLike responses params vTheta =
    let logLik   = zipWith (\a b -> logLikelihood a b vTheta) responses params
        bmePrior = 1
    in sum logLik + log bmePrior

max_stars_repo_name: argiopetech/irt

max_issues_repo_name: argiopetech/irt

max_forks_repo_name: argiopetech/irt


text:
{-# LANGUAGE TypeApplications #-}
module ArrayFire.BLASSpec where

import ArrayFire    hiding (not)

import Data.Complex
import Test.Hspec

spec :: Spec
spec =
  describe "BLAS spec" $ do
    it "Should matmul two matrices" $ do
      (matrix @Double (2,2) [[2,2],[2,2]] `matmul` matrix @Double (2,2) [[2,2],[2,2]]) None None
        `shouldBe` matrix @Double (2,2) [[8,8],[8,8]]
    it "Should dot product two vectors" $ do
      dot (vector @Double 2 (repeat 2)) (vector @Double 2 (repeat 2)) None None
        `shouldBe`
           scalar @Double 8
    it "Should produce scalar dot product between two vectors as a Complex number" $ do
      dotAll (vector @Double 2 (repeat 2)) (vector @Double 2 (repeat 2)) None None
        `shouldBe`
           8.0 :+ 0.0
    it "Should take the transpose of a matrix" $ do
      transpose (matrix @Double (2,2) [[1,1],[2,2]]) False
        `shouldBe`
           matrix @Double (2,2) [[1,2],[1,2]]
    it "Should take the transpose of a matrix in place" $ do
      let m = matrix @Double (2,2) [[1,1],[2,2]]
      transposeInPlace m False
      m `shouldBe` matrix @Double (2,2) [[1,2],[1,2]]

max_stars_repo_name: danielkroeni/arrayfire-haskell

max_issues_repo_name: danielkroeni/arrayfire-haskell

max_forks_repo_name: danielkroeni/arrayfire-haskell