text:
{-# LANGUAGE QuasiQuotes #-}
module Myelin.PyNN.PyNNSpec where

import Control.Lens
import Control.Monad.State.Lazy
import Control.Monad.Trans.Except
import Data.Aeson
import Data.ByteString.Lazy.Char8
import Data.Either
import qualified Data.Map.Strict as Map
import Data.String.Interpolate
import Numeric.LinearAlgebra
import Text.Regex as Regex

import Myelin.Model
import Myelin.Neuron
import Myelin.PyNN.PyNN
import Myelin.SNN

import Test.Hspec

main :: IO ()
main = hspec spec

eval :: PyNNState a -> Either String a
eval s = evalState (runExceptT s) emptyPyNNModel

exec :: PyNNState a -> PyNNModel
exec s = execState (runExceptT s) emptyPyNNModel

spec :: Spec
spec = do
  describe "PyNN backend" $ do
    it "can return the reference for a node" $ do
      let tpe = if_cond_exp
      let pop = Population 2 tpe "p0" 0
      let expected = "p0"
      eval (pyNNNode pop) `shouldBe` Right expected
    it "can translate an IF_exp_cond population to Python" $ do
      let tpe = if_cond_exp
      let dictWithType = unpack $ encode tpe
      let dict = Regex.subRegex typeRegex dictWithType ""
      let pop = Population 2 tpe "p0" 0
      let expectedPopulation = "p0 = pynn.Population(2, pynn.IF_cond_exp(**" ++ dict ++ "))"
      let actualModel = exec (pyNNNode pop)
      actualModel ^. populations `shouldBe` Map.singleton "p0" expectedPopulation
    it "can return the reference for a layer" $ do
      let effect = Static Excitatory (AllToAll (BiasGenerator (Constant 0)) (WeightGenerator (Constant 1.0)))
      let nodeIn = SpikeSourcePoisson 0.1 0 0
      let nodeOut = Population 2 if_cond_exp "p1" 1
      let edge = DenseProjection effect nodeIn nodeOut
      let expected = "layer0"
      eval (pyNNEdge edge) `shouldBe` Right expected
    it "can translate a dense edge to Python" $ do
      let effect = Static Excitatory (AllToAll (BiasGenerator (Constant 0)) (WeightGenerator (Constant 1.0)))
      let nodeIn = SpikeSourcePoisson 0.1 0 0
      let nodeOut = Population 2 if_cond_exp "p1" 1
      let edge = DenseProjection effect nodeIn nodeOut
      let expected = "layer0 = v.Dense(p0, p1, weights=1.0, biases=0.0)"
      exec (pyNNEdge edge) ^. layers `shouldBe` (Map.singleton "layer0" expected)
    it "can translate a replicate edge to Python" $ do
      let effect = Static Excitatory (AllToAll (BiasGenerator (Constant 0)) (WeightGenerator (Constant 1.0)))
      let nodeIn = SpikeSourcePoisson 0.1 0 0
      let nodeOut1 = Population 2 if_cond_exp "p1" 1
      let nodeOut2 = Population 2 if_cond_exp "p2" 2
      let edge = ReplicateProjection effect nodeIn (nodeOut1, nodeOut2)
      let expected = "layer0 = v.Replicate(p0, (p1, p2), weights=(1.0, 1.0), biases=0.0)"
      exec (pyNNEdge edge) ^. layers `shouldBe` (Map.singleton "layer0" expected)
    it "can translate a merge edge to Python" $ do
      let effect = Static Excitatory (AllToAll (BiasGenerator (Constant 0)) (WeightGenerator (Constant 1.0)))
      let nodeIn1 = Population 2 if_cond_exp "p0" 0
      let nodeIn2 = Population 2 if_cond_exp "p1" 1
      let nodeOut = Population 4 if_cond_exp "p2" 2
      let edge = MergeProjection effect (nodeIn1, nodeIn2) nodeOut
      let expected = "layer0 = v.Merge((p0, p1), p2)"
      exec (pyNNEdge edge) ^. layers `shouldBe` (Map.singleton "layer0" expected)
    it "can connect edges to learning nodes with constant weights" $ do
      let effect = Static Excitatory (AllToAll (BiasGenerator (Constant 4)) (WeightGenerator (Constant 12.0)))
      let nodeIn = SpikeSourcePoisson 0.1 0 0
      let nodeOut = Population 2 if_cond_exp "p1" 1
      let edge = DenseProjection effect nodeIn nodeOut
      let expected = "layer0 = v.Dense(p0, p1, weights=12.0, biases=4.0)"
      exec (pyNNEdge edge) ^. layers `shouldBe` Map.singleton "layer0" expected
    it "can correctly translate two nodes into a dictionary" $ do
      let node1 = Population 1 if_cond_exp "p0" 0
      let node2 = Population 2 if_cond_exp "p1" 1
      let network = Network 0 [node1] [node2] [] []
      let code1 = "p0 = " ++ (fromRight "" (pyNNPopulationString if_cond_exp 1 0))
      let code2 = "p1 = " ++ (fromRight "" (pyNNPopulationString if_cond_exp 2 1))
      let expected = Map.fromList [("p0", code1), ("p1", code2)]
      exec (translate' network) ^. populations `shouldBe` expected
    it "can create a model with two inputs and two outputs" $ do
      let effect1 = Static Excitatory (AllToAll (BiasGenerator (Constant 0)) (WeightGenerator (GaussianRandom 1 2)))
      let effect2 = Static Excitatory (AllToAll (BiasGenerator (Constant 0)) (WeightGenerator (Constant 1.0)))
      let nodeIn1 = Population 1 if_cond_exp "p0" 0
      let nodeIn2 = Population 2 if_cond_exp "p1" 1
      let nodeOut1 = Population 3 if_cond_exp "p2" 2
      let nodeOut2 = Population 4 if_cond_exp "p3" 3
      let edge1 = DenseProjection effect1 nodeIn1 nodeOut1
      let edge2 = DenseProjection effect2 nodeIn2 nodeOut2
      let expectedDeclarations = ["l_decode = v.Decode(p3)", "model = v.Model(layer0, layer1, l_decode)"]
      let network = Network 0 [nodeIn1, nodeIn2] [] [edge1, edge2] [nodeOut1, nodeOut2]
      let model = exec (translate' network)
      model ^. declarations `shouldBe` expectedDeclarations
    it "can translate SNN nodes and edges to a PyNN model" $ do
      let input = Population 2 if_cond_exp "p1" 1
      let hidden = Population 4 if_cond_exp "p2" 2
      let output = Population 3 if_cond_exp "p3" 3
      let effect = Static Excitatory (AllToAll (BiasGenerator (Constant 0)) (WeightGenerator (Constant 1)))
      let edges = [ DenseProjection effect input hidden, DenseProjection effect hidden output ]
      let network = Network 0 [input] [hidden] edges [output]
      let expectedDeclarations = ["l_decode = v.Decode(p3)", "model = v.Model(layer0, layer1, l_decode)"]
      let model = exec (translate' network)
      model ^. declarations `shouldBe` expectedDeclarations
      Map.size (model ^. populations) `shouldBe` 3
      Map.keys (model ^. populations) `shouldBe` ["p1", "p2", "p3"]
      Map.size (model ^. layers) `shouldBe` 2
    it "can translate SNN nodes and edges with given weights and biases to a PyNN model" $ do
      let input = Population 2 if_cond_exp "p1" 1
      let output = Population 3 if_cond_exp "p2" 2
      let biases = Biases [0..2]
      let weights=Weights ((3><2) [1..])
      let effect = Static Excitatory (AllToAll biases weights)
      let edges = [ DenseProjection effect input output ]
      let network = Network 0 [input] [] edges [output]
      let expectedDeclarations = ["l_decode = v.Decode(p2)", "model = v.Model(layer0, l_decode)"]
      let expectedLayer = "layer0 = v.Dense(p1, p2, weights=np.array([[1.0,2.0],[3.0,4.0],[5.0,6.0]]), biases=np.array([0.0,1.0,2.0]))"
      let model = exec (translate' network)
      model ^. declarations `shouldBe` expectedDeclarations
      Map.size (model ^. populations) `shouldBe` 2
      model ^. layers `shouldBe` Map.singleton "layer0" expectedLayer
    it "can fail to translate SNN nodes and edges with wrong weight/bias config" $ do
      let input = Population 2 if_cond_exp "p1" 1
      let output = Population 3 if_cond_exp "p2" 2
      let biases = Biases [0..2]
      let weights=Weights ((2><3) [1..])
      let effect = Static Excitatory (AllToAll biases weights)
      let edges = [ DenseProjection effect input output ]
      let network = Network 0 [input] [] edges [output]
      (isLeft $ eval (translate' network)) `shouldBe` True

    it "can translate an entire SNN model to a Python script" $ do
      let input = Population 2 if_cond_exp "p0" 0
      let hidden = Population 4 if_cond_exp "p1" 1
      let output = Population 3 if_cond_exp "p2" 2
      let dictWithType = unpack $ encode if_cond_exp
      let dict = Regex.subRegex typeRegex dictWithType ""
      let effect1 = Static Excitatory (AllToAll (BiasGenerator (Constant 0)) (WeightGenerator (Constant 1)))
      let effect2 = Static Excitatory (AllToAll (Biases [0.0, 1.1, 2.2]) (WeightGenerator (GaussianRandom 2 1)))
      let edges = [ DenseProjection effect1 input hidden, DenseProjection effect2 hidden output ]
      let network = Network 0 [input] [hidden] edges [output]
      let task = Task (Nest 0 0) network 50
      let preample = PyNNPreample "# some config"
      let code = [i|import numpy as np
import volrpynn.nest as v
import pyNN.nest as pynn

# some config

p0 = pynn.Population(2, pynn.IF_cond_exp(**#{dict}))
p1 = pynn.Population(4, pynn.IF_cond_exp(**#{dict}))
p2 = pynn.Population(3, pynn.IF_cond_exp(**#{dict}))
layer0 = v.Dense(p0, p1, weights=1.0, biases=0.0)
layer1 = v.Dense(p1, p2, weights=np.random.normal(2.0, 1.0, (4, 3)), biases=np.array([0.0,1.1,2.2]))
l_decode = v.Decode(p2)
model = v.Model(layer0, layer1, l_decode)

optimiser = v.GradientDescentOptimiser(0.1, simulation_time=50.0)
if __name__ == "__main__":
    v.Main(model).train(optimiser)
|]
      translate task preample `shouldBe` Right code
    it "can translate an parallel SNN model to a Python script" $ do
      let input = Population 100 if_cond_exp "p0" 0
      let h1 = Population 20 if_cond_exp "p1" 1
      let h2 = Population 20 if_cond_exp "p2" 2
      let h3 = Population 10 if_cond_exp "p3" 3
      let h4 = Population 20 if_cond_exp "p4" 4
      let h5 = Population 10 if_cond_exp "p5" 5
      let h6 = Population 20 if_cond_exp "p6" 6
      let output = Population 10 if_cond_exp "p7" 7
      let dictWithType = unpack $ encode if_cond_exp
      let dict = Regex.subRegex typeRegex dictWithType ""
      let effect = Static Excitatory (AllToAll (BiasGenerator (Constant 0)) (WeightGenerator (Constant 1)))
      let edges = [ DenseProjection effect input h1
		  , ReplicateProjection effect h1 (h2, h4)
		  , DenseProjection effect h2 h3
		  , DenseProjection effect h4 h5
		  , MergeProjection effect (h3, h5) h6
		  , DenseProjection effect h6 output
		  ]
      let network = Network 0 [input] [h1, h2, h3, h4, h5, h6] edges [output]
      let task = Task (Nest 0 0) network 50
      let preample = PyNNPreample "# some config"
      let code = [i|import numpy as np
import volrpynn.nest as v
import pyNN.nest as pynn

# some config

p0 = pynn.Population(100, pynn.IF_cond_exp(**#{dict}))
p1 = pynn.Population(20, pynn.IF_cond_exp(**#{dict}))
p2 = pynn.Population(20, pynn.IF_cond_exp(**#{dict}))
p3 = pynn.Population(10, pynn.IF_cond_exp(**#{dict}))
p4 = pynn.Population(20, pynn.IF_cond_exp(**#{dict}))
p5 = pynn.Population(10, pynn.IF_cond_exp(**#{dict}))
p6 = pynn.Population(20, pynn.IF_cond_exp(**#{dict}))
p7 = pynn.Population(10, pynn.IF_cond_exp(**#{dict}))
layer0 = v.Dense(p0, p1, weights=1.0, biases=0.0)
layer1 = v.Replicate(p1, (p2, p4), weights=(1.0, 1.0), biases=0.0)
layer2 = v.Dense(p2, p3, weights=1.0, biases=0.0)
layer3 = v.Dense(p4, p5, weights=1.0, biases=0.0)
layer4 = v.Merge((p3, p5), p6)
layer5 = v.Dense(p6, p7, weights=1.0, biases=0.0)
l_decode = v.Decode(p7)
model = v.Model(layer0, layer1, layer2, layer3, layer4, layer5, l_decode)

optimiser = v.GradientDescentOptimiser(0.1, simulation_time=50.0)
if __name__ == "__main__":
    v.Main(model).train(optimiser)
|]
      translate task preample `shouldBe` Right code


text:
{-# OPTIONS_GHC -fno-warn-unused-do-bind #-}

module Scheme48.Types (
  LispVal(..),
  LispError(..),
  PrimFunc(..),
  Env,
  ThrowsError,
  IOThrowsError,
  IOPrimFunc(..),
  unwordsList,
  makeFunc,
  makeNormalFunc,
  makeVarArgs) where

import System.IO (Handle)
import Numeric ()
import Data.IORef
import Data.Ratio ()
import Data.Complex
import Data.Array
import Control.Monad.Except

import Text.ParserCombinators.Parsec

data IOPrimFunc = IOPrimFunc ([LispVal] -> IOThrowsError LispVal)

instance Eq IOPrimFunc where
   _ ==  _ = False
   _ /=  _ = True

data PrimFunc = PrimFunc ([LispVal] -> ThrowsError LispVal)

instance Eq PrimFunc where
   _ ==  _ = False
   _ /=  _ = True

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

instance Show LispVal where
  show = showVal

showVal :: LispVal -> String
showVal (Atom name)  = name
showVal (String c)   = "\"" ++ c ++ "\""
showVal (Number c)   = show c
showVal (Bool True)  = "#t"
showVal (Bool False) = "#f"
showVal (List c) = "(" ++ unwordsList c ++ ")"
showVal (DottedList x xs) = "(" ++ unwordsList x ++ " . " ++ showVal xs ++ ")"
showVal (Character c) = show c
showVal (Float n) = show n
showVal (Ratio r) = show r
showVal (Complex c) = (show $ realPart c) ++ " + " ++ (show $ imagPart c) ++ "i"
showVal (Vector v) = "#(" ++ unwordsList (elems v) ++ ")"
showVal (PrimitiveFunc _) = "<primitive>"
showVal (Port _)   = "<IO port>"
showVal (IOFunc _) = "<IO primitive>"
showVal (Func {params = args, vararg = varargs, body = _, closure = _}) =
   "(lambda (" ++ unwords (map show args) ++
      (case varargs of
         Nothing -> ""
         Just arg -> " . " ++ arg) ++ ") ...)"

makeFunc :: Maybe String -> Env -> [LispVal] -> [LispVal] -> IOThrowsError LispVal
makeFunc varargs' env' params' body' = return $ Func (map showVal params') varargs' body' env'

makeNormalFunc :: Env -> [LispVal] -> [LispVal] -> IOThrowsError LispVal
makeNormalFunc = makeFunc Nothing

makeVarArgs :: LispVal -> Env -> [LispVal] -> [LispVal] -> IOThrowsError LispVal
makeVarArgs = makeFunc . Just . showVal

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

-- Errors --

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

instance Show LispError where
  show = showError

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

type ThrowsError = Either LispError

-- Env --

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

type IOThrowsError = ExceptT LispError IO

ext: hs


text:
module Main where

import           Data.Either
import           DistanceGeometry
import           Generators
import           Numeric.LinearAlgebra
import           Numeric.LinearAlgebra.Data
import           ReadWrite
import           System.Directory

idir = "test/molecules"
odir = "test/results"

main :: IO ()
main = do
  putStrLn ""
  putStrLn "generateOneMolecule: Begin"
  generateOneMolecule idir odir "example" 16
  putStrLn "generateOneMolecule: End"
  putStrLn ""
  putStrLn "generateManyMolecules: Begin"
  generateManyMolecules idir odir "methan" 16 100
  putStrLn "generateManyMolecules: End"

max_stars_repo_name: wurthel/distance-geometry

max_issues_repo_name: wurthel/distance-geometry

max_forks_repo_name: wurthel/distance-geometry


text:
{-# LANGUAGE DeriveDataTypeable, DeriveGeneric #-}
-- |
-- Module    : Statistics.Distribution.Uniform
-- Copyright : (c) 2011 Aleksey Khudyakov
-- License   : BSD3
--
-- Maintainer  : bos@serpentine.com
-- Stability   : experimental
-- Portability : portable
--
-- Variate distributed uniformly in the interval.
module Statistics.Distribution.Uniform
    (
      UniformDistribution
    -- * Constructors
    , uniformDistr
    -- ** Accessors
    , uniformA
    , uniformB
    ) where

import Data.Aeson (FromJSON, ToJSON)
import Data.Binary (Binary)
import Data.Data (Data, Typeable)
import GHC.Generics (Generic)
import qualified Statistics.Distribution as D
import qualified System.Random.MWC       as MWC
import Data.Binary (put, get)
import Control.Applicative ((<$>), (<*>))

-- | Uniform distribution from A to B
data UniformDistribution = UniformDistribution {
      uniformA :: {-# UNPACK #-} !Double -- ^ Low boundary of distribution
    , uniformB :: {-# UNPACK #-} !Double -- ^ Upper boundary of distribution
    } deriving (Eq, Read, Show, Typeable, Data, Generic)

instance FromJSON UniformDistribution
instance ToJSON UniformDistribution

instance Binary UniformDistribution where
    put (UniformDistribution x y) = put x >> put y
    get = UniformDistribution <$> get <*> get

-- | Create uniform distribution.
uniformDistr :: Double -> Double -> UniformDistribution
uniformDistr a b
  | b < a     = uniformDistr b a
  | a < b     = UniformDistribution a b
  | otherwise = error "Statistics.Distribution.Uniform.uniform: wrong parameters"
-- NOTE: failure is in default branch to guard againist NaNs.

instance D.Distribution UniformDistribution where
  cumulative (UniformDistribution a b) x
    | x < a     = 0
    | x > b     = 1
    | otherwise = (x - a) / (b - a)

instance D.ContDistr UniformDistribution where
  density (UniformDistribution a b) x
    | x < a     = 0
    | x > b     = 0
    | otherwise = 1 / (b - a)
  quantile (UniformDistribution a b) p
    | p >= 0 && p <= 1 = a + (b - a) * p
    | otherwise        =
      error $ "Statistics.Distribution.Uniform.quantile: p must be in [0,1] range. Got: "++show p

instance D.Mean UniformDistribution where
  mean (UniformDistribution a b) = 0.5 * (a + b)

instance D.Variance UniformDistribution where
  -- NOTE: 1/sqrt 12 is not constant folded (#4101) so it's written as
  --       numerical constant. (Also FIXME!)
  stdDev   (UniformDistribution a b) = 0.2886751345948129 * (b - a)
  variance (UniformDistribution a b) = d * d / 12 where d = b - a

instance D.MaybeMean UniformDistribution where
    maybeMean = Just . D.mean

instance D.MaybeVariance UniformDistribution where
    maybeStdDev   = Just . D.stdDev

instance D.Entropy UniformDistribution where
  entropy (UniformDistribution a b) = log (b - a)

instance D.MaybeEntropy UniformDistribution where
  maybeEntropy = Just . D.entropy

instance D.ContGen UniformDistribution where
    genContVar (UniformDistribution a b) gen = MWC.uniformR (a,b) gen

max_stars_repo_name: StefanHubner/statistics

max_issues_repo_name: StefanHubner/statistics

max_forks_repo_name: StefanHubner/statistics


text:
module State where

import           Protolude

import           Data.Complex (Complex (..))
import qualified Data.Vector  as V

import           Instruction  (Instruction)

data MachineState =
  MachineState
    { quantumState :: QuantumState
    , instructions :: [Instruction]
    , pc           :: Integer
    , running      :: Bool
    }

type QuantumState = V.Vector (Complex Double)

-- Construct a quantum state given number of qubits
makeQuantumState :: Int -> QuantumState
makeQuantumState len =
  V.generate
    (2 ^ len)
    (\i ->
       if i == 0
         then (1 :+ 0)
         else (0 :+ 0))

max_stars_repo_name: colescott/quantum-interpreter

max_issues_repo_name: colescott/quantum-interpreter

max_forks_repo_name: colescott/quantum-interpreter


text:
module Data.TimeSeriesSpec (spec) where

import Prelude
import Test.Hspec
import Test.QuickCheck
import Control.Arrow (first)
import Data.Time (UTCTime)
import Data.Time.Clock.POSIX (posixSecondsToUTCTime, utcTimeToPOSIXSeconds)
import qualified Statistics.Sample as S
import qualified Data.Vector as V

import qualified Data.TimeSeries as TS
import Data.TimeSeries.Time (seconds)

sampleSeries :: TS.TimeSeries Double
sampleSeries = TS.tsSeries [1..] [10.0, 1.2, 32.4, 0.6, 11.0]

spec :: Spec
spec = do

  describe "Properties" $ do

    it "toList" $ property $
        \ts -> (TS.series . TS.toList) ts == (ts :: TS.TimeSeries Double)

    it "fold over series" $ property $
        \ts -> (foldr (+) 0.0 ts) - ((sum . (map snd) . TS.toList) (ts :: TS.TimeSeries Double)) < 10**(-6)

    it "resampling the same" $ property $
        \ts -> TS.resample (posixSecondsToUTCTime 1) (seconds 1) ts == (ts :: TS.TimeSeries Double)

    it "groupBy the same" $ property $
        \ts -> TS.groupBy (seconds 1) sum ts == (ts :: TS.TimeSeries Double)

  describe "Selecting sub Series" $ do

    it "return data point value at given index" $ do
        let pos = posixSecondsToUTCTime 2
        TS.valueAt pos sampleSeries `shouldBe` Just 1.2

    it "return Nothing value if wrong index" $ do
        let idx = [10, 20, 30]
        let values = [10.0, 12.0, 32.4]
        let pos = posixSecondsToUTCTime 1234
        TS.valueAt pos (TS.tsSeries idx values) `shouldBe` Nothing

    it "return subset" $ do
        let xs = TS.tsSeries [1..5] [10.0, 1.2, 32.4, 0.65, 11.0]
        let start = utcFromSeconds 2
        let end = utcFromSeconds 3
        TS.size (TS.slice start end xs) `shouldBe` 2

    it "return range" $ do
        let xs = TS.tsSeries [0, 5..] [1.0, 1.2, 32.4, 0.65, 11.0]
        let start = utcFromSeconds 2
        let end = utcFromSeconds 130
        TS.size (TS.slice start end xs) `shouldBe` 4

    it "return empty subset" $ do
        let xs = TS.tsSeries [1..5] [10.0, 1.2, 32.4, 0.65, 11.0]
        let start = utcFromSeconds 12
        let end = utcFromSeconds 34
        TS.slice start end xs `shouldBe` TS.emptySeries

  describe "Basic operations" $ do

    it "return series size" $ do
        let xs = TS.tsSeries [1..] [10.0, 1.2, 32.4, 0.65, 11.0]
        TS.size xs `shouldBe` 5

    it "first element" $ do
        let xs = TS.tsSeries [1..] [10.0, 1.2, 32.4, 0.65, 11.0]
        fmap TS.dpValue (TS.firstElem xs) `shouldBe` Just 10.0

    it "last element" $ do
        let xs = TS.tsSeries [1..] [10.0, 1.2, 32.4, 0.65, 11.0]
        fmap TS.dpValue (TS.lastElem xs) `shouldBe` Just 11.0

    it "element by index" $ do
        let xs = TS.tsSeries [1..] [10.0, 1.2, 32.4, 0.65, 11.0]
        fmap TS.dpValue (TS.elemAt 2 xs) `shouldBe` Just 32.4

    it "map over series" $ do
        let xs = TS.tsSeries [1..] [10.0, 1.2, 32.4, 0.65, 11.0]
        let ys = fmap (+ 2) xs
        TS.values ys `shouldBe` [12.0, 3.2, 34.4, 2.65, 13.0]

    it "maximum value" $ do
        let xs = TS.tsSeries [1..] [10.0, 1.2, 32.4, 0.65, 11.0]
        maximum xs `shouldBe` 32.4

    it "minimum value" $ do
        let xs = TS.tsSeries [1..5] [10.0, 1.2, 32.4, 0.65, 11.0]
        minimum xs `shouldBe` 0.65

  describe "Rolling" $ do

    it "sum" $ do
        let xs = TS.tsSeries [1..] [1.0, 2.0, 3.0, 4.0, 5.0]
        TS.rolling (seconds 2) sum xs `shouldBe` TS.tsSeries [2..] [3.0, 5.0, 7.0, 9.0]

  describe "Smoothing" $ do

    it "mean" $ do
        let xs = TS.tsSeries [1..] [1.0, 2.0, 3.0, 4.0, 5.0]
        TS.rolling (seconds 3) (S.mean . V.fromList) xs `shouldBe` TS.tsSeries [3..] [2.0, 3.0, 4.0]

  describe "Resampling" $ do

    it "missing values" $ do
        let xs = TS.tsSeries [1, 2, 3, 5] [1.0, 2.0, 3.0, 5.0]
        let startTime = posixSecondsToUTCTime 1
        TS.resample startTime (seconds 1) xs `shouldBe` TS.tsSeries [1, 2, 3, 4, 5] [1.0, 2.0, 3.0, 4.0, 5.0]

    it "missing lots of values" $ do
        let xs = TS.tsSeries [1, 5] [1.0, 5.0]
        let startTime = posixSecondsToUTCTime 1
        TS.resample startTime (seconds 1) xs `shouldBe` TS.tsSeries [1, 2, 3, 4, 5] [1.0, 2.0, 3.0, 4.0, 5.0]

    it "middle values" $ do
        let xs = TS.tsSeries [1, 3, 5] [1.0, 3.0, 5.0]
        let startTime = posixSecondsToUTCTime 0
        TS.resample startTime (seconds 2) xs `shouldBe` TS.tsSeries [0, 2, 4] [1.0, 2.0, 4.0]

  describe "Group by" $ do

    it "sum" $ do
        let xs = TS.tsSeries [1..] [1.0, 2.0, 3.0, 4.0, 5.0, 6.0]
        TS.groupBy (seconds 2) sum xs `shouldBe` TS.tsSeries [1, 3, 5] [3.0, 7.0, 11.0]

  describe "zip" $ do

    it "the same index" $ do
        let xs = TS.tsSeries [1..] [1.0, 2.0, 3.0, 4.0, 5.0, 6.0]
        let ys = TS.tsSeries [1..] [2.0, 3.0, 4.0, 5.0, 6.0, 7.0]
        fmap (\(x,y) -> x+y) (TS.zip xs ys) `shouldBe` TS.tsSeries [1..] [3.0, 5.0, 7.0, 9.0, 11.0, 13.0]

-- Convert to UTC Time from timestamp
utcFromSeconds :: Integer -> UTCTime
utcFromSeconds = posixSecondsToUTCTime . fromIntegral

-- Generate Time Series
instance (Arbitrary a) => Arbitrary (TS.TimeSeries a) where
    arbitrary = do
        vs <- arbitrary
        return $ TS.tsSeries [1..] vs

max_stars_repo_name: klangner/timeseries

max_issues_repo_name: klangner/timeseries

max_forks_repo_name: klangner/time-series-lib


text:
{-# LANGUAGE RankNTypes #-}

module Main where

import System.Random.MWC.Probability (Gen)
import qualified System.Random.MWC.Probability as MWC
import Control.Monad.Trans.State.Strict (execStateT)

import Statistics.Distribution.Poisson
import Statistics.Distribution.LogNormal
import Statistics.Distribution

import Conduit

import Data.Model
import Numeric.MCMC

-- ~from https://hackage.haskell.org/package/declarative-0.2.1/docs/src/Numeric-MCMC.html#mcmc
-- A Markov chain driven by an arbitrary transition operator.
-- now using Conduit instead of Pipes
chain
  :: PrimMonad m
  => Transition m b
  -> b
  -> Gen (PrimState m)
  -> Producer m b
chain transition = loop where
  loop state prng = do
    next <- lift (MWC.sample (execStateT transition state) prng)
    yield next
    loop next prng

main :: IO ()
main = withSystemRandom . asGenIO $
        \g -> chain (metropolis 10) c g =$ takeC 99999 =$ mapC (filter (`notElem` "|\"") . show) $$ mapM_C putStrLn

    where c = Chain t (testLH initPred) initPred Nothing
          t = Target testLH Nothing
          testData = [1, 100, 25, 50]
          initPred = [1, 1, 1, 1]
          ln = logNormalDistr' 50 20
          testLH xs = if any (<= 0) xs
                         then log 0.0
                         else predHistLLH testData (fmap poisson xs) + (sum $ map (logDensity ln) xs)

max_stars_repo_name: cspollard/HMCMC

max_issues_repo_name: cspollard/HMCMC

max_forks_repo_name: cspollard/HMCMC


text:
{-# LANGUAGE BangPatterns #-}
module FourierSeriesOfPinwheels where

import           Control.Monad.Parallel          as MP
import           Data.Array.IArray               as IA
import           Data.Array.Repa                 as R
import           Data.Array.Repa.IO.Binary
import           Data.Array.Repa.Repr.ForeignPtr
import           Data.Complex
import           Data.List                       as L
import           Data.Vector.Storable            as VS
import           Data.Vector.Unboxed             as VU
import           FokkerPlanck.GreensFunction
import           FokkerPlanck.Histogram
import           Image.IO
import           Pinwheel.Base
import           Pinwheel.FourierSeries2D
import           STC
import           System.Directory
import           System.Environment
import           System.FilePath
import           Text.Printf
import           Utils.Parallel
import           Utils.Time
import FourierMethod.FourierSeries2D

main = do
  args@(gpuIDsStr:lenStr:deltaStr:inverseDeltaStr:numR2FreqsStr:periodStr:phiFreqsStr:rhoFreqsStr:thetaFreqsStr:scaleFreqsStr:coefFilePath:stdStr:batchSizePolarFreqsStr:batchSizeR2FreqsStr:batchSizeR2Str:numThreadStr:_) <-
    getArgs
  let gpuIDs = read gpuIDsStr :: [Int]
      len = read lenStr :: Double
      delta = read deltaStr :: Double
      inverseDelta = read inverseDeltaStr :: Double
      numR2Freqs = read numR2FreqsStr :: Int
      period = read periodStr :: Double
      phiFreq = read phiFreqsStr :: Int
      phiFreqs = L.map fromIntegral [-phiFreq .. phiFreq]
      rhoFreq = read rhoFreqsStr :: Int
      rhoFreqs = L.map fromIntegral [-rhoFreq .. rhoFreq]
      thetaFreq = read thetaFreqsStr :: Int
      thetaFreqs = L.map fromIntegral [-thetaFreq .. thetaFreq]
      scaleFreq = read scaleFreqsStr :: Int
      scaleFreqs = L.map fromIntegral [-scaleFreq .. scaleFreq]
      numThread = read numThreadStr :: Int
      folderPath = "output/test/FourierSeriesOfPinwheels"
      std = read stdStr :: Double
      batchSizePolarFreqs = read batchSizePolarFreqsStr :: Int
      batchSizeR2Freqs = read batchSizeR2FreqsStr :: Int
      batchSizeR2 = read batchSizeR2Str :: Int
  removePathForcibly (folderPath </> "Recons")
  createDirectoryIfMissing True folderPath
  let !numAngularFreqs = 2 * (phiFreq + thetaFreq) + 1
      !angularFreqsCenter = div numAngularFreqs 2
      !numRadialFreqs = 2 * (rhoFreq + scaleFreq) + 1
      !numPoints' = round $ len / delta :: Int
      !numPoints =
        if odd numPoints'
          then numPoints'
          else numPoints' + 1
  coefficients <-
    computeFourierCoefficients
      gpuIDs
      numR2Freqs
      numPoints
      period
      delta
      batchSizeR2Freqs
      batchSizePolarFreqs
      numAngularFreqs
      numRadialFreqs
  writeArrayToStorableFile coefFilePath coefficients
  coefficientsArr' <-
    readArrayFromStorableFile
      coefFilePath
      (Z :. numRadialFreqs :. numAngularFreqs :. numR2Freqs :. numR2Freqs) :: IO (R.Array F DIM4 (Complex Double))
  coefficientsArr <- applyGaussian 1 std coefficientsArr'
  let coefficientsVecs =
        L.map
          (\(!rf, !af) ->
             VU.convert . toUnboxed . computeS . R.slice coefficientsArr $
             (Z :. rf :. af :. All :. All))
          [ (rf, af)
          | rf <- [0 .. numRadialFreqs - 1]
          , af <- [0 .. numAngularFreqs - 1]
          ]
  createDirectoryIfMissing True (folderPath </> "Coefficients")
  MP.mapM_
    (\(rf, af) ->
       plotImageRepaComplex
         (folderPath </> "Coefficients" </>
          printf
            "Coef(%d,%d).png"
            (rf - div numRadialFreqs 2)
            (af - div numAngularFreqs 2)) .
       ImageRepa 8 .
       computeS . extend (Z :. (1 :: Int) :. All :. All) . R.slice coefficients $
       (Z :. rf :. af :. All :. All)) .
    L.take 10 $
    [ (rf, af)
    | rf <- [0 .. numRadialFreqs - 1]
    , af <- [0 .. numAngularFreqs - 1]
    ]
  recons <-
    computeFourierSeries
      gpuIDs
      numR2Freqs
      (round len)
      period
      inverseDelta
      batchSizePolarFreqs
      batchSizeR2
      numAngularFreqs
      numRadialFreqs
      coefficientsVecs
  createDirectoryIfMissing True (folderPath </> "Recons")
  MP.mapM_
    (\(!rf, !af) ->
       plotImageRepaComplex
         (folderPath </> "Recons" </>
          printf
            "Recon(%d,%d).png"
            (rf - div numRadialFreqs 2)
            (af - div numAngularFreqs 2)) .
       ImageRepa 8 .
       computeS . extend (Z :. (1 :: Int) :. All :. All) . R.slice recons $
       (Z :. rf :. af :. All :. All)) .
    L.take 10 $
    [ (rf, af)
    | rf <- [0 .. numRadialFreqs - 1]
    , af <- [0 .. numAngularFreqs - 1]
    ]

max_stars_repo_name: XinhuaZhang/Stochastic-Completion-Field

max_issues_repo_name: XinhuaZhang/Stochastic-Completion-Field

max_forks_repo_name: XinhuaZhang/Stochastic-Completion-Field


text:
module Task (
writeTask,
readTask
) where
import System.IO.Strict as Strict
import System.IO
import InputParser
import Statistics

writeTask :: String -> [String] -> [Line] -> IO ()
writeTask filename info lns =
    let
        statistics = collectStat lns
        errorLine (Line position (Error text))
            = [concat ["ERROR: ", recordSourceText position, " - ", text]]
        errorLine _ = []
    in do
        outFile <- openFile filename WriteMode
        _ <- mapM (hPutStrLn outFile) info
        _ <- sequence $ hPutStrLn outFile <$> statistics
        _ <- sequence $ hPutStrLn outFile <$> (lns >>= errorLine)
        hClose outFile

readTask :: String -> IO [String]
readTask path = do
        result <- lines <$> Strict.readFile path
        return $ length result `seq` result --read whole file to allow its overwriting

max_stars_repo_name: sbenstewart/macOS-brew-Fortress

max_issues_repo_name: sbenstewart/macOS-brew-Fortress

max_forks_repo_name: sbenstewart/macOS-brew-Fortress


text:
{-# LANGUAGE DataKinds                                #-}

module Main where
import Control.Monad.State                           (execState)
import Data.List.Split                               (splitOn)
import           Control.DeepSeq
import           Control.Exception.Base
import           Control.Monad
import           Control.Monad.IO.Class
import           Control.Monad.Trans.Maybe
import           Control.Monad.Trans.State
import           Data.Bitraversable
import           Data.IDX
import           Data.Maybe
import           Data.Traversable
import           Data.Tuple
import qualified Data.Vector                         as V
import qualified Data.Vector.Generic                 as VG
import qualified Data.Vector.Unboxed                 as VU
import           GHC.TypeLits
import qualified Numeric.LinearAlgebra               as HM  -- non-backprop hmatrix
import qualified Numeric.LinearAlgebra.Data          as HMD
import qualified Numeric.LinearAlgebra.Static        as HMS -- hmatrix with type-checked operations :D
import System.IO

import MNIST
import Neural
import Regression                                    (activate, forward, loss, runTrain)
import Titanic                                       (getAge, processFeatures, processTargets)

readCSV :: [Char] -> [[[Char]]]
readCSV file = [splitOn "," x | x <- (lines file)]

accuracy :: [Double] -> [Double] -> Int -> Double
accuracy targets predictions batchSize = (fromIntegral correct) / (fromIntegral batchSize) where
    correct = length $ filter (\idx -> targets !! idx == (activate predictions) !! idx) [0..batchSize-1]

main = do
    -- Logistic Regression
    --
    -- contents <- readFile "data/titanic/train.csv"
    -- let weights = matrix 1 [0.0 | _ <- [0..7]]
    -- let inputs = processFeatures $ tail (readCSV contents)
    -- let targets = processTargets $ tail (readCSV contents)

    -- -- subtract the header row
    -- let batchSize = (length $ toLists inputs) - 1
    -- let finalWeights = (execState (runTrain inputs targets batchSize) weights)
    -- let finalPredictions = forward inputs finalWeights

    -- print $ accuracy (concat $ toLists targets) (concat $ toLists finalPredictions) batchSize
    -- print $ loss targets finalPredictions batchSize

    -- MNIST Multilayer Perceptron
    --
    -- Just train <- loadMNIST "data/mnist/train-images-idx3-ubyte" "data/mnist/train-labels-idx1-ubyte"
    -- Just test  <- loadMNIST "data/mnist/t10k-images-idx3-ubyte"  "data/mnist/t10k-labels-idx1-ubyte"
    weightsInit <- HM.rand 728 300
    inputsInit <- HM.rand 1 728
    biasInit <- HM.rand 1 300
    yInit <- HM.rand 1 10
    let weights = HMS.matrix (concat $ HM.toLists weightsInit) :: HMS.L 728 300
    let inputs = HMS.matrix (concat $ HM.toLists inputsInit) :: HMS.L 1 728
    let bias = HMS.matrix (concat $ HM.toLists biasInit) :: HMS.L 1 300
    let y = HMS.vector (concat $ HM.toLists yInit) :: HMS.R 10
    print $ relu $ linear weights bias inputs

max_stars_repo_name: staturecrane/haskell-ml

max_issues_repo_name: staturecrane/haskell-ml

max_forks_repo_name: staturecrane/haskell-ml


text:
{-# LANGUAGE ForeignFunctionInterface #-}
module Grenade.Layers.Internal.Pooling (
    poolForward
  , poolBackward
  ) where

import qualified Data.Vector.Storable as U ( unsafeToForeignPtr0, unsafeFromForeignPtr0 )

import           Foreign ( mallocForeignPtrArray, withForeignPtr )
import           Foreign.Ptr ( Ptr )

import           Numeric.LinearAlgebra ( Matrix , flatten )
import qualified Numeric.LinearAlgebra.Devel as U

import           System.IO.Unsafe ( unsafePerformIO )

poolForward :: Int -> Int -> Int -> Int -> Int -> Int -> Int -> Matrix Double -> Matrix Double
poolForward channels height width kernelRows kernelColumns strideRows strideColumns dataIm =
  let vec             = flatten dataIm
      rowOut          = (height - kernelRows) `div` strideRows + 1
      colOut          = (width - kernelColumns) `div` strideColumns + 1
      numberOfPatches = rowOut * colOut
  in unsafePerformIO $ do
    outPtr <- mallocForeignPtrArray (numberOfPatches * channels)
    let (inPtr, _) = U.unsafeToForeignPtr0 vec

    withForeignPtr inPtr $ \inPtr' ->
      withForeignPtr outPtr $ \outPtr' ->
        pool_forwards_cpu inPtr' channels height width kernelRows kernelColumns strideRows strideColumns outPtr'

    let matVec = U.unsafeFromForeignPtr0 outPtr (numberOfPatches * channels)
    return $ U.matrixFromVector U.RowMajor (rowOut * channels) colOut matVec

foreign import ccall unsafe
    pool_forwards_cpu
      :: Ptr Double -> Int -> Int -> Int -> Int -> Int -> Int -> Int -> Ptr Double -> IO ()

poolBackward :: Int -> Int -> Int -> Int -> Int -> Int -> Int -> Matrix Double -> Matrix Double -> Matrix Double
poolBackward channels height width kernelRows kernelColumns strideRows strideColumns dataIm dataGrad =
  let vecIm     = flatten dataIm
      vecGrad   = flatten dataGrad
  in unsafePerformIO $ do
    outPtr <- mallocForeignPtrArray (height * width * channels)
    let (imPtr, _) = U.unsafeToForeignPtr0 vecIm
    let (gradPtr, _) = U.unsafeToForeignPtr0 vecGrad

    withForeignPtr imPtr $ \imPtr' ->
      withForeignPtr gradPtr $ \gradPtr' ->
        withForeignPtr outPtr $ \outPtr' ->
          pool_backwards_cpu imPtr' gradPtr' channels height width kernelRows kernelColumns strideRows strideColumns outPtr'

    let matVec = U.unsafeFromForeignPtr0 outPtr (height * width * channels)
    return $ U.matrixFromVector U.RowMajor (height * channels) width matVec

foreign import ccall unsafe
    pool_backwards_cpu
      :: Ptr Double -> Ptr Double -> Int -> Int -> Int -> Int -> Int -> Int -> Int -> Ptr Double -> IO ()

max_issues_repo_name: Alien-Inc/grenade

max_forks_repo_name: Alien-Inc/grenade


text:
{-# LANGUAGE BangPatterns #-}
module Util where
import Numeric.LinearAlgebra
import Numeric.LinearAlgebra.Util
import Data.List
import Data.Ord (comparing)
import Control.Applicative
import System.Random
import System.IO.Unsafe
import Data.IORef
import qualified Data.Map as Map

partitionList :: [Double] -> [a] -> [[a]]
partitionList ps xs = iterate (proportional ps) xs
  where
    sumP = sum ps
    lenX = length xs
    proportional = map (/sumP)
    iterate [] leftover = [leftover]
    iterate (p:rest) li = (take proportion li):iterate rest (drop proportion li)
      where
        proportion = floor $ p * (fromIntegral lenX)

{- shuffle the input and partition into training, test and validation sets in the
   given proportions
   -}
shuffleAndPartitionInput :: RandomGen g => g -> (Matrix Double,Matrix Double) ->
                                           (Double,Double,Double) ->
                                           ((Matrix Double,Matrix Double),
                                            (Matrix Double,Matrix Double),
                                            (Matrix Double,Matrix Double))
shuffleAndPartitionInput gen (inputX,inputY) (trainProp,valProp,testProp) = (train,val,test)
  where
   m = rows inputY
   (trainIndices:valIndices:testIndices:rest) = partitionList [trainProp,valProp,testProp] (shuffle gen m [0..])
   train   = get trainIndices
   val     = get valIndices
   test    = get testIndices
   get ind = (extractRows ind inputX,extractRows ind inputY)

-- |shuffle a list
shuffle :: RandomGen g => g -> Int -> [a] -> [a]
shuffle gen 0   _  = []
shuffle gen len xs =
        let
                n = fst $ randomR (0, len - 1) gen
                (y, ys) =  choose n xs
                ys' = shuffle gen (len - 1) ys
        in y:ys'

choose _ [] = error "choose: index out of range"
choose 0 (x:xs) = (x, xs)
choose i (x:xs) = let (y, ys) = choose (i - 1) xs in (y, x:ys)

-- | choose a random sample from x and y
randomSample :: RandomGen g => g -> (Matrix Double,Matrix Double) -> Int -> (Matrix Double, Matrix Double)
randomSample gen (x,y) size = (newX,newY)
  where
    m = rows y
    indices = sort $ uniqueRandoms size m gen
    newX = extractRows indices x
    newY = extractRows indices y

uniqueRandoms :: RandomGen g => Int -> Int -> g -> [Int]
uniqueRandoms size max gen | max < size = error "size greater than max"
                           | otherwise = take size $ shuffle gen max [0..(max-1)]

-- fatureNormalize(X)
-- Normalizes the features in X
-- FEATURENORMALIZE(X) returns a normalized version of X where
-- the mean value of each feature is 0 and the standard deviation
-- is 1. This is often a good preprocessing step to do when
-- working with learning algorithms.
featureNormalize :: Matrix Double -> (Matrix Double, Matrix Double, Matrix Double)
featureNormalize x = ((doNormalize x mu sigma), mu, sigma)
                      where
                        mu = fromRows $ [fromList $ map mean $ map toList $ toColumns x]
                        sigma = fromRows $ [fromList $ map stddev $ map toList $ toColumns x]

-- normalize an input matrix given sigma and mu
doNormalize :: Matrix Double -> Matrix Double -> Matrix Double -> Matrix Double
doNormalize x mu sigma = (x - mu) / sigma

-- split a list into words by a function
wordsBy :: (a -> Bool) -> [a] -> [[a]]
wordsBy p s = case dropWhile p s of
    []      -> []
    s':rest -> (s':w) : wordsBy p (drop 1 s'')
        where (w, s'') = break p rest

-- read a csv file of numbers into a matrix
readMatrixCSV :: FilePath -> IO (Matrix Double)
readMatrixCSV path = do
    s <- map (wordsBy (\ x -> x == ',' )) <$> lines <$> readFile path
    let d = map (map read) $ s :: [[Double]]
    let cols = length $ head d
    return $ reshape cols $ fromList (concat d)

-- sum up the rows in a matrix, returning a column vector of sums
sumRows m = fromColumns $ [ sum $ toColumns m ]

-- sum up the columns in a matrix, returning a row vector of sums
sumColumns m = fromRows $ [ sum $ toRows m ]

-- pad a matrix on the left with a value
padLeft :: Matrix Double -> Double -> Matrix Double
padLeft !m num = fromBlocks [[scalar num, m]]

-- convenience function similar to above
addBias !m = fromBlocks [[1,m]]

-- |Numerically stable mean
mean :: Floating a => [a] -> a
mean x = fst $ foldl' (\(!m, !n) x -> (m+(x-m)/(n+1),n+1)) (0,0) x

-- |Same as 'mean'
average :: Floating a => [a] -> a
average = mean

-- |Median
median :: (Floating a, Ord a) => [a] -> a
median x | odd n  = head  $ drop (n `div` 2) x'
         | even n = mean $ take 2 $ drop i x'
                  where i = (length x' `div` 2) - 1
                        x' = sort x
                        n  = length x
-- |Range
range :: (Num a, Ord a) => [a] -> a
range xs = maximum xs - minimum xs

-- |Standard deviation of sample
stddev :: (Floating a) => [a] -> a
stddev xs = sqrt $ var xs

-- |Sample variance
var xs = (var' 0 0 0 xs) / (fromIntegral $ length xs - 1)
    where
      var' _ _ s [] = s
      var' m n s (x:xs) = var' nm (n + 1) (s + delta * (x - nm)) xs
         where
           delta = x - m
           nm = m + delta/(fromIntegral $ n + 1)

-- |Replace items in a list
replace :: Eq a => a -> a -> [a] -> [a]
replace _ _ [] = []
replace target replacement (x:xs) | target == x = replacement:(replace target replacement xs)
                                  | otherwise   = x:(replace target replacement xs)

-- The Sigmoid function
sigmoid x = 1.0 / (1 + exp (-x))

sigmoidGradient z = (sigmoid z) * (1 - (sigmoid z))

-- |Find the maximum value of each row and return a column vector of their indices
rowMaxIndex :: Matrix Double -> Matrix Double
rowMaxIndex m = asColumn $ fromList $ map maxIndex $ toLists m
  where
    maxIndex xs = snd $ maximum $ zip xs [0,1..]

-- |Turn a list of matrices into a vector
flattenConcat ms = fromList $ concat $ map (toList . flatten) ms

-- |Memoizes a function with unsafePerformIO
{-# NOINLINE unsafeMemoize #-}
unsafeMemoize :: Ord a => (a -> b) -> (a -> b)
unsafeMemoize f = unsafePerformIO $ do
    r <- newIORef Map.empty
    return $ \ x -> unsafePerformIO $ do
        m <- readIORef r
        case Map.lookup x m of
            Just y  -> return y
            Nothing -> do
                    let y = f x
                    writeIORef r (Map.insert x y m)
                    return y

{- Generates the train and cross validation set errors needed
  to plot a learning curve
   [error_train, error_val] = ...
       LEARNINGCURVE(X, y, Xval, yval, lambda) returns the train and
       cross validation set errors for a learning curve. In particular,
       it returns two vectors of the same length - error_train and
       error_val. Then, error_train(i) contains the training error for
       i examples (and similarly for error_val(i)).

   In this function, you will compute the train and test errors for
   dataset sizes from 1 up to m. In practice, when working with larger
   datasets, you might want to do this in larger intervals.
-}

learningCurves :: (Matrix Double -> Matrix Double -> Matrix Double) ->
                  (Matrix Double -> Matrix Double -> Matrix Double -> Double) ->
                   Matrix Double ->
                   Matrix Double ->
                   Matrix Double ->
                   Matrix Double ->
                   Int ->
                   [(Double,Double,Int)]
learningCurves train cost x y xval yval step_size = iterate 1
  where
    m = rows y
    iterate i | i >= m    = []
              | otherwise = (train_cost,cval_cost,i):iterate (i+step_size)
      where
        x' = takeRows i x
        y' = takeRows i y
        theta = train x' y'
        train_cost = cost x' y' theta
        cval_cost = cost xval yval theta

{- POLYFEATURES Maps X (1D vector) into the p-th power
   [X_poly] = POLYFEATURES(X, p) takes a data matrix X (size m x 1) and
   maps each example into its polynomial features where
   X_poly(i, :) = [X(i) X(i).^2 X(i).^3 ...  X(i).^p];
-}
oneDPolyFeatures :: Matrix Double -> Int -> Matrix Double
oneDPolyFeatures x p = fromColumns $ iterate 1
  where
    xVect = (toColumns x) !! 0
    iterate i | i > p     = []
              | otherwise = (xVect ^ i):iterate (i+1)

{- VALIDATIONCURVE Generate the train and validation errors needed to
   plot a validation curve that we can use to select lambda
   [lambda_vec, error_train, error_val] = ...
   VALIDATIONCURVE(X, y, Xval, yval) returns the train
   and validation errors (in error_train, error_val)
   for different values of lambda. You are given the training set (X,
   y) and validation set (Xval, yval).
-}
validationCurve :: (Double -> Matrix Double) -> (Matrix Double -> Double) -> (Matrix Double -> Double) -> [Double] -> [(Double,Double)]
validationCurve trainF trainErrorF valErrorF candidates = iterate candidates
  where
    iterate [] = []
    iterate (c:cs) = (trainEr,valEr):(iterate cs)
      where
        theta = trainF c
        trainEr = trainErrorF theta
        valEr = valErrorF theta

{- mapOutputsVectorToTargetMatrix
   Takes an m-dimension column vector of outputs in space {0:(k-1)} and maps it to an m><k target
   matrix in space {0,1} where each row represents an output and a 1 is used to indicate which
   class (column) an output is in.  Eg if there are 3 output classes and the input vector looks
   like this:

   [ 0,
     0,
     1,
     2 ]

  The output matrix will look like this:

  [ 1, 0, 0;
    1, 0, 0;
    0, 1, 0;
    0, 0, 1 ]
-}
mapOutputsVectorToTargetMatrix y = (targetYsMatrix,targets)
  where
    targets = nub $ toList $ flatten y
    targetYsList = map (\t -> cond y t 0 1 0) (map scalar targets)
    targetYsMatrix = fromBlocks [targetYsList]

max_stars_repo_name: thade/haskellml

max_issues_repo_name: thade/haskellml

max_forks_repo_name: thade/haskellml


text:
{-# LANGUAGE TypeOperators #-}
{-# LANGUAGE FlexibleContexts #-}
{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE ExistentialQuantification #-}
{-# LANGUAGE DefaultSignatures #-}
{-# LANGUAGE OverlappingInstances #-}
{-# LANGUAGE ScopedTypeVariables #-}

module Test.SmartCheck.Types
  ( SubT(..)
  , subT
  , Result(..)
  , SubTypes(..)
  , Idx(..)
  , Subst(..)
  , Replace(..)
  , errorMsg
  -- ** For constructing new instances of `SubTypes`
  , gst
  , grc
  , gtc
  , gsf
--  , gsz
  ) where

import GHC.Generics
import Data.Tree
import Data.Typeable
import Control.Applicative
import Control.Monad (ap)

-- For instances
import Data.Word
import Data.Int
import Data.Ratio
import Data.Complex

import qualified Test.QuickCheck as Q

-------------------------------------------------------------------------------

-- | Nominally, a list for value generalization indexes and existential
-- generalization.
data Replace a = Replace { unVals :: [a], unConstrs :: [a] }
  deriving (Show, Read, Eq)

--------------------------------------------------------------------------------
-- Result type
--------------------------------------------------------------------------------

-- | Possible results of iterateArb.
data Result a =
    BaseType      -- ^ Base type. Won't analyze.
  | FailedPreCond -- ^ Couldn't satisfy the precondition of a QuickCheck
                  --   property
  | FailedProp    -- ^ Failed the property---either we expect failure and it
                  --   passes or we expect to pass it and we fail.
  | Result a      -- ^ Satisfied it, with the satisfying value.
  deriving (Show, Read, Eq)

instance Functor Result where
  fmap _ BaseType      = BaseType
  fmap _ FailedPreCond = FailedPreCond
  fmap _ FailedProp    = FailedProp
  fmap f (Result a)    = Result (f a)

instance Monad Result where
  return a            = Result a
  BaseType      >>= _ = BaseType
  FailedPreCond >>= _ = FailedPreCond
  FailedProp    >>= _ = FailedProp
  Result a      >>= f = f a

instance Applicative Result where
 pure  = return
 (<*>) = ap

-------------------------------------------------------------------------------
-- Indexing
-------------------------------------------------------------------------------

-- | Index into a Tree/Forest, where level is the depth from the root and column
-- is the distance d is the dth value on the same level.  Thus, all left-most
-- nodes are in column 0.  This is a "matrix view" of tree-structured data.
data Idx = Idx { level :: Int, column :: Int }
  deriving Eq

instance Show Idx where
  show (Idx l c) = foldr1 (++) ["Idx ", show l, " ", show c]

-- | Keep or substitue a value in the tree.
data Subst = Keep | Subst
  deriving (Show, Eq, Read)

-- | Sort in order of depth first then left to right.
instance Ord Idx where
  compare (Idx l0 c0) (Idx l1 c1) | l0 < l1 = LT
                                  | l0 > l1 = GT
                                  | c0 < c1 = LT
                                  | c0 > c1 = GT
                                  | True    = EQ

-------------------------------------------------------------------------------
-- User-defined subtypes of data
-------------------------------------------------------------------------------

data SubT = forall a. (Q.Arbitrary a, SubTypes a)
          => SubT { unSubT :: a }

subT :: (Q.Arbitrary a, SubTypes a) => a -> SubT
subT = SubT

-- Would require SubT to derive Eq.
-- instance Eq SubT where
--   SubT a == SubT b = cast a == Just b

instance Show SubT where
  show (SubT t) = show t

-- | This class covers algebraic datatypes that can be transformed into Trees.
-- subTypes is the main method, placing values into trees.
--
-- for a datatype with constructors A and C,
--
-- > subTypes (A (C 0) 1)
-- > [Node {rootLabel = C 0, subForest = []}]
--
class (Q.Arbitrary a, Show a, Typeable a) => SubTypes a where
  -----------------------------------------------------------
  -- | Turns algebraic data into a forest representation.
  subTypes :: a -> Forest SubT
  default subTypes :: (Generic a, GST (Rep a))
                   => a -> Forest SubT
  subTypes = gst . from
  -----------------------------------------------------------
  -- | Base types (e.g., Int, Char) aren't analyzed.
  baseType :: a -> Bool
  baseType _ = False
  -----------------------------------------------------------
  -- | Generically replace child i in m with value s.  A total function: returns
  -- Nothing if you try to replace a child with an ill-typed child s.  (Returns
  -- Just (the original data) if your index is out of bounds).
  replaceChild :: Typeable b => a -> Forest Subst -> b -> Maybe a
  default replaceChild :: (Generic a, GST (Rep a), Typeable b)
                       => a -> Forest Subst -> b -> Maybe a
  replaceChild a forest b = fmap to $ grc (from a) forest b
  -----------------------------------------------------------
  -- | Get the string representation of the constructor.
  toConstr :: a -> String
  default toConstr :: (Generic a, GST (Rep a)) => a -> String
  toConstr = gtc . from
  -----------------------------------------------------------
  -- | showForest generically shows a value while preserving its structure (in a
  -- Tree).  Always returns either a singleton list containing the tree (a
  -- degenerate forest) or an empty list for baseTypes.  An invariant is that
  -- the shape of the tree produced by showForest is the same as the one
  -- produced by subTypes.
  showForest :: a -> Forest String
  default showForest :: (Generic a, GST (Rep a))
                     => a -> Forest String
  showForest = gsf . from
  -----------------------------------------------------------

-------------------------------------------------------------------------------
-- Generic representation
-------------------------------------------------------------------------------

class GST f where
  -- Names are abbreviations of the corresponding method names above.
  gst :: f a -> Forest SubT
  grc :: Typeable b => f a -> Forest Subst -> b -> Maybe (f a)
  gtc :: f a -> String
  gsf :: f a -> Forest String
  gsz :: f a -> Int

instance GST U1 where
  gst U1 = []
  grc _ _ _ = Nothing
  gtc U1 = ""
  gsf U1 = []
  gsz U1 = 0

instance (GST a, GST b) => GST (a :*: b) where
  gst (a :*: b) = gst a ++ gst b

  grc (a :*: b) forest c =
    case forest of
      []    -> Just (a :*: b)
      ls    -> do let (x,y) = splitAt (gsz a) ls
                  left  <- grc a x c
                  right <- grc b y c
                  return $ left :*: right

  gtc (a :*: b) = gtc a ++ gtc b
  gsf (a :*: b) = gsf a ++ gsf b
  gsz (a :*: b) = gsz a + gsz b

instance (GST a, GST b) => GST (a :+: b) where
  gst (L1 a) = gst a
  gst (R1 b) = gst b

  grc (L1 a) forest c = grc a forest c >>= return . L1
  grc (R1 a) forest c = grc a forest c >>= return . R1

  gtc (L1 a) = gtc a
  gtc (R1 a) = gtc a

  gsf (L1 a) = gsf a
  gsf (R1 a) = gsf a

  gsz (L1 a) = gsz a
  gsz (R1 a) = gsz a

-- Constructor meta-information
instance (Constructor c, GST a) => GST (M1 C c a) where
  gst (M1 a) = gst a
  grc (M1 a) forest c = grc a forest c >>= return . M1
  gtc = conName

  gsf m@(M1 a) = [ Node (conName m) (gsf a) ]

  gsz (M1 a) = gsz a

-- All the other meta-information (selector, module, etc.)
instance GST a => GST (M1 i k a) where
  gst (M1 a) = gst a
  grc (M1 a) forest c = grc a forest c >>= return . M1
  gtc (M1 a) = gtc a
  gsf (M1 a) = gsf a
  gsz (M1 a) = gsz a

instance (Show a, Q.Arbitrary a, SubTypes a, Typeable a) => GST (K1 i a) where
  gst (K1 a) = if baseType a
                 then [ Node (subT a) [] ]
                 else [ Node (subT a) (subTypes a) ]

  grc (K1 a) forest c =
    case forest of
      []                  -> Just (K1 a)
      (Node Keep  _  : _) -> Just (K1 a)
      (Node Subst [] : _) -> fmap K1 (cast c)
      (Node Subst ls : _) -> replaceChild a ls c >>= return . K1

  gtc _ = ""

  gsf (K1 a) = if baseType a then [ Node (show a) [] ] else showForest a

  gsz _ = 1

-------------------------------------------------------------------------------
-- We cover the instances supported by QuickCheck:
-- http://hackage.haskell.org/packages/archive/QuickCheck/2.4.2/doc/html/Test-QuickCheck-Arbitrary.html

instance SubTypes Bool    where baseType _    = True
instance SubTypes Char    where baseType _    = True
instance SubTypes Double  where baseType _    = True
instance SubTypes Float   where baseType _    = True
instance SubTypes Int     where baseType _    = True
instance SubTypes Integer where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance SubTypes Int8    where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance SubTypes Int16   where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance SubTypes Int32   where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance SubTypes Int64   where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance SubTypes Word    where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance SubTypes Word8   where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance SubTypes Word16  where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance SubTypes Word32  where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance SubTypes Word64  where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance SubTypes ()      where baseType _    = True
instance ( Q.Arbitrary a, SubTypes a, Typeable a) => SubTypes [a]
instance (Integral a, Q.Arbitrary a, SubTypes a, Typeable a)
  => SubTypes (Ratio a)   where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance (RealFloat a, Q.Arbitrary a, SubTypes a, Typeable a)
  => SubTypes (Complex a) where
  subTypes _    = []
  baseType _    = True
  replaceChild  = replaceChild'
  toConstr      = toConstr'
  showForest    = showForest'
instance (Q.Arbitrary a, SubTypes a, Typeable a) => SubTypes (Maybe a)
instance ( Q.Arbitrary a, SubTypes a, Typeable a
         , Q.Arbitrary b, SubTypes b, Typeable b)
         => SubTypes (Either a b)
instance ( Q.Arbitrary a, SubTypes a, Typeable a
         , Q.Arbitrary b, SubTypes b, Typeable b)
         => SubTypes (a, b)
instance ( Q.Arbitrary a, SubTypes a, Typeable a
         , Q.Arbitrary b, SubTypes b, Typeable b
         , Q.Arbitrary c, SubTypes c, Typeable c)
         => SubTypes (a, b, c)
instance ( Q.Arbitrary a, SubTypes a, Typeable a
         , Q.Arbitrary b, SubTypes b, Typeable b
         , Q.Arbitrary c, SubTypes c, Typeable c
         , Q.Arbitrary d, SubTypes d, Typeable d)
         => SubTypes (a, b, c, d)
instance ( Q.Arbitrary a, SubTypes a, Typeable a
         , Q.Arbitrary b, SubTypes b, Typeable b
         , Q.Arbitrary c, SubTypes c, Typeable c
         , Q.Arbitrary d, SubTypes d, Typeable d
         , Q.Arbitrary e, SubTypes e, Typeable e)
         => SubTypes (a, b, c, d, e)

-------------------------------------------------------------------------------
-- Helpers

-- These should never be directly called.  We provide compatible instances anyway.
toConstr' :: Show a => a -> String
toConstr' = show

replaceChild' :: (Typeable a, Typeable b)
              => a -> Forest Subst -> b -> Maybe a
replaceChild' a []                 _ = Just a
replaceChild' a (Node Keep  _ : _) _ = Just a
replaceChild' _ (Node Subst _ : _) b = cast b

showForest' :: Show a => a -> Forest String
showForest' _ = []

-------------------------------------------------------------------------------

errorMsg :: String -> a
errorMsg loc = error $ "SmartCheck error: unexpected error in " ++ loc
    ++ ".  Please file a bug report at "
    ++ "<https://github.com/leepike/SmartCheck/issues>."

-------------------------------------------------------------------------------


text:
{-# LANGUAGE CPP #-}
{-# LANGUAGE GADTs #-}
{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE FlexibleContexts #-}
{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE DefaultSignatures #-}
{-# LANGUAGE MultiParamTypeClasses #-}
{-# LANGUAGE UndecidableInstances #-}
{-# LANGUAGE FunctionalDependencies #-}

#ifndef MIN_VERSION_base
#define MIN_VERSION_base(x,y,z) 1
#endif
-----------------------------------------------------------------------------
-- |
-- Module      :  Control.Lens.Each
-- Copyright   :  (C) 2012-13 Edward Kmett
-- License     :  BSD-style (see the file LICENSE)
-- Maintainer  :  Edward Kmett <ekmett@gmail.com>
-- Stability   :  experimental
-- Portability :  non-portable
--
-----------------------------------------------------------------------------
module Control.Lens.Each
  (
  -- * Indices
    Index
  -- * Each
  , Each(..)
  ) where

import Control.Applicative
import Control.Lens.Cons as Lens
import Control.Lens.Internal.Deque
import Control.Lens.Internal.Setter
import Control.Lens.Indexed as Lens
import Control.Lens.Iso
import Control.Lens.Type
import Control.Lens.Traversal
import Data.Array.Unboxed as Unboxed
import Data.Array.IArray as IArray
import Data.ByteString as StrictB
import Data.ByteString.Lazy as LazyB
import Data.Complex
import Data.Foldable as Foldable
import Data.Functor.Identity
import Data.HashMap.Lazy as HashMap
import Data.HashSet
import Data.Int
import Data.IntMap as IntMap
import Data.IntSet
import Data.List.NonEmpty as NonEmpty
import Data.Map as Map
import Data.Set
import Data.Sequence as Seq
import Data.Text as StrictT
import Data.Text.Lazy as LazyT
import Data.Traversable
import Data.Tree as Tree
import qualified Data.Vector as Vector
import qualified Data.Vector.Primitive as Prim
import Data.Vector.Primitive (Prim)
import qualified Data.Vector.Storable as Storable
import Data.Vector.Storable (Storable)
import qualified Data.Vector.Unboxed as Unboxed
import Data.Vector.Unboxed (Unbox)
import Data.Word

-- | This is a common 'Index' type shared by 'Each', 'Control.Lens.At.At', 'Control.Lens.At.Contains' and 'Control.Lens.At.Ixed'.
type family Index (s :: *) :: *
type instance Index (e -> a) = e
type instance Index IntSet = Int
type instance Index (Set a) = a
type instance Index (HashSet a) = a
type instance Index [a] = Int
type instance Index (NonEmpty a) = Int
type instance Index (Seq a) = Int
type instance Index (a,b) = Int
type instance Index (a,b,c) = Int
type instance Index (a,b,c,d) = Int
type instance Index (a,b,c,d,e) = Int
type instance Index (a,b,c,d,e,f) = Int
type instance Index (a,b,c,d,e,f,g) = Int
type instance Index (a,b,c,d,e,f,g,h) = Int
type instance Index (a,b,c,d,e,f,g,h,i) = Int
type instance Index (IntMap a) = Int
type instance Index (Map k a) = k
type instance Index (HashMap k a) = k
type instance Index (Array i e) = i
type instance Index (UArray i e) = i
type instance Index (Vector.Vector a) = Int
type instance Index (Prim.Vector a) = Int
type instance Index (Storable.Vector a) = Int
type instance Index (Unboxed.Vector a) = Int
type instance Index (Complex a) = Int
type instance Index (Identity a) = ()
type instance Index (Maybe a) = ()
type instance Index (Tree a) = [Int]
type instance Index StrictT.Text = Int
type instance Index LazyT.Text = Int64
type instance Index StrictB.ByteString = Int
type instance Index LazyB.ByteString = Int64

-- $setup
-- >>> :set -XNoOverloadedStrings
-- >>> import Control.Lens
-- >>> import Data.Text.Strict.Lens as Text
-- >>> import Data.Char as Char

-- | Extract 'each' element of a (potentially monomorphic) container.
--
-- Notably, when applied to a tuple, this generalizes 'Control.Lens.Traversal.both' to arbitrary homogeneous tuples.
--
-- >>> (1,2,3) & each *~ 10
-- (10,20,30)
--
-- It can also be used on monomorphic containers like 'StrictT.Text' or 'StrictB.ByteString'.
--
-- >>> over each Char.toUpper ("hello"^.Text.packed)
-- "HELLO"
--
-- 'each' is an 'IndexedTraversal', so it can be used to access keys in many containers:
--
-- >>> itoListOf each $ Map.fromList [("hello",2),("world",4)]
-- [("hello",2),("world",4)]
--
-- >>> ("hello","world") & each.each %~ Char.toUpper
-- ("HELLO","WORLD")
class (Functor f, Index s ~ Index t) => Each f s t a b | s -> a, t -> b, s b -> t, t a -> s where
  each :: IndexedLensLike (Index s) f s t a b
#ifndef HLINT
  default each :: (Applicative f, Traversable g, s ~ g a, t ~ g b, Index s ~ Int, Index t ~ Int) => IndexedLensLike Int f s t a b
  each = traversed
  {-# INLINE each #-}
#endif

-- | @'each' :: 'IndexedTraversal' 'Int' (a,a) (b,b) a b@
instance (Applicative f, a~a', b~b') => Each f (a,a') (b,b') a b where
  each f ~(a,b) = (,) <$> f' (0 :: Int) a <*> f' 1 b
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' (a,a,a) (b,b,b) a b@
instance (Applicative f, a~a2, a~a3, b~b2, b~b3) => Each f (a,a2,a3) (b,b2,b3) a b where
  each f ~(a,b,c) = (,,) <$> f' (0 :: Int) a <*> f' 1 b <*> f' 2 c
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' (a,a,a,a) (b,b,b,b) a b@
instance (Applicative f, a~a2, a~a3, a~a4, b~b2, b~b3, b~b4) => Each f (a,a2,a3,a4) (b,b2,b3,b4) a b where
  each f ~(a,b,c,d) = (,,,) <$> f' (0 :: Int) a <*> f' 1 b <*> f' 2 c <*> f' 3 d
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' (a,a,a,a,a) (b,b,b,b,b) a b@
instance (Applicative f, a~a2, a~a3, a~a4, a~a5, b~b2, b~b3, b~b4, b~b5) => Each f (a,a2,a3,a4,a5) (b,b2,b3,b4,b5) a b where
  each f ~(a,b,c,d,e) = (,,,,) <$> f' (0 :: Int) a <*> f' 1 b <*> f' 2 c <*> f' 3 d <*> f' 4 e
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' (a,a,a,a,a,a) (b,b,b,b,b,b) a b@
instance (Applicative f, a~a2, a~a3, a~a4, a~a5, a~a6, b~b2, b~b3, b~b4, b~b5, b~b6) => Each f (a,a2,a3,a4,a5,a6) (b,b2,b3,b4,b5,b6) a b where
  each f ~(a,b,c,d,e,g) = (,,,,,) <$> f' (0 :: Int) a <*> f' 1 b <*> f' 2 c <*> f' 3 d <*> f' 4 e <*> f' 5 g
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' (a,a,a,a,a,a,a) (b,b,b,b,b,b,b) a b@
instance (Applicative f, a~a2, a~a3, a~a4, a~a5, a~a6, a~a7, b~b2, b~b3, b~b4, b~b5, b~b6, b~b7) => Each f (a,a2,a3,a4,a5,a6,a7) (b,b2,b3,b4,b5,b6,b7) a b where
  each f ~(a,b,c,d,e,g,h) = (,,,,,,) <$> f' (0 :: Int) a <*> f' 1 b <*> f' 2 c <*> f' 3 d <*> f' 4 e <*> f' 5 g <*> f' 6 h
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' (a,a,a,a,a,a,a,a) (b,b,b,b,b,b,b,b) a b@
instance (Applicative f, a~a2, a~a3, a~a4, a~a5, a~a6, a~a7, a~a8, b~b2, b~b3, b~b4, b~b5, b~b6, b~b7, b~b8) => Each f (a,a2,a3,a4,a5,a6,a7,a8) (b,b2,b3,b4,b5,b6,b7,b8) a b where
  each f ~(a,b,c,d,e,g,h,i) = (,,,,,,,) <$> f' (0 :: Int) a <*> f' 1 b <*> f' 2 c <*> f' 3 d <*> f' 4 e <*> f' 5 g <*> f' 6 h <*> f' 7 i
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' (a,a,a,a,a,a,a,a,a) (b,b,b,b,b,b,b,b,b) a b@
instance (Applicative f, a~a2, a~a3, a~a4, a~a5, a~a6, a~a7, a~a8, a~a9, b~b2, b~b3, b~b4, b~b5, b~b6, b~b7, b~b8, b~b9) => Each f (a,a2,a3,a4,a5,a6,a7,a8,a9) (b,b2,b3,b4,b5,b6,b7,b8,b9) a b where
  each f ~(a,b,c,d,e,g,h,i,j) = (,,,,,,,,) <$> f' (0 :: Int) a <*> f' 1 b <*> f' 2 c <*> f' 3 d <*> f' 4 e <*> f' 5 g <*> f' 6 h <*> f' 7 i <*> f' 8 j
    where f' = Lens.indexed f
  {-# INLINE each #-}

#if MIN_VERSION_base(4,4,0)
-- | @'each' :: ('RealFloat' a, 'RealFloat' b) => 'IndexedTraversal' 'Int' ('Complex' a) ('Complex' b) a b@
instance Applicative f => Each f (Complex a) (Complex b) a b where
  each f (a :+ b) = (:+) <$> f' (0 :: Int) a <*> f' (1 :: Int) b
    where f' = Lens.indexed f
  {-# INLINE each #-}
#else
-- | @'each' :: 'IndexedTraversal' 'Int' ('Complex' a) ('Complex' b) a b@
instance (Applicative f, RealFloat a, RealFloat b) => Each f (Complex a) (Complex b) a b where
  each f (a :+ b) = (:+) <$> f' (0 :: Int) a <*> f' 1 b
    where f' = Lens.indexed f
  {-# INLINE each #-}
#endif

-- | @'each' :: 'IndexedTraversal' c ('Map' c a) ('Map' c b) a b@
instance Applicative f => Each f (Map c a) (Map c b) a b where
  each f m = sequenceA $ Map.mapWithKey f' m
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' c ('Map' c a) ('Map' c b) a b@
instance Applicative f => Each f (IntMap a) (IntMap b) a b where
  each f m = sequenceA $ IntMap.mapWithKey f' m
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' c ('HashMap' c a) ('HashMap' c b) a b@
instance Applicative f => Each f (HashMap c a) (HashMap c b) a b where
  each = HashMap.traverseWithKey . Lens.indexed
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' [a] [b] a b@
instance Applicative f => Each f [a] [b] a b
  {- default each -}

-- | @'each' :: 'IndexedTraversal' () ('Identity' a) ('Identity' b) a b@
instance Functor f => Each f (Identity a) (Identity b) a b where
  each f (Identity a) = Identity <$> Lens.indexed f () a
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' () ('Maybe' a) ('Maybe' b) a b@
instance Applicative f => Each f (Maybe a) (Maybe b) a b where
  each f (Just a) = Just <$> Lens.indexed f () a
  each _ Nothing  = pure Nothing
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' ('Seq' a) ('Seq' b) a b@
instance Applicative f => Each f (Seq a) (Seq b) a b
  {- default each -}

-- | @'each' :: 'IndexedTraversal' ['Int'] ('Tree' a) ('Tree' b) a b@
instance Applicative f => Each f (Tree a) (Tree b) a b where
  each pafb = go (BD 0 [] 0 []) where
    go dq (Node a as) = Node <$> Lens.indexed pafb (Foldable.toList dq) a <*> itraverse (\i n -> go (Lens.snoc dq i) n) as
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' ('Vector.Vector' a) ('Vector.Vector' b) a b@
instance Applicative f => Each f (Vector.Vector a) (Vector.Vector b) a b
  {- default each -}

-- | @'each' :: ('Prim' a, 'Prim' b) => 'IndexedTraversal' 'Int' ('Prim.Vector' a) ('Prim.Vector' b) a b@
instance (Applicative f, Prim a, Prim b) => Each f (Prim.Vector a) (Prim.Vector b) a b where
  each f v = Prim.fromListN (Prim.length v) <$> traversed (Indexed f') (Prim.toList v)
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: ('Storable' a, 'Storable' b) => 'IndexedTraversal' 'Int' ('Storable.Vector' a) ('Storable.Vector' b) a b@
instance (Applicative f, Storable a, Storable b) => Each f (Storable.Vector a) (Storable.Vector b) a b where
  each f v = Storable.fromListN (Storable.length v) <$> traversed (Indexed f') (Storable.toList v)
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: ('Unbox' a, 'Unbox' b) => 'IndexedTraversal' 'Int' ('Unboxed.Vector' a) ('Unboxed.Vector' b) a b@
instance (Applicative f, Unbox a, Unbox b) => Each f (Unboxed.Vector a) (Unboxed.Vector b) a b where
  each f v = Unboxed.fromListN (Unboxed.length v) <$> traversed (Indexed f') (Unboxed.toList v)
    where f' = Lens.indexed f
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' 'StrictT.Text' 'StrictT.Text' 'Char' 'Char'@
instance Applicative f =>  Each f StrictT.Text StrictT.Text Char Char where
  each = iso StrictT.unpack StrictT.pack . traversed
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int64' 'LazyT.Text' 'LazyT.Text' 'Char' 'Char'@
instance Applicative f => Each f LazyT.Text LazyT.Text Char Char where
  each = iso LazyT.unpack LazyT.pack . traversed64
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int' 'StrictB.ByteString' 'StrictB.ByteString' 'Word8' 'Word8'@
instance Applicative f => Each f StrictB.ByteString StrictB.ByteString Word8 Word8 where
  each = iso StrictB.unpack StrictB.pack . traversed
  {-# INLINE each #-}

-- | @'each' :: 'IndexedTraversal' 'Int64' 'LazyB.ByteString' 'LazyB.ByteString' 'Word8' 'Word8'@
instance Applicative f => Each f LazyB.ByteString LazyB.ByteString Word8 Word8 where
  each = iso LazyB.unpack LazyB.pack . traversed64
  {-# INLINE each #-}

-- | @'each' :: 'Ix' i => 'IndexedTraversal' i ('Array' i a) ('Array' i b) a b@
instance (Applicative f, Ix i) => Each f (Array i a) (Array i b) a b where
  each f arr = array (bounds arr) <$> traverse (\(i,a) -> (,) i <$> Lens.indexed f i a) (IArray.assocs arr)
  {-# INLINE each #-}

-- | @'each' :: ('Ix' i, 'IArray' 'UArray' a, 'IArray' 'UArray' b) => 'IndexedTraversal' i ('Array' i a) ('Array' i b) a b@
instance (Applicative f, Ix i, IArray UArray a, IArray UArray b) => Each f (UArray i a) (UArray i b) a b where
  each f arr = array (bounds arr) <$> traverse (\(i,a) -> (,) i <$> Lens.indexed f i a) (IArray.assocs arr)
  {-# INLINE each #-}

-- | @'each' :: 'Control.Lens.IndexedSetter' i (i -> a) (i -> b) a b@
instance Settable f => Each f (i -> a) (i -> b) a b where
  each f g = pure (\i -> untaintedDot (Lens.indexed f i) (g i))
  {-# INLINE each #-}

max_stars_repo_name: haasn/lens

max_issues_repo_name: haasn/lens

max_forks_repo_name: haasn/lens


text:
module Mllib.Types
    ( module Reexport
    ) where

import Data.Tree as Reexport (Tree, Forest)

import Numeric.LinearAlgebra.Data    as Reexport (fromList, toList, vector, size)
import Numeric.LinearAlgebra.Data    as Reexport (Vector, R)

-- TODO: replace it
import Numeric.LinearAlgebra.HMatrix as Reexport (norm_2)

import System.Random as Reexport (StdGen, mkStdGen)

-- ======== ANOTHER WAY TO IMPORT ========

{-
module Mllib.Types
    ( norm_2
    , fromList
    , Vector
    , R
    ) where

import Numeric.LinearAlgebra.HMatrix (norm_2)
import Numeric.LinearAlgebra.Data (fromList)
import Numeric.LinearAlgebra.Data (Vector, R)
-}


text:
module Main where

import Numeric.LinearAlgebra
import Control.Arrow
import Common
import Forward
import BackProp
import AutoEncoder
import ActivFunc
import Other
import Mnist

main :: IO ()
main = do
  tm <- parseCsvToMatrixR "data/mnist_test_100.csv"
  let (ty, tx) = tr *** tr $ mnistRead 10 tm
  m <- parseCsvToMatrixR "data/mnist_train_600.csv"
  let (y, x) = tr *** tr $ mnistRead 10 m
  ws <- genWeights [784, 10]

  -- let pws = preTrains 0.5 1000 ramps x ws
  nws <- last . take 1000 $ iterateM (sgdMethod 100 (x, y) $ backPropClassification 0.5 ramps) ws

  -- putStrLn "not trained outputs"
  -- print . tr . classesToLabels . tr . forwardClassification rampC ws $ tx
  -- TODO: use function
  let o = tr . classesToLabels . tr . forwardClassification rampC nws $ tx
  let t = tr . classesToLabels . tr $ ty
  let l = head . toLists $ t - o
  let a = length l
  let s = length . filter (==0) $ l
  print $ fromIntegral s / fromIntegral a

  -- putStrLn "pretrainined outputs"
  -- print . classesToLabels . tr . forwardClassification rampC npws $ x

max_stars_repo_name: pupuu/deep-neuralnet

max_issues_repo_name: pupuu/deep-neuralnet

max_forks_repo_name: pupuu/deep-neuralnet


text:
-- Test of the computational correctness of my current FFT candidate.
--
-- Original author: David Banas <capn.freako@gmail.com>
-- Original date:   October 3, 2015
--
-- Copyright (c) 2015 David Banas; all rights reserved World wide.
--
-- I'm waiting for Conal to find time to look into my latest non-termination
-- failure through his compiler. While I do, I'm attempting, here, to verify
-- the computational correctness of my current FFT candidate.

{-# LANGUAGE GADTs #-}
{-# LANGUAGE TemplateHaskell #-}
{-# LANGUAGE TypeSynonymInstances #-}
{-# LANGUAGE MultiParamTypeClasses #-}
{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE FlexibleContexts #-}

module Main where

import Prelude hiding ({- id,(.), -}foldl,foldr,sum,product,zipWith,reverse,and,or,scanl,minimum,maximum)

import Control.Applicative
import Control.Arrow
import Control.Monad (forM_, unless)
import Data.Complex
import Data.Foldable (Foldable, sum, foldl')
import Data.Newtypes.PrettyDouble
import System.Exit (exitFailure)
import TypeUnary.Nat (IsNat, Nat(..), nat, N2, N3, N4, N5)  -- , N6)

-- import Test.QuickCheck (choose, vectorOf, elements, collect)
import Test.QuickCheck (choose, vectorOf)
import Test.QuickCheck.Arbitrary
import Test.QuickCheck.All (quickCheckAll)

import Circat.Scan (lproducts, LScan)
import qualified Circat.Pair as P
import qualified Circat.RTree as RT
import Circat.RTree (bottomSplit)

type RTree = RT.Tree

-- FFT, as a class
-- (The LScan constraint comes from the use of 'lproducts', in 'addPhase'.)
class (LScan f) => FFT f a where
    fft  :: f a -> f a  -- Computes the FFT of a functor.

-- Note that this definition of the FFT instance for Pair assumes DIT.
-- How can we eliminate this assumption and make this more general?
instance (RealFloat a, Applicative f, Foldable f, Num (f (Complex a)), FFT f (Complex a)) => FFT P.Pair (f (Complex a)) where
    fft = P.inP (uncurry (+) &&& uncurry (-)) . P.secondP addPhase . fmap fft

instance (IsNat n, RealFloat a) => FFT (RTree n) (Complex a) where
    fft = fft' nat
        where   fft' :: (RealFloat a) => Nat n -> RTree n (Complex a) -> RTree n (Complex a)
                fft' Zero     = id
                fft' (Succ _) = inDIT fft
                    where   inDIT g  = RT.toB . g . bottomSplit

-- Adds the proper phase adjustments to a functor containing Complex RealFloats,
-- and instancing Num.
addPhase :: (Applicative f, Foldable f, LScan f, RealFloat a, Num (f (Complex a))) => f (Complex a) -> f (Complex a)
addPhase = liftA2 (*) id phasor
  where phasor f = fst $ lproducts (pure phaseDelta)
          where phaseDelta = cis ((-pi) / fromIntegral n)
                n          = flen f

-- Gives the "length" (i.e. - number of elements in) of a Foldable.
-- (Soon, to be provided by the Foldable class, as "length".)
flen :: (Foldable f) => f a -> Int
flen = foldl' (flip ((+) . const 1)) 0

-- Test config.
realData :: [[PrettyDouble]]
realData = [  [1.0,   0.0,   0.0,   0.0]  -- Delta
            , [1.0,   1.0,   1.0,   1.0]  -- Constant
            , [1.0,  -1.0,   1.0,  -1.0]  -- Nyquist
            , [1.0,   0.0,  -1.0,   0.0]  -- Fundamental
            , [0.0,   1.0,   0.0,  -1.0]  -- Fundamental w/ 90-deg. phase lag
           ]
complexData :: [[Complex PrettyDouble]]
complexData = map (map (:+ 0.0)) realData

myTree2 :: [a] -> RTree N2 a
myTree2 [w, x, y, z] = RT.tree2 w x y z
myTree2 _            = error "Something went horribly wrong!"

myTree3 :: [a] -> RTree N3 a
myTree3 [a, b, c, d, e, f, g, h] = RT.tree3 a b c d e f g h
myTree3 _            = error "Something went horribly wrong!"

myTree4 :: [a] -> RTree N4 a
myTree4 [a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p] = RT.tree4 a b c d e f g h i j k l m n o p
myTree4 _            = error "Something went horribly wrong!"

myTree5 :: [a] -> RTree N5 a
myTree5 [a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p,
         a', b', c', d', e', f', g', h', i', j', k', l', m', n', o', p'] =
            RT.tree5 a b c d e f g h i j k l m n o p a' b' c' d' e' f' g' h' i' j' k' l' m' n' o' p'
myTree5 _            = error "Something went horribly wrong!"

-- myTree6 :: [a] -> RTree N6 a
-- myTree6 [a1, b1, c1, d1, e1, f1, g1, h1, i1, j1, k1, l1, m1, n1, o1, p1,
--          a2, b2, c2, d2, e2, f2, g2, h2, i2, j2, k2, l2, m2, n2, o2, p2,
--          a3, b3, c3, d3, e3, f3, g3, h3, i3, j3, k3, l3, m3, n3, o3, p3,
--          a4, b4, c4, d4, e4, f4, g4, h4, i4, j4, k4, l4, m4, n4, o4, p4] =
--             RT.tree6 a1 b1 c1 d1 e1 f1 g1 h1 i1 j1 k1 l1 m1 n1 o1 p1
--                      a2 b2 c2 d2 e2 f2 g2 h2 i2 j2 k2 l2 m2 n2 o2 p2
--                      a3 b3 c3 d3 e3 f3 g3 h3 i3 j3 k3 l3 m3 n3 o3 p3
--                      a4 b4 c4 d4 e4 f4 g4 h4 i4 j4 k4 l4 m4 n4 o4 p4
-- myTree6 _            = error "Something went horribly wrong!"

-- Discrete Fourier Transform (DFT) (our "truth" reference)
-- O(n^2)
--
dft :: RealFloat a => [Complex a] -> [Complex a]
dft xs = [ sum [ x * exp((0.0 :+ (-1.0)) * 2 * pi / lenXs * fromIntegral(k * n))
                 | (x, n) <- Prelude.zip xs [0..]
               ]
           | k <- [0..(length xs - 1)]
         ]
    where lenXs = fromIntegral $ length xs

-- QuickCheck types & propositions
newtype FFTTestVal = FFTTestVal {
    getVal :: [Complex PrettyDouble]
} deriving (Show)
instance Arbitrary FFTTestVal where
    arbitrary = do
        xs <- vectorOf 32 $ choose (-1.0::Double, 1.0)
        let zs = map ((:+ 0) . PrettyDouble) xs
        return $ FFTTestVal zs

prop_fft_test_N2 :: FFTTestVal -> Bool
prop_fft_test_N2 testVal = fft (myTree2 zs) == RT.fromList (dft zs)
    where zs = take 4 $ getVal testVal

prop_fft_test_N3 :: FFTTestVal -> Bool
prop_fft_test_N3 testVal = fft (myTree3 zs) == RT.fromList (dft zs)
    where zs = take 8 $ getVal testVal

prop_fft_test_N4 :: FFTTestVal -> Bool
prop_fft_test_N4 testVal = fft (myTree4 zs) == RT.fromList (dft zs)
    where zs = take 16 $ getVal testVal

prop_fft_test_N5 :: FFTTestVal -> Bool
prop_fft_test_N5 testVal = fft (myTree5 zs) == RT.fromList (dft zs)
    where zs = take 32 $ getVal testVal

-- Test definitions & choice
basicTest :: IO ()
basicTest = forM_ complexData (\x -> do
                putStr "\nTesting input: "
                print x
                putStr "Expected output: "
                print $ dft x
                putStr "Actual output:   "
                print $ fft $ myTree2 x
                )

-- This weirdness is required, as of GHC 7.8.
return []

runTests :: IO Bool
runTests = $quickCheckAll
-- End weirdness.

advancedTest :: IO ()
advancedTest = do
    allPass <- runTests -- Run QuickCheck on all prop_ functions
    unless allPass exitFailure

main :: IO ()
-- main = basicTest
main = advancedTest

max_stars_repo_name: capn-freako/fft-ccc

max_issues_repo_name: capn-freako/fft-ccc

max_forks_repo_name: capn-freako/fft-ccc


text:
module Common (
    genWeights,
    initWeights,
    inputWithBias,
    weightWithoutBias,
    pickupSets
  ) where

import Numeric.LinearAlgebra
import System.Random.Shuffle

genWeights :: [Int] -> IO [Matrix R]
genWeights ns = mapM genWeight ts
  where ts = zip ns (tail ns)

genWeight :: (Int, Int) -> IO (Matrix R)
genWeight (i, o) = randn o (i + 1) -- +1: bias

initWeights :: [Int] -> [Matrix R]
initWeights ns = fmap initWeight ts
  where ts = zip ns (tail ns)

initWeight :: (Int, Int) -> Matrix R
initWeight (i, o) = matrix (i + 1) . replicate (o * (i + 1)) $ 0 -- +1: bias

inputWithBias :: Matrix R -> Matrix R
inputWithBias v = v === konst 1 (1, cols v) -- [1,1,..1]: bias

weightWithoutBias :: Matrix R -> Matrix R
weightWithoutBias w = w ?? (All, DropLast 1)

pickupSets :: Int -> (Matrix R, Matrix R) -> IO (Matrix R, Matrix R)
pickupSets n (x, y) = do
  s <- shuffleM [0..cols y - 1]
  let ps = take n s in
    return (x ¿ ps, y ¿ ps)


text:
module STCR2Z1T0_1 where

import           Control.Monad             as M
import           Data.Array.Repa           as R
import           Data.Binary               (decodeFile)
import           Data.Complex
import           Data.List                 as L
import           DFT.Plan
import           FokkerPlanck.DomainChange
import           FokkerPlanck.MonteCarlo
import           FokkerPlanck.Pinwheel
import           Image.IO
import           Image.Transform
import           STC.CompletionField
import           System.Directory
import           System.Environment
import           System.FilePath
import           Types
import           Utils.Array

main = do
  args@(numPointStr:numOrientationStr:sigmaStr:taoStr:lenStr:initStr:numTrailStr:maxTrailStr:theta0FreqsStr:thetaFreqsStr:initDistStr:histFilePath:alphaStr:pinwheelFlagStr:numThreadStr:_) <-
    getArgs
  print args
  let numPoint = read numPointStr :: Int
      numOrientation = read numOrientationStr :: Int
      sigma = read sigmaStr :: Double
      tao = read taoStr :: Double
      len = read lenStr :: Int
      init = read initStr :: (Double, Double, Double, Double, Double, Double)
      numTrail = read numTrailStr :: Int
      maxTrail = read maxTrailStr :: Int
      theta0Freq = read theta0FreqsStr :: Double
      theta0Freqs = [-theta0Freq .. theta0Freq]
      thetaFreq = read thetaFreqsStr :: Double
      thetaFreqs = [-thetaFreq .. thetaFreq]
      initDist = read initDistStr :: [R2S1RPPoint]
      alpha = read alphaStr :: Double
      pinwheelFlag = read pinwheelFlagStr :: Bool
      numThread = read numThreadStr :: Int
      sourceDist = L.take 1 initDist
      sinkDist = L.drop 1 initDist
      folderPath = "output/test/STCR2Z1T0_1"
  flag <- doesFileExist histFilePath
  arrR2Z1T0' <-
    if pinwheelFlag
      then computeR2Z1T0Array numPoint numPoint alpha thetaFreqs theta0Freqs
      else if flag
             then getNormalizedHistogramArr <$> decodeFile histFilePath
             else solveMonteCarloR2Z1T0
                    numThread
                    numTrail
                    maxTrail
                    numPoint
                    numPoint
                    sigma
                    tao
                    len
                    theta0Freqs
                    thetaFreqs
                    histFilePath
                    (emptyHistogram
                       [ numPoint
                       , numPoint
                       , L.length theta0Freqs
                       , L.length thetaFreqs
                       ]
                       0)
  let arrR2Z1T0 =
        computeUnboxedS .
        pad [numPoint, numPoint, L.length theta0Freqs, L.length thetaFreqs] 0 .
        downsample [2, 2, 1, 1] $
        arrR2Z1T0'
  createDirectoryIfMissing True folderPath
  plan <- makeR2Z1T0Plan emptyPlan arrR2Z1T0
  sourceDistArr <-
    computeInitialDistributionR2Z1T0
      plan
      numPoint
      numPoint
      thetaFreqs
      theta0Freqs
      sourceDist
  sinkDistArr <-
    computeInitialDistributionR2Z1T0
      plan
      numPoint
      numPoint
      thetaFreqs
      theta0Freqs
      sinkDist
  arrR2Z1T0F <- dftR2Z1T0 plan . makeFilterR2Z1T0 $ arrR2Z1T0
  arrR2Z1T0TRF <-
    dftR2Z1T0 plan . makeFilterR2Z1T0 . timeReverseR2Z1T0 thetaFreqs theta0Freqs $
    arrR2Z1T0
  -- Source field
  sourceArr <- convolveR2Z1T0 plan arrR2Z1T0F sourceDistArr
  sourceR2 <- R.sumP . R.sumS . rotate4D . rotate4D $ sourceArr
  let sourceField =
        computeS . R.extend (Z :. (1 :: Int) :. All :. All) $ sourceR2
  plotImageRepaComplex (folderPath </> "Source.png") . ImageRepa 8 $ sourceField
  -- Sink field
  sinkArr <- convolveR2Z1T0 plan arrR2Z1T0TRF sinkDistArr
  sinkR2 <- R.sumP . R.sumS . rotate4D . rotate4D $ sinkArr
  let sinkField = computeS . R.extend (Z :. (1 :: Int) :. All :. All) $ sinkR2
  plotImageRepaComplex (folderPath </> "Sink.png") . ImageRepa 8 $ sinkField
  -- Completion Filed
  completionFiled <- convolveR2Z1' plan thetaFreqs theta0Freqs sourceArr sinkArr
  completionFiledR2 <- R.sumP . R.sumS . rotate4D . rotate4D $ completionFiled
  plotImageRepa (folderPath </> "Completion.png") .
    ImageRepa 8 .
    computeS . R.extend (Z :. (1 :: Int) :. All :. All) . R.map magnitude $
    completionFiledR2


text:
module Main where
import Plot
import Graphics.GD
import Data.Complex

imageDim :: Size
imageDim = (600,600)

-- How many iterations before  we declare the point doesn't diverge
maxIter = 96

-- How close successive tests must be before we declare they have converged
minDiff = 0.0001

-- formula to iterate
f z = 3 * z^4 - 1
-- derivative
f' z = 12 * ( z^3)

-- Coordinates in the complex plane to graph
-- First coordinate must be top left.
mbWindow :: Window
mbWindow = ( (-2.0, 2.0), (2.0, -2.0) ) -- Full view

main :: IO ()
main = drawPlot drawNewton imageDim mbWindow "newtons-method.png"

drawNewton :: Coordinate -> Color
drawNewton (x,y) = colorIterations $ newton (x :+ y) 0

-- Count number of iterations for the given point to diverge
-- Non-diverging points return 0
newton :: Complex Double -- Coordinate to test
       -> Int -- Current iteration
       -> Int -- Iterations before diverging
newton z iter
    | iter > maxIter = 0
    | ( abs $ magnitude r - magnitude z) > minDiff = newton r (iter + 1)
    | otherwise = iter
    where r = z - f z / f' z

colorIterations :: Int -> Color
colorIterations x
    | x > maxIter = rgb 255 255 255
    | otherwise = rgb (x^2) x (x*3)

max_stars_repo_name: typeclasses/fractal-sets

max_issues_repo_name: typeclasses/fractal-sets

max_forks_repo_name: typeclasses/fractal-sets