text:
{-# LANGUAGE DeriveDataTypeable #-}
{-# LANGUAGE OverloadedStrings  #-}
{-# LANGUAGE RecordWildCards    #-}
{-# LANGUAGE TemplateHaskell    #-}

-----------------------------------------------------------------------------
-- |
-- Module      :
-- Copyright   :  (C) 2012 Soostone Inc
-- License     :  All Rights Reserved
--
-- Maintainer  :  Ozgun Ataman <oz@soostone.com>
-- Stability   :  experimental
--
-- Defines a standard set of statistics computed for collections of
-- data. For theory around combining cells monoidally, see:
--
-- http://stackoverflow.com/questions/7753002/adding-combining-standard-deviations
-- http://en.wikipedia.org/wiki/Algorithms_for_calculating_variance#On-line_algorithm
----------------------------------------------------------------------------

module Data.StatVar
    ( StatVar (..)
    , svStdDev
    , svVariance

    -- * Operating on StatVars
    , newStatVar
    , updateSV
    , updateSVRemExt
    , trimExt
    , isExtremeVal
    , listToStatVar
    , expandStatVar
    ) where

-------------------------------------------------------------------------------
import           Control.Applicative
import           Control.DeepSeq
import           Data.Aeson
import           Data.DeriveTH
import           Data.Generics
import           Data.List           (foldl')
import           Data.Monoid
import           Data.SafeCopy
import           Data.Serialize
import           Data.Text           (Text)
import qualified Data.Vector         as V
import           Statistics.Sample
-------------------------------------------------------------------------------

-------------------------------------------------------------------------------
-- | A set of standard stats over a collection of numbers
data StatVar = StatVar {
      svMax   :: !Double
    , svMin   :: !Double
    , svSum   :: !Double
    , svCount :: !Int
    , svMean  :: !Double
    , svSOS   :: !Double
    -- ^ sum of squares of differences from the (current) mean. Needed
    -- for Welford's method of variance estimation.
    } deriving (Eq,Show,Read,Ord,Typeable,Data)

-------------------------------------------------------------------------------
$(deriveSafeCopy 0 'base ''StatVar)
$(derives [makeSerialize] [''StatVar])
$(derives [makeNFData] [''StatVar])

-------------------------------------------------------------------------------
instance FromJSON StatVar where
    parseJSON (Object v) = StatVar
      <$> v .: "mx"
      <*> v .: "mn"
      <*> v .: "sm"
      <*> v .: "cnt"
      <*> v .: "avg"
      <*> v .: "sos"

-------------------------------------------------------------------------------
instance ToJSON StatVar where
    toJSON StatVar{..} = object
      [ "mx" .= svMax
      , "mn" .= svMin
      , "sm" .= svSum
      , "cnt" .= svCount
      , "avg" .= svMean
      , "sos" .= svSOS
      ]

-------------------------------------------------------------------------------
-- | Extract components out of a 'StatVar' for convenience in human
-- readable contexts.
expandStatVar :: StatVar -> [(String, Double)]
expandStatVar st =
    [ ("max", svMax st)
    , ("min", svMin st)
    , ("mean", svMean st)
    , ("count", fromIntegral $ svCount st)
    , ("sum", svSum st)
    , ("stdDev", svStdDev st)
    , ("var", svVariance st)
    ]

-------------------------------------------------------------------------------
-- | Convert a list of numerical data into a 'StatVar'
listToStatVar [] = Nothing
listToStatVar (x:xs) = Just $ foldl' updateSV (newStatVar x) xs

-------------------------------------------------------------------------------
newStatVar :: Double -> StatVar
newStatVar x = StatVar { svMax = x
                       , svMin = x
                       , svSum = x
                       , svCount = 1
                       , svMean = x
                       , svSOS = 0 }

-------------------------------------------------------------------------------
-- | Update the stats variable with a new datapoint
updateSV :: StatVar -> Double -> StatVar
updateSV StatVar{..} x =
    StatVar { svMax = mx
            , svMin = mn
            , svSum = sm
            , svCount = cnt
            , svMean = mean
            , svSOS = sos }
    where
      mx = max svMax x
      mn = min svMin x
      sm = svSum + x
      cnt = svCount + 1
      delta = x - svMean
      mean = sm / fromIntegral cnt
      sos = svSOS + delta * (x - mean)

-- | Add new value while trimming extreme cases.
--
-- Require that at least n samples have been collected so far. If new
-- value is outside of the 3-sigma range, it will be replaced by the
-- mean value so far observed.
updateSVRemExt n sv x = updateSV sv $ trimExt n sv x

-- | Maybe trim the value based on given StatVar with historical
-- stats.
trimExt :: Int
        -- ^ Min # of samples required for trim
        -> StatVar
        -- ^ Historical stats
        -> Double
        -- ^ New observed value
        -> Double
        -- ^ (Possibly) trimmed value
trimExt n sv x = if isExtremeVal n sv x
                 then svMean sv
                 else x

-- | Is the value extreme (+- 3 sigmas)?
isExtremeVal :: Int -> StatVar -> Double -> Bool
isExtremeVal n sv x =
    svCount sv >= n && (x > ulimit || x < llimit)
    where
      sdev = svStdDev sv
      mean = svMean sv
      ulimit = mean + 3 * sdev
      llimit = mean - 3 * sdev

-------------------------------------------------------------------------------
svStdDev :: StatVar -> Double
svStdDev = sqrt . svVariance

-------------------------------------------------------------------------------
svVariance :: StatVar -> Double
svVariance st
    | svCount st >= 2 = svSOS st / (fromIntegral $ svCount st - 1)
    | otherwise = 0

-------------------------------------------------------------------------------
instance Monoid StatVar where
    mempty = StatVar 0 0 0 0 0 0

    mappend a b = StatVar {
                        svMax = max (svMax a) (svMax b)
                      , svMin = min (svMin a) (svMin b)
                      , svSum = sum'
                      , svCount = cnt
                      , svMean = mean'
                      , svSOS = sos
                      }
        where
          mean' = sum' / cnt'
          cnt = svCount a + svCount b
          cnt' = fromIntegral cnt
          sum' = svSum a + svSum b

          -- variance calculation; back-compute sos from combined variance
          sos = s1 + s2
          var = (s1 + s2) / (cnt' - 1)
          s1 = (fromIntegral $ svCount a - 1) * svVariance a +
               (fromIntegral $ svCount b - 1) * svVariance b
          s2 = fromIntegral (svCount a) * ((svMean a - mean') ** 2) +
               fromIntegral (svCount b) * ((svMean b - mean') ** 2)

                                 -------------
                                 -- Testing --
                                 -------------

va = V.fromList [3.0,7,2,4,5,7,8]

vb = V.fromList [1.0,6,9,5,3,11,24]

v = V.concat [va, vb]

varVa = varianceUnbiased va
varVb = varianceUnbiased vb
varV = varianceUnbiased v

stA = mkStat va
stB = mkStat vb

-------------------------------------------------------------------------------
-- | stAll should be equal to stComb if the method to merge variances is right.
stAll = mkStat v

-------------------------------------------------------------------------------
stComb = stA `mappend` stB

-------------------------------------------------------------------------------
mkStat va = StatVar mx mn sm (floor cnt) avg sos
    where
      sos = var * (cnt - 1)
      [mx, mn, sm, cnt, avg, var] = map ($ va)
          [V.maximum, V.minimum, V.sum, fromIntegral . V.length, mean, varianceUnbiased]

buildSV = foldl updateSV (newStatVar 175) (replicate 32 133)

max_stars_repo_name: cartazio/analytics

max_issues_repo_name: cartazio/analytics

max_forks_repo_name: cartazio/analytics


text:
module Stepik where

import Prelude hiding (Monoid, mempty, mappend, mconcat, Functor, fmap)
import Data.Char
import Data.Function
import Data.List (transpose)
import Data.Functor
-- import Data.Semigroup ((<>))
import Control.Monad (ap, liftM)

import Data.Time.Clock
import Data.Time.Format
-- import System.Locale

main = undefined

{- Create own operator -}
infixl 6 *+*
(*+*) a b = a ^ 2 + b ^ 2

infixl 6 |-|
x |-| y = abs (x - y)

twoDigits2Int :: Char -> Char -> Int
twoDigits2Int x y = if isDigit x && isDigit y then read [x, y] :: Int else 100

doubleFact :: Integer -> Integer
doubleFact 1 = 1
doubleFact 2 = 2
doubleFact n = n * doubleFact(n-2)

{- Negative Factorial -}
fibonacci :: Integer -> Integer
fibonacci 0 = 0
fibonacci 1 = 1
fibonacci n
    | n < 0 = (-1)^(abs n + 1) * fibonacci (abs n)
    | otherwise = fibonacci (n - 1) + fibonacci (n - 2)

{- Tail recursion -}
factorial :: Integer -> Integer
factorial n
    | n >= 0 = helper 1 n
    | otherwise = error "Arg must be >= 0!"

helper :: Integer -> Integer -> Integer
helper acc 0 = acc
helper acc n = helper (acc * n) (n - 1)

{- Negative Fibonacci tail recursion -}
fibonacciTail :: Integer -> Integer
fibonacciTail n
    | n < 0 = (-1)^(abs n + 1) * fibonacciTail (abs n)
    | otherwise = helperFib 0 1 n

helperFib :: Integer -> Integer -> Integer -> Integer
helperFib p c 1 = c
helperFib p c n = helperFib c (c+p) (n-1)

{- Negative Fibonacci tail recursion `let ... in ...`. Task -}
fibonacciTail' :: Integer -> Integer
fibonacciTail' n
    | n < 0 = (-1)^(abs n + 1) * fibonacciTail' (abs n)
    | otherwise = let
        helperFib p c 0 = 0
        helperFib p c 1 = c
        helperFib p c n = helperFib c (c + p) (n - 1)
    in helperFib 0 1 n

fibonacciTail'' :: Integer -> Integer
fibonacciTail'' n
    | n < 0 = (-1)^(abs n + 1) * fibonacciTail'' (abs n)
    | otherwise = helper 0 1 n
    where
        helper p c 0 = 0
        helper p c 1 = c
        helper p c n = helper c (c + p) (n - 1)

{- Task -}
seqA :: Integer -> Integer
seqA 0 = 1
seqA 1 = 2
seqA 2 = 3
seqA n = seqA (n-1) + seqA (n-2) - 2 * seqA (n-3)

seqA' :: Integer -> Integer
seqA' n = helper 3 3 2 n
    where helper curr prev multip n
            | n == 0 = 1
            | n == 1 = 2
            | n == 2 = 3
            | n == 3 = curr
            | n > 3 = helper (curr + prev - 2 * multip) curr prev (n - 1)

seqA'' :: Integer -> Integer
seqA'' n = let
    helper curr prev multip 0 = 1
    helper curr prev multip 1 = 2
    helper curr prev multip 2 = 3
    helper curr prev multip 3 = curr
    helper curr prev multip n = helper (curr + prev - 2 * multip) curr prev (n - 1)
    in helper 3 3 2 n

seqA''' :: Integer -> Integer
seqA''' n = helper 3 3 2 n
    where
        helper curr prev multip 0 = 1
        helper curr prev multip 1 = 2
        helper curr prev multip 2 = 3
        helper curr prev multip 3 = curr
        helper curr prev multip n = helper (curr + prev - 2 * multip) curr prev (n - 1)

{- Task -}
sum'n'count :: Integer -> (Integer, Integer)
sum'n'count x = (s, l)
        where
            xs = show $ abs x
            l = toInteger $ length xs
            s = toInteger $ sum $ map digitToInt xs

{- Task -}
integration :: (Double -> Double) -> Double -> Double -> Double
integration f a b = let
                        xs = [a + j * h | j <- [0..1000]]
                        h = (b - a) / 1000
                    in sum [( ( f (xs !! i) + f (xs !! (i + 1)) ) / 2 ) * ( xs !! (i + 1) - xs !! i ) | i <- [0..999]]

integration' :: (Double -> Double) -> Double -> Double -> Double
integration' f a b = h * ( (f a + f b) / 2 + (sum $ map f xs) )
    where
        n = 1000
        h = (b - a) / n
        xs = [a + h * x | x <- [0..n-1]]

{- Task -}
getSecondFrom :: t1 -> t2 -> t3 -> t2
getSecondFrom x y z = y

{- Task -}
multSecond = g `on` h

g :: (Num a) => a -> a -> a
g x y = x * y

h :: (a, b) -> b
h (x, y) = y

{- Task -}
on3 :: (b -> b -> b -> c) -> (a -> b) -> a -> a -> a -> c
on3 op f x y z = op (f x) (f y) (f z)

{- Task -}
doItYourself = doItYourselfF . doItYourselfG . doItYourselfH

doItYourselfF = logBase 2
doItYourselfG = (^ 3)
doItYourselfH = max 42

{- Task -}
swap' = uncurry (flip (,))

{- Task -}
class Printable a where
    toString :: a -> String

instance Printable Bool where
  toString True = "true"
  toString False = "false"

instance Printable () where
  toString _ = "unit type"

instance (Printable a, Printable b) => Printable (a, b) where
  toString (x,y) = "(" ++ toString x ++ "," ++ toString y ++ ")"

{- Task -}
class KnownToGork a where
    stomp :: a -> a
    doesEnrageGork :: a -> Bool

instance KnownToGork Int where
  stomp _ = 4
  doesEnrageGork = (< 10)

class KnownToMork a where
    stab :: a -> a
    doesEnrageMork :: a -> Bool

instance KnownToMork Int where
  stab _ = 3
  doesEnrageMork = (> 5)

class (KnownToGork a, KnownToMork a) => KnownToGorkAndMork a where
    stompOrStab :: a -> a

instance KnownToGorkAndMork Int where
    stompOrStab x
        | doesEnrageGork x && doesEnrageMork x = stomp $ stab x
        | doesEnrageGork x = stab x
        | doesEnrageMork x = stomp x
        | otherwise = x

{- Task -}
class (Enum a, Bounded a, Eq a) => SafeEnum a where
  ssucc :: a -> a
  ssucc x | maxBound == x = minBound | otherwise = succ x
--  ssucc maxBound = minBound
--  ssucc x = succ x

  spred :: a -> a
  spred x | minBound == x = maxBound | otherwise = pred x
--  spred minBound = maxBound
--  spred x = pred x

instance SafeEnum Bool
instance SafeEnum Int

{- Task -}
avg :: Int -> Int -> Int -> Double
avg x y z = sum(map fromIntegral [x,y,z]) / 3

{- Task -}
--bar x y z = x + y
--foo a b = bar a a (a + b)
--value = foo (3 * 10) (5 - 2)

{-
value ~> foo (3 * 10) (5 - 2) ~> bar (3 * 10) (3 * 10) ((3* 10) + (5 - 2)) ~>
(3 * 10) + (3 * 10) ~> 30 + (3 * 10) ~> 30 + 30 ~> 60
-}

{- Task -}
addTwoElements :: a -> a -> [a] -> [a]
addTwoElements x y lst = x : y : lst

{- Task -}
nTimes:: a -> Int -> [a]
nTimes x 0 = []
nTimes x n = x : nTimes x (n-1)

{- Task -}
oddsOnly :: Integral a => [a] -> [a]
oddsOnly xs = [x | x <- xs, odd x]

oddsOnly' :: Integral a => [a] -> [a]
oddsOnly' [] = []
oddsOnly' (x:xs)
  | odd x = x : oddsOnly' xs
  | otherwise = oddsOnly' xs

reverse' :: Integral a => [a] -> [a]
reverse' l = rev l [] where
  rev [] a = a
  rev (x:xs) a = rev xs (x:a)

{- Task -}
isPalindrome :: Eq a => [a] -> Bool
isPalindrome xs = helper xs True where
  helper _ False = False
  helper [] b = b
  helper [_] b = b
  helper xs b = helper (tail $ init xs) (head xs == last xs)

{- Task -}
sum3 :: Num a => [a] -> [a] -> [a] -> [a]
sum3 as bs cs = map sum $ transpose [as, bs, cs]

{- Task -}
groupElems :: Eq a => [a] -> [[a]]
groupElems xs = let
                  helper [] res = res
                  helper (x:xs) [] = helper xs [[x]]
                  helper (x:xs) res
                    | (last $ lres) == x = helper xs (init res ++ [x:lres])
                    | otherwise = helper xs (res ++ [[x]])
                    where lres = last res
                in helper xs []

filter' :: (a -> Bool) -> [a] -> [a]
filter' p [] = []
filter' p (x:xs)
  | p x = x : filter' p xs
  | otherwise = filter' p xs

{- Task -}
readDigits :: String -> (String, String)
readDigits x = span isDigit x

{- Task -}
filterDisj :: (a -> Bool) -> (a -> Bool) -> [a] -> [a]
filterDisj p1 p2 xs = filter (\x -> p1 x || p2 x) xs

{- Task -}
qsort :: Ord a => [a] -> [a]
qsort xs | length xs < 2 = xs
qsort xs = let
              pivot = head xs
              less = [x | x <- tail xs, x <= pivot]
              greater = [x | x <- tail xs, x > pivot]
           in (qsort less) ++ [pivot] ++ (qsort greater)

qsort' :: Ord a => [a] -> [a]
qsort' xs | length xs < 2 = xs
qsort' xs = let
              pivot = head xs
              less = filter (<= pivot) (tail xs)
              greater = filter (> pivot) (tail xs)
           in (qsort' less) ++ [pivot] ++ (qsort' greater)

qsort'' :: Ord a => [a] -> [a]
qsort'' [] = []
qsort'' (x:xs) = (qsort'' $ filter (<=x) xs) ++ [x]  ++ (qsort'' $ filter (>x) xs)

{- Task -}
squares'n'cubes :: Num a => [a] -> [a]
squares'n'cubes xs = concat $ map (\(x,y) -> [x,y]) (zip (map (^2) xs) (map (^3) xs))

squares'n'cubes' :: Num a => [a] -> [a]
squares'n'cubes' = concat . map (\x -> [x^2,x^3])

squares'n'cubes'' :: Num a => [a] -> [a]
squares'n'cubes'' = concatMap (\x -> [x^2,x^3])

{- Task -}
rotations :: Int -> [a] -> [[a]]
rotations len xs = take len (iterate (\(y:ys) -> ys ++ [y]) xs)

perms :: [a] -> [[a]]
perms []        = [[]]
perms il@(x:xs) = concatMap ((rotations len).(x:)) (perms xs)
                  where len = length il

{- Task -}
delAllUpper :: String -> String
delAllUpper = unwords . filter (any isLower) . words

{- Task -}
max3 :: Ord a => [a] -> [a] -> [a] -> [a]
max3 xs ys zs = map (\(x,y,z) -> max x $ max y z) (zip3 xs ys zs)

max3' :: Ord a => [a] -> [a] -> [a] -> [a]
max3' = zipWith3 (\a b c -> a `max` b `max` c)

max3'' :: Ord a => [a] -> [a] -> [a] -> [a]
max3'' = zipWith3 (\x y z -> maximum [x,y,z])

{- Task -}
fibStream :: [Integer]
fibStream = 0 : 1 : zipWith (+) fibStream (tail fibStream)

{- Task -}
data Odd = Odd Integer
  deriving (Eq, Show)

instance Enum Odd where
  succ (Odd x) = Odd (x + 2)
  pred (Odd x) = Odd (x - 2)

  fromEnum (Odd x) = fromInteger $ x `div` 2
  toEnum x = Odd (toInteger x * 2 + 1)

  enumFrom = iterate succ
  enumFromThen (Odd x) (Odd y) = map Odd [x, y ..]
  enumFromTo (Odd x) (Odd y) = map Odd [x, x + 2 .. y]
  enumFromThenTo (Odd x) (Odd y) (Odd z) = map Odd [x , y .. z]

{- Folding: Remember -}
{- foldr (-) x [2,1,5] ~> (2 - (1 - (5 - x))) -}
{- foldl (-) x [2,1,5] ~> (((x - 5) - 1) - 2) -}

{- Task -}
meanList :: [Double] -> Double
meanList = (\(x, y) -> x / y) . foldr (\x (s, l) -> (x + s, 1 + l)) (0, 0)

{- Task -}
evenOnly :: [a] -> [a]
evenOnly = fst . foldl (\(s, p) x -> if even p then (s ++ [x], p + 1) else (s, p + 1))  ([], 1)

-- Endless list
evenOnly' :: [a] -> [a]
evenOnly' = foldr (\(n,x) xs -> if even n then x:xs else xs) [] . zip [1, 2..]

maximum' :: (Ord a) => [a] -> a
maximum' = foldl1 max

{- Task -}
lastElem :: [a] -> a
lastElem = foldl1 (\_ y -> y)

lastElem' :: [a] -> a
lastElem' = foldl1 seq

{- Task -}
unfoldr :: (b -> Maybe (a, b)) -> b -> [a]
unfoldr f ini = helper (f ini) where
  helper (Just (x, ini')) = x : unfoldr f ini'
  helper Nothing          = []

revRange :: (Char,Char) -> [Char]
revRange (a, b) = unfoldr g b
  where g = \x -> if x >= a then Just (x, pred x) else Nothing

{- Task -}
data Color = Red | Green | Blue

instance Show Color where
  show Red = "Red"
  show Green = "Green"
  show Blue = "Blue"

{- Task -}
change :: (Ord a, Num a) => a -> [[a]]
change 0 = [[]]
change amount = [c:cs | c <- [2,3,7], amount >= c, cs <- change (amount - c) ]

{- Task -}
data Point = Point Double Double

origin :: Point
origin = Point 0.0 0.0

distanceToOrigin :: Point -> Double
distanceToOrigin (Point x y) = sqrt (x ^ 2 + y ^ 2)

distance :: Point -> Point -> Double
distance (Point x1 y1) (Point x2 y2)= sqrt ((x2 - x1) ^ 2 + (y2 - y1) ^ 2)

{- Task -}
data Shape = Circle Double | Rectangle Double Double

area :: Shape -> Double
area (Circle r) = pi * r ^ 2
area (Rectangle a b) = a * b

{- Task -}
data Result' = Fail' Int | Success'

instance Show Result' where
    show (Fail' c) = "Fail: " ++ show c
    show Success' = "Success"

--doSomeWork' :: SomeData -> Result'
--doSomeWork' d = case doSomeWork d of
--                    (_, 0) -> Success'
--                    (_, c) -> Fail' c

{- Task -}
square :: Double -> Shape
square a = Rectangle a a

isSquare :: Shape -> Bool
isSquare (Rectangle a b) = a == b
isSquare _ = False

{- Constructor -}
data Person = Person String String Int
  deriving Show

firstName :: Person -> String
firstName (Person f _ _) = f

lastName :: Person -> String
lastName (Person _ l _) = l

age :: Person -> Int
age (Person _ _ a) = a

data Person' = Person' { firstName' :: String, lastName' :: String, age' :: Int}
  deriving (Show, Eq)

john = Person' "John" "Smith" 33
johnAge= john & age'

{- Task -}
timeToString :: UTCTime -> String
timeToString = formatTime defaultTimeLocale "%a %d %T"

data LogLevel = Error | Warning | Info

data LogEntry = LogEntry {timestamp :: UTCTime, logLevel :: LogLevel, message :: String}

logLevelToString :: LogLevel -> String
logLevelToString x = case x of
                        Error   -> "Error"
                        Warning -> "Warning"
                        _    -> "Info"

logEntryToString :: LogEntry -> String
logEntryToString l = timeToString (timestamp l) ++ ": " ++ logLevelToString (logLevel l)  ++ ": " ++ message l

{- Task -}
updateLastName :: Person' -> Person' -> Person'
updateLastName p1 p2 = p2 {lastName' = lastName' p1}

name :: Person' -> String
name person = firstName' person ++ " " ++ lastName' person

name' :: Person' -> String
name' (Person' fn ln _) = fn ++ " " ++ ln

name'' :: Person' -> String
name'' (Person' {lastName' = ln, firstName' = fn}) = fn ++ " " ++ ln

{- Task -}
abbrFirstName :: Person' -> Person'
abbrFirstName person@(Person' {firstName' = fn})
    | length fn > 2 = person {firstName' = (++ ".") . take 1 $ fn}
    | otherwise = person

{- Task -}
data Coord a = Coord a a
  deriving Show

getCenter :: Double -> Coord Int -> Coord Double
getCenter size (Coord x y) = Coord (center x) (center y)
  where center c = fromIntegral c * size + size / 2

getCell :: Double -> Coord Double -> Coord Int
getCell size (Coord x y)  = Coord (point x) (point y)
  where point c = floor $ c / size

roots :: Double -> Double -> Double -> Either [Char] (Double, Double)
roots a b c
  | discr >= 0 = Right (x1, x2)
  | otherwise = Left "Negative discriminant"
  where
      x1 = helper . negate $ d
      x2 = helper d
      helper x = (-b + x) / (2 * a)
      d = sqrt discr
      discr = b ^ 2 - 4 * a * c

{- Task -}
findDigit :: [Char] -> Maybe Char
findDigit xs
  | c /= "" = Just (head c)
  | otherwise = Nothing
  where
    c = filter isDigit xs

findDigit' :: [Char] -> Maybe Char
findDigit' [] = Nothing
findDigit' (x:xs) | isDigit x = Just x
                 | otherwise = findDigit' xs

{- Task -}
findDigitOrX :: [Char] -> Char
findDigitOrX xs = case findDigit xs of
                    Just x -> x
                    Nothing -> 'X'

{- Task -}
maybeToList :: Maybe a -> [a]
maybeToList (Just x) = [x]
maybeToList Nothing = []

maybeToList' :: Maybe a -> [a]
maybeToList' x = case x of
                  Just x -> [x]
                  _ -> []

listToMaybe :: [a] -> Maybe a
listToMaybe [] = Nothing
listToMaybe (x:[]) = Just x

listToMaybe' :: [a] -> Maybe a
listToMaybe' xs = case xs of
                    (x:_) -> Just x
                    _ -> Nothing

eitherToMaybe :: Either a a -> Maybe a
eitherToMaybe (Left a) = Just a
eitherToMaybe (Right _) = Nothing

{- strict and lazy calculations -}
data CoordLazy a = CoordLazy a a
  deriving Show

data CoordStrict a = CoordStrict !a !a
  deriving Show

getXLazy :: CoordLazy a -> a
getXLazy (CoordLazy x _) = x

getXStrict :: CoordStrict a -> a
getXStrict (CoordStrict x _) = x

-- getXLazy 3 undefined - Success
-- getXStrict 3 undefined - Fail

{-
import Data.Complex
import Data.Ratio

data Complex a = !a :+ !a
data Ratio a = !a :% !a
-}

{- Recursive data types -}

-- data [] a = [] | a : ([] a)
data List a = Nil | Cons a (List a)
  deriving Show

yz = Cons 'y' (Cons 'z' Nil)

{- Task -}
fromList :: List a -> [a]
fromList Nil = []
fromList (Cons x y) = x : fromList y

toList :: [a] -> List a
toList [] = Nil
toList (x:xs) = Cons x (toList xs)

{- Task -}
data Nat = Zero | Suc Nat
  deriving Show

nat1 = Suc (Suc Zero)
nat2 = Suc (Suc (Suc Zero))

fromNat :: Nat -> Integer
fromNat Zero = 0
fromNat (Suc n) = fromNat n + 1

toNat :: Integer -> Nat
toNat 0 = Zero
toNat n = Suc (toNat (n-1))

add :: Nat -> Nat -> Nat
add x y = toNat $ fromNat x + fromNat y

mul :: Nat -> Nat -> Nat
mul x y = toNat $ fromNat x * fromNat y

fac :: Nat -> Nat
fac x
  | n >= 0 = toNat $ helper 1 n
  | otherwise = error "Argumnet must be >= 0!"
  where
    n = fromNat x
    helper acc 0 = acc
    helper acc n = helper (acc * n) (n-1)

{- Task -}
data Tree a = Leaf a | Node (Tree a) (Tree a)
  deriving Show

height :: Tree a -> Int
height (Leaf _) = 0
height (Node n1 n2) = 1 + max (height n1) (height n2)

size :: Tree a -> Int
size (Leaf _)  = 1
size (Node n1 n2) = 1 + size n1 + size n2

{- Task -}
avg' :: Tree Int -> Int
avg' t =
    let (s,c) = go t
    in s `div` c
  where
    go :: Tree Int -> (Int,Int)
    go (Leaf n) = (n, 1)
    go (Node n1 n2) = (sum1 + sum2, count1 + count2)
      where
        (sum1, count1) = go n1
        (sum2, count2) = go n2

type AssociativeList k v = [(k, v)]

lookUp :: Eq k => k -> AssociativeList k v -> Maybe v
lookUp _ [] = Nothing
lookUp key ((x,y):xys)
  | key == x = Just y
  | otherwise = lookUp key xys

newtype IntList = IList [Int]
  deriving Show

newtype Identity a = Identity {runIdentity :: a}
  deriving (Eq, Ord, Show)

-- Monoid --
class Monoid a where
  mempty :: a             -- Neutral element
  mappend :: a -> a -> a  -- Operation that introduced on the elements of the monoid

  mconcat :: [a] -> a
  mconcat = foldr mappend mempty

{- Example with List -}
instance Monoid [a] where
  mempty = []
  mappend = (++)

{- Addition implementation for monoid -}
newtype Sum a = Sum {getSum :: a}
  deriving (Eq, Show, Read, Bounded, Ord)

instance Num a => Monoid (Sum a) where
  mempty = Sum 0
  (Sum x) `mappend` (Sum y) = Sum (x + y)

{- Product implementation for monoid -}
newtype Product a = Product {getProduct :: a}
  deriving (Eq, Show, Read, Bounded, Ord)

instance Num a => Monoid (Product a) where
  mempty = Product 1
  (Product x) `mappend` (Product y) = Product (x * y)

{- Task -}
xor :: Bool -> Bool -> Bool
x `xor` y = (x || y) && (not (x && y))

newtype Xor = Xor { getXor :: Bool }
  deriving (Eq, Show)

instance Monoid Xor where
  mempty = Xor False
  (Xor x) `mappend` (Xor y) = Xor (x `xor` y)

{- Monoid Pairs (tuple) -}
instance (Monoid a, Monoid b) => Monoid (a,b) where
  mempty = (mempty, mempty)
  (x1,x2) `mappend` (y1,y2) = (x1 `mappend` y1, x2 `mappend` y2)

{- Example: ("ABC", Product 1) `mappend` ("CDE", Product 3) -}

{- Monoid Maybe -}
instance Monoid a => Monoid (Maybe a) where
  mempty = Nothing
  Nothing `mappend` m = m
  m `mappend` Nothing = m
  Just m1 `mappend` Just m2 = Just (m1 `mappend` m2)

newtype First a = First { getFirst :: Maybe a }
  deriving (Eq, Ord, Read, Show)

instance Monoid (First a) where
  mempty = First Nothing
  First Nothing `mappend` r = r
  l `mappend` _             = l

{- mconcat $ map First [Nothing, Just 3, Just 5] -}
firstConcat = getFirst . mconcat . map First

-- MONAD --

{- Uncomment to check (Data.Functor) -}
--class Functor f where
--  fmap :: (a -> b) -> (f a -> f b)
--
--instance Functor [] where
--  fmap = map
--
--instance Functor Maybe where
--  fmap _ Nothing = Nothing
--  fmap f (Just a) = Just (f a)

{- Task -}
data Point3D a = Point3D a a a
  deriving Show

instance Functor Point3D where
  fmap f (Point3D x1 x2 x3)= Point3D (f x1) (f x2) (f x3)

{- Task -}
data GeomPrimitive a = Point' (Point3D a) | LineSegment (Point3D a) (Point3D a)

instance Functor GeomPrimitive where
  fmap f (Point' p) = Point' (fmap f p)
  fmap f (LineSegment p1 p2) = LineSegment (fmap f p1) (fmap f p2)

{- Binary trees -}
-- SOOO BEAUTIFULL !!!
data Tree' a = Leaf' a | Branch (Tree' a) a (Tree' a)
  deriving Show

testTree = Branch (Leaf' 2) 3 (Branch (Leaf' 4) 5 (Leaf' 6))

instance Functor Tree' where
  fmap g (Leaf' x)       = Leaf' (g x)
  fmap g (Branch l x r) = Branch (fmap g l) (g x) (fmap g r)

{- fmap == (<$>) -}
{- (^2) <$> testTree -}
{- (^2) `fmap` testTree -}

-- Operator (<$) replace all values in data structure

{- Task -}
data Tree'' a = Leaf'' (Maybe a) | Branch' (Tree'' a) (Maybe a) (Tree'' a)
  deriving Show

instance Functor Tree'' where
  fmap g (Leaf'' x)       = case x of
                                Just x  ->  Leaf'' (Just (g x))
                                _       ->  Leaf'' (Nothing)
  fmap g (Branch' l x r)  = case x of
                                Just x  ->  Branch' (fmap g l) (Just (g x)) (fmap g r)
                                _       ->  Branch' (fmap g l) (Nothing) (fmap g r)

--instance Functor ((,), s) where
--  fmap g (x,y) = (x, g y)

-- (a -> b) -> (,) s a -> (,) s b
-- (a -> b) -> (s,a) -> (s,b)
-- fmap succ (1, 'A')

--instance Functor (Either e) where
--  fmap _ (Left x) = Left x
--  fmap g (Right y) = Right (g y)

-- (a -> b) -> Either e a -> Either e b
-- fmap (+3) $ Right 1  ->  Right 4
-- fmap (+3) $ Left 2   ->  Left 2

{- Task -}
data Entry k1 k2 v = Entry (k1, k2) v
  deriving Show

data Map k1 k2 v = Map [Entry k1 k2 v]
  deriving Show

instance Functor (Entry k1 k2) where
  fmap g (Entry (k1, k2) v) = Entry (k1, k2) (g v)

instance Functor (Map k1 k2) where
  fmap g (Map xs) = Map (map (fmap g) xs)

{- Rules for Functors -}
-- 1. fmap id = id
-- 2. fmap (f . g) = fmap f . fmap g

{-
 -
 - Functions with effects
 -
 - f :: a -> Maybe b       -- can sometimes fail (Just b `or` Nothing)
 - f :: a -> [b]           -- can return many results
 - f :: a -> (Either s) b  -- can sometimes end with a typed exception
 - f :: a -> (s,b)         -- can write to the log
 - f :: a -> ((->) e) b    -- can read from external environments
 - f :: a -> (State s) b   -- work with a named state
 - f :: a -> IO a          -- carry out inputs/outputs (files, console ...)
 -
-}

{- Task -}
data Log a = Log [String] a
  deriving Show

add1Log = toLogger (+1) "added one"
mult2Log = toLogger (*2) "multiplied by 2"

toLogger :: (a -> b) -> String -> (a -> Log b)
toLogger f msg = Log [msg] . f

execLoggers :: a -> (a -> Log b) -> (b -> Log c) -> Log c
execLoggers x f g = Log (msg ++ msg') res'
  where
    Log msg res = f x
    Log msg' res' = g res

{-
 - class Monad m where
 -    return :: a -> m a
 -    (>>=) :: m a -> (a -> m b) -> m b   -- bind
 -    ...
 -
 - infixl 1 >>=
 - return True :: [] Bool       ~> [True]
 - return True :: Maybe Bool    ~> [Just True]
-}

toKleisli :: Monad m => (a -> b) -> (a -> m b)
toKleisli f x = return (f x)
-- toKleisli f = \x -> return (f x)
-- toKleisli f = return .

-- toKleisli cos 0 :: [Double]    ~> [1.0]
-- toKleisli cos 0 :: IO Double   ~> 1.0

{- Task -}
returnLog :: a -> Log a
returnLog x = Log [] x

{- Task -}
bindLog :: Log a -> (a -> Log b) -> Log b
bindLog (Log msg x) f =  Log (msg ++ msg') res
  where Log msg' res = f x

{- Task -}
instance Functor Log where
  fmap = liftM

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

instance Monad Log where
    return = returnLog
    (>>=) = bindLog

execLoggersList :: a -> [a -> Log a] -> Log a
execLoggersList x ys = foldl (>>=) (return x) ys

{-
 - class Monad m where
 -    return :: a -> m a
 -    (>>=) :: m a -> (a -> m b) -> m b   -- bind
 -    (>>) :: m a -> m b -> m b   -- bind
 -    x >> y = x >>= \_ -> y
 -    fail :: String -> m a
 -    fail s = error s
-}

-- Control.Monad
(=<<) :: Monad m => (a -> m b) -> m a -> m b
(=<<) = flip (>>=)

-- Fish
(<=<) :: Monad m => (b -> m c) -> (a -> m b) -> (a -> m c)
f <=< g = \x -> g x >>= f

{- Task -}
-- instance Functor SomeType where
--   fmap f x = x >>= (return . f)

-- `do` notation --

-- do { e1 ; e2 }           == e1 >> e2
-- do { p <- e1 ; e2 }      == e1 >>= \p -> e2
-- do { let v = e1 ; e2 }   == let v = e1 in do e2

-- Monad Identity
instance Monad Identity where
  return x = Identity x
  Identity x >>= k = k x

instance Functor Identity where
  fmap  f (Identity x) = Identity (f x)

instance Applicative Identity where
  pure x = Identity x
  Identity f <*> Identity v = Identity (f v)

wrap'n'succ :: Integer -> Identity Integer
wrap'n'succ x = Identity (succ x)

-- Monad Maybe
type Name = String
type Database = [(Name, Name)]

fathers, mothers :: Database
fathers = [("Bill", "John"),
           ("Ann", "John"),
           ("John", "Piter")]

mothers = [("Bill", "Jane"),
           ("Ann", "Jane"),
           ("John", "Alice"),
           ("Jane", "Dorothy"),
           ("Alice", "Mary")]

getM, getF :: Name -> Maybe Name
getM person = lookup person mothers
getF person = lookup person fathers

-- getF "Bill" >>= getM >>= getM
-- do {f <- getF "Bill"; gm <- getM f; getM gm}

granmas :: Name -> Maybe (Name, Name)
granmas person = do
  m <- getM person
  gmm <- getM m
  f <- getF person
  gmf <- getM f
  return (gmm, gmf)

{- Task -}
data Token = Number Int | Plus | Minus | LeftBrace | RightBrace
  deriving (Eq, Show)

asToken :: String -> Maybe Token
asToken x | all isDigit x = Just (Number (read x :: Int))
asToken x = case x of
                "+"       ->  Just Plus
                "-"       ->  Just Minus
                "("       ->  Just LeftBrace
                ")"       ->  Just RightBrace
                _         ->  Nothing

tokenize :: String -> Maybe [Token]
tokenize input = sequence (map asToken (words input))

-- Monad List

-- return 4 :: [Int] ~> [4]
-- [1,2] >>= (\x -> [x,x,x]) ~> [1,1,1,2,2,2]
-- concat $ map (\x -> [x,x,x]) [1,2] ~> [1,1,1,2,2,2]

list = [(x,y) | x <- [1,2,3], y <- [4,5,6]]

list' = do
  x <- [1,2,3]
  y <- [4,5,6]
  return (x,y)

list'' =
  [1,2,3] >>= (\x ->
  [4,5,6] >>= (\y ->
  return (x,y)))

lst = [(x,y) | x <- [1,2,3], y <- [1,2], x /= y]

{- Task -}
-- nextPositionsN :: Board -> Int -> (Board -> Bool) -> [Board]
-- nextPositionsN b 0 pred | pred b = [b]
-- nextPositionsN _ n _    | n <= 0  = []
-- nextPositionsN b n pred = do
--   p <- nextPositions b
--     nextPositionsN p (n-1) pred

lst' = do
  x <- [1,2,3]
  y <- [1,2]
  True <- return (x /= y)
  return (x,y)

lst'' =
  [1,2,3]         >>= (\x ->
  [1,2]           >>= (\y ->
  return (x /= y) >>= (\b ->
  case b of True -> return (x,y)
            _    -> fail "...")))

lst''' = do
  x <- [1,2,3]
  y <- [1,2]
  if x /= y then "E" else []
  return (x,y)

{- Task -}
pythagoreanTriple :: Int -> [(Int, Int, Int)]
pythagoreanTriple x = do
  a <- [1..x]
  b <- [1..x]
  c <- [1..x]
  if a^2 + b^2 == c^2 && a < b then "E" else []
  return (a,b,c)

-- Monad IO

-- type IO a = IO (RealWorld -> (RealWorld, a))
-- return :: a -> IO a ~> return a = \w -> (w,a)
-- (>>=) :: IO a -> (a -> IO b) -> IO b

-- instance Monad IO where
--  return a = \w -> (w,a)
--  (>>=) m k = \w -> case m w of (w',a) -> k a w'

getLine' :: IO String
getLine' = do
  c <- getChar
  if c == '\n' then
    return []
  else do
    cs <- getLine'
    return (c:cs)

putStr' :: String -> IO ()
putStr' [] = return ()
putStr' (x:xs) = putChar x >> putStr' xs

{-
import Control.Monad

sequence_ :: Monad m => [m a] -> m ()
sequence_ = foldr (>>) (return ())

-}

putStr'' :: String -> IO ()
putStr'' = sequence_ . map putChar

{-
mapM_ :: Monad m => (a -> m b) -> [a] -> m ()
mapM_ f = sequence_ , map f
-}

putStr''' :: String -> IO ()
putStr''' = mapM_ putChar

ext: hs

max_stars_repo_name: endygamedev/learn-haskell

max_issues_repo_name: endygamedev/learn-haskell

max_forks_repo_name: endygamedev/learn-haskell


text:
{-# LANGUAGE BangPatterns    #-}
{-# LANGUAGE RecordWildCards #-}
{-# LANGUAGE TupleSections   #-}

module Main where

import Data.List
import Data.Maybe
import Data.Ord
import Numeric.GSL.SimulatedAnnealing
import qualified Numeric.LinearAlgebra as H
import Text.Parsec hiding (parse)
import qualified Text.Parsec as Parsec

main :: IO ()
main = do
  !input <- parse <$> readFile "inputs/day-17.txt"

  putStrLn $ "Target: " <> show input

  putStrLn "\nPart 1:"
  -- print . hitsTarget input $! fancySolve input
  print $! bruteforce input

  putStrLn "\nPart 2:"
  print $! bruteforceBy length input

-- | The target area.
data Target = Target { left :: Int, right :: Int, bottom :: Int, top :: Int }
  deriving (Show)

-- | Find the velocity with the highest y-coordinate that passes through the
-- target.
--
-- This doesn't actually work
fancySolve :: Target -> (Int, Int)
fancySolve target =
  simanSolve 420 2 solverParams (20, 20)
    (solverPenalty target)
    (\(oldX, oldY) (newX, newY) -> fromIntegral $ abs (newX - oldX) + abs (newY - oldY))
    step
    (Just $ \velo -> show velo <> ": dist = " <> show (overshootDistance target velo))
  where
    solverParams = SimulatedAnnealingParams
      1000   -- Tries per step
      1000   -- Tries per temperature
      1      -- Maximum step size in random walk
      1.0    -- Boltzman constant for random walks
      2.0    -- Initial temperature
      1.69   -- Cooling rate
      0.0069 -- Final temperature

    step :: H.Vector Double -> Double -> (Int, Int) -> (Int, Int)
    step rands _ velo@(veloX, veloY) =
      ( veloX + round (((rands H.! 0) - 0.5) * overshootDistance target velo * 10)
      , veloY + round (((rands H.! 1) - 0.5) * overshootDistance target velo * 10)
      )

-- | Just exhaustively try all options to find starting velocity that results in
-- highest reached Y-coordinate.
bruteforce :: Target -> ((Int, Int), Int)
bruteforce = bruteforceBy (maximumBy (comparing snd))

bruteforceBy :: ([((Int, Int), Int)] -> a) -> Target -> a
bruteforceBy f target@Target{..}
  = f
  $ mapMaybe (\p -> (p,) <$> hitsTarget target p)
    [(x, y) | x <- [1 .. right], y <- [bottom .. abs bottom]]

-- | The penalty for the solver. This is negative the maximum Y-coordinate
-- unless the trajectory would not hit the target area, in which case it is the
-- distance to the center of the target area.
solverPenalty :: Target -> (Int, Int) -> Double
solverPenalty target initialVelo
  | Just maxYPos <- hitsTarget target initialVelo = fromIntegral (negate maxYPos)
  | otherwise                                     = (overshootDistance target initialVelo + 100) ** 2

-- | Whether the given velocity causes the target to be hit. Returns the highest
-- Y-position if it does.
hitsTarget :: Target -> (Int, Int) -> Maybe Int
hitsTarget Target{..} = go (0, 0)
  where
    go :: (Int, Int) -> (Int, Int) -> Maybe Int
    go (posX, posY) (veloX, veloY)
      | posX > right || posY < bottom = Nothing
      | otherwise =
          let newPos  = (posX + veloX, posY + veloY)
              newVelo = (veloX - signum veloX, veloY - 1)
           in if posX >= left && posX <= right && posY >= bottom && posY <= top
                then Just posY
                else max posY <$> go newPos newVelo

-- | The minimal distance the target's center is overshot by while following the
-- trajectory, bounded by double the target's bottom right corner.
overshootDistance :: Target -> (Int, Int) -> Double
overshootDistance Target{..} = lowestDistance (0, 0)
  where
    centerX = (left + right) `div` 2
    centerY = (top + bottom) `div` 2
    boundX = right * 2
    boundY = bottom * 2

    lowestDistance :: (Int, Int) -> (Int, Int) -> Double
    lowestDistance (posX, posY) (veloX, veloY) =
      let distance                  = sqrt $ fromIntegral (centerX - posX) ** 2
                                           + fromIntegral (centerY - posY) ** 2
          newPos@(newPosX, newPosY) = (posX + veloX, posY + veloY)
          newVelo                   = (veloX - signum veloX, veloY - 1)
       in if newPosX > boundX || newPosY < boundY
            then distance
            else min distance (lowestDistance newPos newVelo)

parse :: String -> Target
parse = fromRight' . Parsec.parse pTarget ""
  where
    fromRight' (Right x) = x
    fromRight' _         = error "This wasn't in our agreement!"

    pTarget :: Parsec String () Target
    pTarget = Target <$> (string "target area: x=" *> pInt)     <*> (string ".." *> pInt)
                     <*> (string ", y=" *> pInt <* string "..") <*> pInt

    pInt :: Parsec String () Int
    pInt = read <$> many (digit <|> char '-')


text:
{-# LANGUAGE BangPatterns      #-}
{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE Strict            #-}
{-# LANGUAGE StrictData        #-}
module Utils.BLAS where

import           Data.Complex
import           Data.Vector.Storable            as VS
import           Foreign.CUDA.BLAS               as BLAS
import           Foreign.CUDA.Driver             as CUDA
import           Foreign.Marshal
import           Numerical.HBLAS.BLAS.FFI
import           Numerical.HBLAS.BLAS.FFI.Level1
import           Numerical.HBLAS.BLAS.FFI.Level3

type BLASMMT e
   = Int -> Int -> Int -> VS.Vector e -> VS.Vector e -> IO (VS.Vector e)

type BLASDOTU e = VS.Vector e -> VS.Vector e -> IO e

class BLAS a where
  gemmBLAS :: BLASMMT a
  dotuBLAS :: BLASDOTU a

instance BLAS Double where
  {-# INLINE gemmBLAS #-}
  gemmBLAS m' n' k' a b = do
    let output = VS.replicate (m' * n') 0
        !m = fromIntegral m'
        !n = fromIntegral n'
        !k = fromIntegral k'
    unsafeWith a $ \aPtr ->
      unsafeWith b $ \bPtr ->
        unsafeWith output $ \cPtr ->
          cblas_dgemm_safe
            (encodeOrder BLASRowMajor)
            (encodeTranspose BlasNoTranspose)
            (encodeTranspose BlasNoTranspose)
            m
            n
            k
            1
            aPtr
            k
            bPtr
            n
            0
            cPtr
            n
    return output
  {-# INLINE dotuBLAS #-}
  dotuBLAS vec1 vec2 =
    unsafeWith vec1 $ \ptr1 ->
      unsafeWith vec2 $ \ptr2 ->
        cblas_ddot_safe (fromIntegral . VS.length $ vec1) ptr1 1 ptr2 1

instance BLAS (Complex Double) where
  {-# INLINE gemmBLAS #-}
  gemmBLAS m' n' k' a b = do
    let output = VS.replicate (m' * n') 0
        !m = fromIntegral m'
        !n = fromIntegral n'
        !k = fromIntegral k'
    with 1 $ \alphaPtr ->
      with 0 $ \betaPtr ->
        unsafeWith a $ \aPtr ->
          unsafeWith b $ \bPtr ->
            unsafeWith output $ \cPtr ->
              cblas_zgemm_safe
                (encodeOrder BLASRowMajor)
                (encodeTranspose BlasNoTranspose)
                (encodeTranspose BlasNoTranspose)
                m
                n
                k
                alphaPtr
                aPtr
                k
                bPtr
                n
                betaPtr
                cPtr
                n
    return output
  {-# INLINE dotuBLAS #-}
  dotuBLAS vec1 vec2 = do
    let output = VS.singleton 0
    unsafeWith vec1 $ \ptr1 ->
      unsafeWith vec2 $ \ptr2 ->
        unsafeWith output $ \outPtr ->
          cblas_zdotu_safe
            (fromIntegral . VS.length $ vec1)
            ptr1
            1
            ptr2
            1
            outPtr
    return . VS.head $ output

type CUBLASMMT a b e
   = Handle -> Int -> Int -> Int -> a e -> b e -> IO (VS.Vector e)

class CUBLAS a where
  gemmCuBLAS :: CUBLASMMT DevicePtr DevicePtr a

instance CUBLAS Float where
  {-# INLINE gemmCuBLAS #-}
  gemmCuBLAS handle !m !n !k matA matB = do
    let !sizeC = m * n
        output = VS.replicate sizeC 0
    with 1 $ \alpha ->
      with 0 $ \beta ->
        CUDA.allocaArray sizeC $ \matC ->
          unsafeWith output $ \outputPtr -> do
            sgemm handle N N n m k alpha matB n matA k beta matC n
            CUDA.peekArray sizeC matC outputPtr
    return output

instance CUBLAS Double where
  {-# INLINE gemmCuBLAS #-}
  gemmCuBLAS handle !m !n !k matA matB = do
    let !sizeC = m * n
        output = VS.replicate sizeC 0
    with 1 $ \alpha ->
      with 0 $ \beta ->
        CUDA.allocaArray sizeC $ \matC ->
          unsafeWith output $ \outputPtr -> do
            dgemm handle N N n m k alpha matB n matA k beta matC n
            CUDA.peekArray sizeC matC outputPtr
    return output

instance CUBLAS (Complex Float) where
  {-# INLINE gemmCuBLAS #-}
  gemmCuBLAS handle !m !n !k matA matB = do
    let !sizeC = m * n
        output = VS.replicate sizeC 0
    with 1 $ \alpha ->
      with 0 $ \beta ->
        CUDA.allocaArray sizeC $ \matC ->
          unsafeWith output $ \outputPtr -> do
            cgemm handle N N n m k alpha matB n matA k beta matC n
            CUDA.peekArray sizeC matC outputPtr
    return output

instance CUBLAS (Complex Double) where
  {-# INLINE gemmCuBLAS #-}
  gemmCuBLAS handle !m !n !k matA matB = do
    let !sizeC = m * n
        output = VS.replicate sizeC 0
    with 1 $ \alpha ->
      with 0 $ \beta ->
        CUDA.allocaArray sizeC $ \matC ->
          unsafeWith output $ \outputPtr -> do
            zgemm handle N N n m k alpha matB n matA k beta matC n
            CUDA.peekArray sizeC matC outputPtr
    return output

{-# INLINE unsafeWithGPU #-}
unsafeWithGPU :: (Storable e) => VS.Vector e -> (DevicePtr e -> IO b) -> IO b
unsafeWithGPU vec f = do
  let !len = VS.length vec
  CUDA.allocaArray len $ \devPtr -> do
    unsafeWith vec $ \ptr -> CUDA.pokeArray len ptr devPtr
    f devPtr

{-# INLINE gemmCuBLAS10 #-}
gemmCuBLAS10 :: (CUBLAS e, Storable e) => CUBLASMMT VS.Vector DevicePtr e
gemmCuBLAS10 handle m n k a matB =
  unsafeWithGPU a $ \matA -> gemmCuBLAS handle m n k matA matB

{-# INLINE gemmCuBLAS01 #-}
gemmCuBLAS01 :: (CUBLAS e, Storable e) => CUBLASMMT DevicePtr VS.Vector e
gemmCuBLAS01 handle m n k matA b =
  unsafeWithGPU b $ \matB -> gemmCuBLAS handle m n k matA matB

{-# INLINE gemmCuBLAS11 #-}
gemmCuBLAS11 :: (CUBLAS e, Storable e) => CUBLASMMT VS.Vector VS.Vector e
gemmCuBLAS11 handle m n k a b =
  unsafeWithGPU a $ \matA ->
    unsafeWithGPU b $ \matB -> gemmCuBLAS handle m n k matA matB

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 LispTypes where

import Data.Ratio
import Data.Complex
import Data.Array
import Data.IORef

import Numeric
import Text.ParserCombinators.Parsec.Error
import Control.Monad.Except
import System.IO

-- Env is a stateful variable (IORef) that holds a map of strings -> IORef LispVal
type Env = IORef[(String, IORef LispVal)]

-- |Lisp Atom data type
data LispVal = Atom String -- Simple Types
    | Number Integer
    | Float Double
    | String String
    | Character Char
    | Bool Bool
    | Ratio Rational -- Composite types
    | Complex (Complex Double)
    | List [LispVal]
    | DottedList [LispVal] LispVal
    | Vector (Array Int LispVal)
    -- Stores a primitive function
    | PrimitiveFunc ([LispVal] -> ThrowsError LispVal)
    -- Generalized function type
    | Func {
        params :: [String], -- Parameters name
        vararg :: Maybe String, -- name of a variable-length argument list
        body :: [LispVal], -- list of expressions
        closure :: Env -- the environment where the function was created
        }
    | IOFunc ([LispVal] -> IOThrowsError LispVal) -- A dirty function that performs IO
    | Port Handle -- Represents input and output devices

-- Defining Equality methods for LispVal
instance Eq LispVal where
    (==) (Number x)         (Number y)  = x == y
    (==) (Float x)          (Float y)   = x == y
    (==) (String x)         (String y)  = x == y
    (==) (Bool x)           (Bool y)    = x == y
    (==) (Ratio x)          (Ratio y)   = x == y
    (==) (Complex x)        (Complex y) = x == y
    (==) (DottedList xs x)  (DottedList ys y) =
        (List $ xs ++ [x]) == (List $ ys ++ [y])
    (==) l1@(List x) l2@(List y) =
        (length x == length y) && all (uncurry (==)) (zip x y)

instance Ord LispVal where
    compare (Number x)      (Number y)  = compare x y
    compare (Float x)       (Float y)   = compare x y
    compare (String x)      (String y)  = compare x y
    compare (Bool x)        (Bool y)    = compare x y
    compare (Ratio x)       (Ratio y)   = compare x y
    -- There is no compare instance for complex numbers since it violates trichotomy
    -- https://math.stackexchange.com/questions/257184/defining-the-complex-numbers/257208#257208
    compare (Complex x)     (Complex y) = error "Complex numbers do not satisfy the axiom of trichotomy"

    -- |Helper function to show a Double in full precision
showFullPrecision :: Double -> String
showFullPrecision x = showFFloat Nothing x ""

-- |Print the content of a LispVal
-- Pattern matching is used to destructure
-- The algebraic data type selecting a clause
-- Based on constructors
showVal :: LispVal -> String
-- Simple Types
showVal (Atom name) = name
showVal (String s) = "\"" ++ s ++ "\""
showVal (Number n) = show n
showVal (Float f) = showFullPrecision f
showVal (Character '\n') = ""
showVal (Character c) = "#\\" ++ [c]
showVal (Bool True) = "#t"
showVal (Bool False) = "#f"

-- Composite types
showVal (Ratio r) =
    (show . numerator) r ++ "/" ++ (show . denominator) r
showVal (Complex c) = (showFullPrecision. realPart) c
    ++ (if imagPart c >= 0 then "+" else "")
    ++ (showFullPrecision . imagPart) c ++ "i"
showVal (List l) = "(" ++ unwordsList l ++ ")"
showVal (DottedList hd tl) = "(" ++ unwordsList hd ++ " . " ++ showVal tl ++ ")"
showVal v@(Vector x) = "#(" ++ unwordsVector v ++ ")"

-- Show for functions and ports
showVal (PrimitiveFunc _) = "<primitive>"
showVal Func { params = args, vararg = varargs, body = body, closure = env} =
    "(lambda (" ++ unwords args ++
    (case varargs of
        Nothing -> ""
        Just arg -> " . " ++ arg) ++ ") ...)"
showVal (IOFunc _) = "<IO primitive>"
showVal (Port _) = "<IO port>"
-- |Helper function mapping showVal over a Lisp List
-- |Basically unwords for LispVal

{-
 unwordsList is defined in point-free style
 no argument is specified and it is written only in
 terms of function composition and partial application:
 map is partially applied to showVal (creates a fun that
 takes a list of LispVals and returns a list of their
 string representation). Then, it is composed to
 Haskell's unwords which merges strings in a list
 into a single string separated by spaces.

 This is an important example of currying
-}
unwordsList :: [LispVal] -> String
unwordsList = unwords . map showVal

unwordsVector :: LispVal -> String
unwordsVector (Vector y) = unwords . map showVal $ elems y

-- |Defining the show method for LispVal
-- See https://www.haskell.org/tutorial/classes.html
-- for details on typeclasses
instance Show LispVal where show = showVal

-- ERRORS

{-
    Some interpreted languages like PHP silently assign
    default values like #f or 0 as a result when errors occur
    This approach means that errors pass silently throughout the
    programs until they can become problematic and need long
    debugging sessions to get rid of the bug completely.

    Here, errors are signaled as soon as they happen and
    immediately break out of execution
-}

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

-- Make LispError an instance of Show
showError :: LispError -> String
showError (Default msg) = show msg
showError DivideByZero = show "Division by zero!"
showError (TypeMismatch expected found) = "Invalid type: expected "
    ++ expected ++ ", found " ++ show found
showError (Parser parseErr) = "Parse error at " ++ show parseErr
showError (BadSpecialForm msg form) = msg ++ ": " ++ show form
showError (NotFunction msg func) = msg ++ ": " ++ show func
showError (UnboundVar msg varname) = msg ++ ": " ++ varname
showError (NumArgs expected found) = "Expected " ++ show expected
    ++ " args: found values " ++ unwordsList found

instance Show LispError where show = showError

-- |Type to represent functions that may throw a LispError or return a value.

-- Either is only partially applied (curried) to LispError so that ThrowsError
-- can be used with any data type
-- Either is a monad and the "extra information" carried by it is whether or not
-- An error has occured. This is why haskell does not need a separate
-- Control-flow construct to handle exceptions.
-- Either also provides throwError and catchError
type ThrowsError = Either LispError

-- Needs an Error Monad to handle errors like unbound variables
-- ErrorT is a monad transformer that layers error handling on top of IO
-- May contain IO actions and throws a LispError
-- Curried type constructor that still accepts an argument: return type of the function
type IOThrowsError = ExceptT LispError IO

-- |Convert errors to their string representation and return it
-- trapError action
-- The result of calling trapError is an Either action which will
-- always have valid (Right) data.
trapError action = catchError action (return . show)

-- |Extract a value from a ThrowsError type returned by trapError
-- extractValue (Right val) = val
-- extractValue is purposedly left undefined for Left values
extractValue :: ThrowsError a -> a
extractValue (Right val) = val

max_stars_repo_name: zfnmxt/yasih

max_issues_repo_name: zfnmxt/yasih

max_forks_repo_name: zfnmxt/yasih


text:
module ComposableSDR.Sink
  ( fileSink
  , stdOutSink
  , audioFileSink
  , constellationPlotSink
  , AudioFormat(..)
  ) where

import           System.FilePath.Posix               (takeBaseName)
import           System.IO                           (stdout)

import           Data.Complex                        (imagPart, realPart)

import           Text.Printf

import qualified Control.Monad.Catch                 as MC

import qualified Sound.File.Sndfile                  as SF

import qualified Streamly.Internal.Data.Fold         as F
import qualified Streamly.Internal.Data.Fold.Types   as FL
import qualified Streamly.Internal.FileSystem.File   as FS
import qualified Streamly.Internal.FileSystem.Handle as FH
import qualified Streamly.Memory.Array               as A

import           ComposableSDR.Common
import           ComposableSDR.Types

{-# INLINABLE fileSink #-}
fileSink ::
     (MonadIO m, MC.MonadCatch m, Storable a)
  => FilePath
  -> FL.Fold m (A.Array a) ()
fileSink = FS.writeChunks

{-# INLINABLE stdOutSink #-}
stdOutSink ::
     (MonadIO m, MC.MonadCatch m, Storable a) => FL.Fold m (A.Array a) ()
stdOutSink = FH.writeChunks stdout

openAudioFile :: AudioFormat -> FilePath -> Int -> Int -> Int -> IO SF.Handle
openAudioFile fmt fp sr sn nch = SF.openFile (fp ++ ext) SF.WriteMode info
  where
    sfFmt AU  = (SF.HeaderFormatAu, ".au")
    sfFmt WAV = (SF.HeaderFormatWav, ".wav")
    (sff, ext) = sfFmt fmt
    info =
      SF.Info
        sn
        sr
        nch
        (SF.Format sff SF.SampleFormatFloat SF.EndianBig)
        1
        False

writeToAudioFile :: (MonadIO m) => SF.Handle -> A.Array Float -> m ()
writeToAudioFile h a = do
  _ <- liftIO $ SF.hPutBuffer h (toArray a)
  return ()

closeAudioFile :: SF.Handle -> IO ()
closeAudioFile = SF.hClose

audioFileSink ::
     (MC.MonadCatch m, MonadIO m)
  => AudioFormat
  -> Int
  -> Int
  -> Int
  -> String
  -> FL.Fold m (A.Array Float) ()
audioFileSink fmt sr sn nch fp =
  let wav = openAudioFile fmt fp sr sn nch
   in bracketIO wav closeAudioFile (F.drainBy . writeToAudioFile)

constellationPlotSink ::
     (MC.MonadCatch m, MonadIO m)
  => FilePath
  -> FL.Fold m (A.Array SamplesIQCF32) ()
constellationPlotSink fp =
  let writeHeader = writeFile fp $ unlines ["clear all; close all;", "v = [];"]
      toFloat (CFloat f) = f
      writeBody as =
        liftIO $
        appendFile fp $
        unlines $
          (\s ->
             printf
               "v(end+1) = %12.4e + j*%12.4e;"
               (toFloat $ realPart s)
               (toFloat $ imagPart s)) <$>
        A.toList as
      writeFooter =
        appendFile fp $
        unlines
          [ "n = length(v);"
          , "figure('color','white','position',[100 100 1200 400]);"
          , "plot(real(v), imag(v), 'x', 'Color',[0 0.2 0.4]);"
          , "xlabel('In-Phase');"
          , "ylabel('Quadrature');"
          , "grid on;"
          , "print -dpng -color \"-S1200,600\" " ++ takeBaseName fp ++ ".png"
          ]
   in bracketIO writeHeader (const writeFooter) (F.drainBy . const writeBody)

max_stars_repo_name: mryndzionek/composable-sdr

max_issues_repo_name: mryndzionek/composable-sdr

max_forks_repo_name: mryndzionek/composable-sdr


text:
module SampleWeight
    ( Weight
    , Bias
    , SampleWeight
    , createBinarySW
    , loadSW
    ) where

import Control.Applicative ((<$>))
import Control.Monad
import Numeric.LinearAlgebra
import qualified Data.ByteString.Lazy as BL
import qualified Codec.Compression.GZip as GZ (compress, decompress)
import Data.Binary (encode, decode)
import Text.Parsec

type Weight = Matrix R
type Bias = Vector R
type SampleWeight = ([Weight], [Bias])

assetsDir = "assets"
weightFiles = [ "sample-weight-w1"
              , "sample-weight-w2"
              , "sample-weight-w3"
              , "sample-weight-b1"
              , "sample-weight-b2"
              , "sample-weight-b3"
              ]

generatePath :: String -> String
generatePath p = assetsDir ++ "/" ++ p

createBinary :: String -> IO ()
createBinary p = do
    let bp = generatePath p

    ws <- readFile $ bp ++ ".csv"
    case parseCSV ws of
      Left e -> print e
      Right w -> do
        putStrLn $ "Creating binary Matrix file: " ++ bp
        let wm = fromLists $ fmap (read :: String -> Double) <$> w
        createPickle (bp ++ ".dat") wm
        putStrLn "Done"

createBinarySW :: IO ()
createBinarySW = forM_ weightFiles createBinary

csvStruct = endBy line eol
line = sepBy cell $ char ','
cell = many $ noneOf ",\n"
eol = char '\n'

parseCSV :: String -> Either ParseError [[String]]
parseCSV = parse csvStruct "* ParseError *"

createPickle :: String -> Matrix R -> IO ()
createPickle p w = BL.writeFile p $ (GZ.compress . encode) w

loadPickle :: String -> IO (Matrix R)
loadPickle p = do
    esw <- BL.readFile $ generatePath p ++ ".dat"
    return $ (decode . GZ.decompress) esw

loadSW :: IO SampleWeight
loadSW = do
    sw <- forM weightFiles loadPickle
    let (w,b) = splitAt 3 sw
    return (w, fmap flatten b)


text:
{-# LANGUAGE BangPatterns, DeriveDataTypeable, DeriveGeneric #-}
-- |
-- Module    : Statistics.Distribution.Normal
-- Copyright : (c) 2009 Bryan O'Sullivan
-- License   : BSD3
--
-- Maintainer  : bos@serpentine.com
-- Stability   : experimental
-- Portability : portable
--
-- The normal distribution.  This is a continuous probability
-- distribution that describes data that cluster around a mean.

module Statistics.Distribution.Normal
    (
      NormalDistribution
    -- * Constructors
    , normalDistr
    , normalFromSample
    , standard
    ) where

import Data.Aeson (FromJSON, ToJSON)
import Control.Applicative ((<$>), (<*>))
import Data.Binary (Binary)
import Data.Binary (put, get)
import Data.Data (Data, Typeable)
import GHC.Generics (Generic)
import Numeric.MathFunctions.Constants (m_sqrt_2, m_sqrt_2_pi)
import Numeric.SpecFunctions (erfc, invErfc)
import qualified Statistics.Distribution as D
import qualified Statistics.Sample as S
import qualified System.Random.MWC.Distributions as MWC

-- | The normal distribution.
data NormalDistribution = ND {
      mean       :: {-# UNPACK #-} !Double
    , stdDev     :: {-# UNPACK #-} !Double
    , ndPdfDenom :: {-# UNPACK #-} !Double
    , ndCdfDenom :: {-# UNPACK #-} !Double
    } deriving (Eq, Read, Show, Typeable, Data, Generic)

instance FromJSON NormalDistribution
instance ToJSON NormalDistribution

instance Binary NormalDistribution where
    put (ND w x y z) = put w >> put x >> put y >> put z
    get = ND <$> get <*> get <*> get <*> get

instance D.Distribution NormalDistribution where
    cumulative      = cumulative
    complCumulative = complCumulative

instance D.ContDistr NormalDistribution where
    logDensity = logDensity
    quantile   = quantile

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

instance D.Mean NormalDistribution where
    mean = mean

instance D.MaybeVariance NormalDistribution where
    maybeStdDev   = Just . D.stdDev
    maybeVariance = Just . D.variance

instance D.Variance NormalDistribution where
    stdDev = stdDev

instance D.Entropy NormalDistribution where
  entropy d = 0.5 * log (2 * pi * exp 1 * D.variance d)

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

instance D.ContGen NormalDistribution where
    genContVar d = MWC.normal (mean d) (stdDev d)

-- | Standard normal distribution with mean equal to 0 and variance equal to 1
standard :: NormalDistribution
standard = ND { mean       = 0.0
              , stdDev     = 1.0
              , ndPdfDenom = log m_sqrt_2_pi
              , ndCdfDenom = m_sqrt_2
              }

-- | Create normal distribution from parameters.
--
-- IMPORTANT: prior to 0.10 release second parameter was variance not
-- standard deviation.
normalDistr :: Double            -- ^ Mean of distribution
            -> Double            -- ^ Standard deviation of distribution
            -> NormalDistribution
normalDistr m sd
  | sd > 0    = ND { mean       = m
                   , stdDev     = sd
                   , ndPdfDenom = log $ m_sqrt_2_pi * sd
                   , ndCdfDenom = m_sqrt_2 * sd
                   }
  | otherwise =
    error $ "Statistics.Distribution.Normal.normalDistr: standard deviation must be positive. Got " ++ show sd

-- | Create distribution using parameters estimated from
--   sample. Variance is estimated using maximum likelihood method
--   (biased estimation).
normalFromSample :: S.Sample -> NormalDistribution
normalFromSample xs
  = normalDistr m (sqrt v)
  where
    (m,v) = S.meanVariance xs

logDensity :: NormalDistribution -> Double -> Double
logDensity d x = (-xm * xm / (2 * sd * sd)) - ndPdfDenom d
    where xm = x - mean d
          sd = stdDev d

cumulative :: NormalDistribution -> Double -> Double
cumulative d x = erfc ((mean d - x) / ndCdfDenom d) / 2

complCumulative :: NormalDistribution -> Double -> Double
complCumulative d x = erfc ((x - mean d) / ndCdfDenom d) / 2

quantile :: NormalDistribution -> Double -> Double
quantile d p
  | p == 0         = -inf
  | p == 1         = inf
  | p == 0.5       = mean d
  | p > 0 && p < 1 = x * ndCdfDenom d + mean d
  | otherwise      =
    error $ "Statistics.Distribution.Normal.quantile: p must be in [0,1] range. Got: "++show p
  where x          = - invErfc (2 * p)
        inf        = 1/0

max_stars_repo_name: StefanHubner/statistics

max_issues_repo_name: StefanHubner/statistics

max_forks_repo_name: StefanHubner/statistics


text:
{-# LANGUAGE DeriveDataTypeable, DeriveGeneric, FlexibleContexts,
    FlexibleInstances, UndecidableInstances #-}
-- |
-- Module    : Statistics.Distribution.Transform
-- Copyright : (c) 2013 John McDonnell;
-- License   : BSD3
--
-- Maintainer  : bos@serpentine.com
-- Stability   : experimental
-- Portability : portable
--
-- Transformations over distributions
module Statistics.Distribution.Transform
    (
      LinearTransform (..)
    , linTransFixedPoint
    , scaleAround
    ) where

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

-- | Linear transformation applied to distribution.
--
-- > LinearTransform μ σ _
-- > x' = μ + σ·x
data LinearTransform d = LinearTransform
  { linTransLocation :: {-# UNPACK #-} !Double
    -- ^ Location parameter.
  , linTransScale    :: {-# UNPACK #-} !Double
    -- ^ Scale parameter.
  , linTransDistr    :: d
    -- ^ Distribution being transformed.
  } deriving (Eq, Show, Read, Typeable, Data, Generic)

instance (FromJSON d) => FromJSON (LinearTransform d)
instance (ToJSON d) => ToJSON (LinearTransform d)

instance (Binary d) => Binary (LinearTransform d) where
    get = LinearTransform <$> get <*> get <*> get
    put (LinearTransform x y z) = put x >> put y >> put z

-- | Apply linear transformation to distribution.
scaleAround :: Double           -- ^ Fixed point
            -> Double           -- ^ Scale parameter
            -> d                -- ^ Distribution
            -> LinearTransform d
scaleAround x0 sc = LinearTransform (x0 * (1 - sc)) sc

-- | Get fixed point of linear transformation
linTransFixedPoint :: LinearTransform d -> Double
linTransFixedPoint (LinearTransform loc sc _) = loc / (1 - sc)

instance Functor LinearTransform where
  fmap f (LinearTransform loc sc dist) = LinearTransform loc sc (f dist)

instance D.Distribution d => D.Distribution (LinearTransform d) where
  cumulative (LinearTransform loc sc dist) x = D.cumulative dist $ (x-loc) / sc

instance D.ContDistr d => D.ContDistr (LinearTransform d) where
  density    (LinearTransform loc sc dist) x = D.density    dist ((x-loc) / sc) / sc
  logDensity (LinearTransform loc sc dist) x = D.logDensity dist ((x-loc) / sc) - log sc
  quantile      (LinearTransform loc sc dist) p = loc + sc * D.quantile      dist p
  complQuantile (LinearTransform loc sc dist) p = loc + sc * D.complQuantile dist p

instance D.MaybeMean d => D.MaybeMean (LinearTransform d) where
  maybeMean (LinearTransform loc _ dist) = (+loc) <$> D.maybeMean dist

instance (D.Mean d) => D.Mean (LinearTransform d) where
  mean (LinearTransform loc _ dist) = loc + D.mean dist

instance D.MaybeVariance  d => D.MaybeVariance (LinearTransform d) where
  maybeVariance (LinearTransform _ sc dist) = (*(sc*sc)) <$> D.maybeVariance dist
  maybeStdDev   (LinearTransform _ sc dist) = (*sc)      <$> D.maybeStdDev dist

instance (D.Variance d) => D.Variance (LinearTransform d) where
  variance (LinearTransform _ sc dist) = sc * sc * D.variance dist
  stdDev   (LinearTransform _ sc dist) = sc * D.stdDev dist

instance (D.MaybeEntropy d) => D.MaybeEntropy (LinearTransform d) where
  maybeEntropy (LinearTransform _ _ dist) = D.maybeEntropy dist

instance (D.Entropy d) => D.Entropy (LinearTransform d) where
  entropy (LinearTransform _ _ dist) = D.entropy dist

instance D.ContGen d => D.ContGen (LinearTransform d) where
  genContVar (LinearTransform loc sc d) g = do
    x <- D.genContVar d g
    return $! loc + sc * x


text:
{-# LANGUAGE Strict, StrictData #-}
{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE DeriveGeneric, DeriveAnyClass #-}
{-# LANGUAGE TypeSynonymInstances #-}
{-# LANGUAGE FlexibleContexts #-}
{-# LANGUAGE StandaloneDeriving #-}
{-# LANGUAGE UndecidableInstances #-}
{-# LANGUAGE TypeFamilies, TypeFamilyDependencies #-}

------------------------------------------------
-- |
-- Module    :  Numeric.Algebra.Homology
-- Copyright :  (c) Jun Yoshida 2019
-- License   :  BSD3
--
-- Chain complexes over integers
--
------------------------------------------------

module Numeric.Algebra.Homology
  (ChainEliminable, Homology(..), homology) where

import GHC.Generics (Generic, Generic1)
import Control.DeepSeq (NFData, force, deepseq)

import Control.Parallel.Strategies

import Control.Monad
import Control.Monad.ST

import Control.Applicative

import Data.STRef

import Data.Bifunctor
import Data.Traversable
import Data.Tuple (swap)
import Data.Maybe
import Data.List

import qualified Numeric.LinearAlgebra as LA
import qualified Numeric.LinearAlgebra.Devel as LA

import qualified Data.Vector as V
import qualified Data.Vector.Mutable as MV

import Numeric.Algebra.FreeModule
import Numeric.Algebra.Presentation
import Numeric.Matrix.Integral
import Numeric.Matrix.Integral.NormalForms
import Numeric.F2 (F2)
import qualified Numeric.Matrix.F2 as F2
import qualified Numeric.Matrix.F2.Mutable as F2

{-- for Debug
import Debug.Trace

traceCond :: Bool -> String -> a -> a
traceCond False _  = id
traceCond True msg = trace msg
--}

--------------
-- * Classes
--------------
class (Normalized a) => ChainEliminable a where
  elimChainHead :: Mat a -> Mat a -> (NormalForm a, Mat a)
  -- ^ If
  -- @ let (a',v) = elimChainHead u a @
  -- Then
  -- @ a' == toNF (a <> u^{-1}) @
  -- and the vectors of v (i.e. toVecs v) induces a basis for the cokernel of a⊗ℚ.

-------------------------------------
-- * Computation of homology groups
-------------------------------------
data Homology a = Homology {
  freeCycs :: [Vec a],
  torsions :: [(a,Vec a)],
  bndries :: [Vec a]
  } deriving (Generic)

deriving instance (Coefficient a, Show a, Show (Vec a)) => Show (Homology a)
deriving instance (Coefficient a) => NFData (Homology a)

-- | Compute the kernel and the image of a matrix from its normal form.
homologyOfNF :: (Normalized a) => NormalForm a -> (Homology a, Homology a)
homologyOfNF nf
  = let invs = invariants nf ++ repeat 0
        cycs = map snd $ filter ((==0) . fst) $ zip invs (toVecs (basisDom nf))
        img = filter ((/=0) . fst) $ zip invs (toVecs (basisCod nf))
        (tors,bnds') = partition ((/=1) . fst) img
        bnds = map (uncurry scale) bnds'
    in (Homology [] tors bnds, Homology cycs [] [])

-- | Decompose complexes into two-term complexes with differential in a Smith normal form.
chainToNormals :: (ChainEliminable a, Traversable t) => t (Mat a) -> (Mat a,t (NormalForm a))
chainToNormals = first (maybe (ident 0) id) . mapAccumL chainToNormals' Nothing
  where
    chainToNormals' mayU mat
      = let (u,mat') = case mayU of
                         Just u -> (u, mat <> u)
                         Nothing -> (ident (rankDom mat), mat)
            !(!matNF, !coker) = force (elimChainHead u mat')
        in (Just coker, matNF)

-- | Compute homology group with integral coefficient.
homology :: (ChainEliminable a, Traversable t, Alternative t, NFData (t (NormalForm a)), NFData (t (Homology a, Homology a))) => t (Mat a) -> t (Homology a)
homology diffs =
  let !(!coker, !nfs) = force (chainToNormals diffs)
      !hs = force (fmap homologyOfNF nfs `using` parTraversable (rparWith rdeepseq))
  in (uncurry ((<|>) . pure)) (mapAccumR pile (Homology (toVecs coker) [] []) hs)
  where
    pile (Homology cycs _ _) (thisH, nextCycH) =
      (nextCycH, thisH {freeCycs = cycs})

-----------------
-- * Instances
-----------------
-- | Enables chain-level Gaussian eliminations on integer matrices.
instance ChainEliminable LA.Z where
  elimChainHead basis matA = runST $ do
    stMatA <- LA.thawMatrix matA
    stBasis <- LA.thawMatrix basis
    stU <- LA.thawMatrix (LA.ident (LA.rows matA))
    -- Compute the Smith normal form
    smithRepST stU stMatA stBasis
    -- Write the new basis of the source
    resBasis <- LA.unsafeFreezeMatrix stBasis
    -- Take the invariant factor.
    invFs <- LA.takeDiag <$!> LA.unsafeFreezeMatrix stMatA
    let rkQ = length (takeWhile (/=0) (LA.toList invFs)) -- rank over the field of fractions; e.g. Q for Z
    -- Free closure of the image of matA & its orthogonal.
    fimage <- if LA.rows matA > 0 && rkQ > 0
              then LA.extractMatrix stU LA.AllRows (LA.ColRange 0 (rkQ-1))
              else return $ (LA.rows matA LA.>< 0) []
    fcoker <- if LA.rows matA > 0 && rkQ < LA.rows matA
              then LA.extractMatrix stU LA.AllRows (LA.FromCol rkQ)
              else return $ (LA.rows matA LA.>< 0) []
    return (SmithNF invFs resBasis fimage, fcoker)

instance ChainEliminable F2 where
  elimChainHead basis matA = runST $ do
    stMatA <- F2.thawMatrixF2 matA
    stBasis <- F2.thawMatrixF2 basis
    stU <- F2.newIdentF2 (F2.rows matA)
    -- Compute the Smith normal form
    rk <- F2.mdiagRepF2 stU stMatA stBasis
    -- Write the new basis of the source
    resBasis <- F2.unsafeFreezeMatrixF2 stBasis
    (fimage,fcoker) <- F2.splitColumnsAt rk <$> F2.unsafeFreezeMatrixF2 stU
    return (DiagNF rk resBasis fimage, fcoker)

max_stars_repo_name: Junology/linkedgram

max_issues_repo_name: Junology/linkedgram

max_forks_repo_name: Junology/linkedgram


text:
module CreateTerrain
where

import Data.Maybe
import Control.Monad.State
import System.Random

import Galaxy
import Math
import CreateGalaxy
import DataFunction
import Utils
import Statistics
import Good
import Terrain

nearsystems :: [String]
nearsystems = ["Alpha Centauri",
               "Barnard's Star",
               "Wolf 359",
               "Lalande 21185",
               "Sirius",
               "Lyuten 726-8",
               "Ross 154",
               "Ross 248",
               "Epsilon Eridani",
               "Lacaille 9352",
               "Ross 128",
               "EZ Aquarii",
               "Procyon",
               "61 Cygni",
               "Tau Ceti",
               "Fomalhaut",
               "Struve 2398",
               "Groombridge 34",
               "Epsilon Indi",
               "DX Cancri",
               "GJ 1061",
               "YZ Ceti",
               "Lyuten's Star",
               "Kapteyn's Star"
              ]

testGalaxy :: Galaxy Terrain
testGalaxy = testRandomGalaxy 20 16

testRandomGalaxy :: Int -> Int -> Galaxy Terrain
testRandomGalaxy v numsys =
  let r = mkStdGen v
  in evalState (createGalaxy (createTerrain stdGoods) "milky way" (take numsys $ nearsystems ++ map show [1..numsys])) r

createTerrain :: [Good] -> Planet () -> Rnd Terrain
createTerrain gs p =
  case planettype p of
    Planetoid         -> createRockyTerrain gs p
    NoAtmosphere      -> createRockyTerrain gs p
    RockyPlanet _     -> createRockyTerrain gs p
    SmallGasGiant     -> return (Terrain [] [])
    MediumGasGiant    -> return (Terrain [] [])
    LargeGasGiant     -> return (Terrain [] [])
    VeryLargeGasGiant -> return (Terrain [] [])

createRockyTerrain :: [Good] -> Planet () -> Rnd Terrain
createRockyTerrain gs p = do
  massmult <- randomRM (0, 1000 * planetMass p)
  gs' <- mapMaybeM (createNaturalGood massmult (planettype p)) gs
  return (Terrain gs' [])

createNaturalGood :: Flt -> PlanetType -> Good -> Rnd (Maybe (Resource, ResourceUnit))
createNaturalGood massmult pt g =
  case natural g of
    Nothing                       -> return Nothing
    Just (Natural atms _ initial) -> do
      let atmNeeded = not $ null atms
      let invAtm = if not atmNeeded
                     then False
                     else case pt of
                       RockyPlanet a -> a `notElem` atms
                       _             -> True
      if invAtm
        then return Nothing
        else do
          mult <- randomRM (0, initial)
          let v = floor $ mult * massmult
          return $ Just ((g, v), v)

max_stars_repo_name: anttisalonen/starrover

max_issues_repo_name: anttisalonen/starrover

max_forks_repo_name: anttisalonen/starrover


text:
-- root finding examples
import Numeric.GSL
import Numeric.LinearAlgebra
import Text.Printf(printf)

rosenbrock a b [x,y] = [ a*(1-x), b*(y-x^2) ]

test method = do
    print method
    let (s,p) = root method 1E-7 30 (rosenbrock 1 10) [-10,-5]
    print s -- solution
    disp p -- evolution of the algorithm

jacobian a b [x,y] = [ [-a    , 0]
                     , [-2*b*x, b] ]

testJ method = do
    print method
    let (s,p) = rootJ method 1E-7 30 (rosenbrock 1 10) (jacobian 1 10) [-10,-5]
    print s
    disp p

disp = putStrLn . format "  " (printf "%.3f")

main = do
    test Hybrids
    test Hybrid
    test DNewton
    test Broyden

    mapM_ testJ [HybridsJ .. GNewton]

max_stars_repo_name: curiousleo/liquidhaskell

max_issues_repo_name: curiousleo/liquidhaskell

max_forks_repo_name: curiousleo/liquidhaskell


text:
--
-- Algebra: algebra library
--

module CNN.Algebra (
  vadd
, vsub
, vsum
, vscale
, vdot
, madd
, msub
, msum
, mscale
, mmul
, mdot
, mavg
, transpose
) where

import Data.List (foldl')
import Numeric.LinearAlgebra
--import Numeric.LinearAlgebra.Data

-- VECTOR FUNCTIONS

vadd :: Vector R -> Vector R -> Vector R
vadd = (+)

vsub :: Vector R -> Vector R -> Vector R
vsub = (-)

{- |
vsum

>>> let vs = [vector [1.0,2.0,3.0], vector [4.0,5.0,6.0]]
>>> vsum vs
[5.0,7.0,9.0]

-}

vsum :: [Vector R] -> Vector R
vsum [] = fromList [0.0]
vsum (v:vs) = v + (vsum vs)

vscale :: Double -> Vector R -> Vector R
vscale s a = scale s a

vdot :: Vector R -> Vector R -> Double
vdot a b = a `dot` b

-- MATRIX FUNCTIONS

{- |
madd

>>> let m1 = fromLists [[1.0,10.0,100.0],[2.0,20.0,200.0]]
>>> let m2 = fromLists [[3.0,30.0,300.0],[4.0,40.0,400.0]]
>>> madd m1 m2
(2><3)
 [ 4.0, 40.0, 400.0
 , 6.0, 60.0, 600.0 ]

-}

madd :: Matrix R -> Matrix R -> Matrix R
madd a b = a + b

{- |
msub

>>> let m1 = fromLists [[3.0,30.0,300.0],[5.0,50.0,500.0]]
>>> let m2 = fromLists [[1.0,10.0,100.0],[2.0,20.0,200.0]]
>>> msub m1 m2
(2><3)
 [ 2.0, 20.0, 200.0
 , 3.0, 30.0, 300.0 ]

-}

msub :: Matrix R -> Matrix R -> Matrix R
msub a b = a - b

{- |
msum

>>> let m1 = fromLists [[3.0,30.0,300.0],[5.0,50.0,500.0]]
>>> let m2 = fromLists [[3.0,30.0,300.0],[4.0,40.0,400.0]]
>>> let m3 = fromLists [[1.0,10.0,100.0],[2.0,20.0,200.0]]
>>> msum [m1,m2,m3]
(2><3)
 [  7.0,  70.0,  700.0
 , 11.0, 110.0, 1100.0 ]

-}

msum :: [Matrix R] -> Matrix R
--msum = foldr madd []
msum = foldl' madd (matrix 1 [0.0])

{- |
mscale

>>> let m0 = fromLists [[0.0,1.0,2.0],[3.0,4.0,5.0]]
>>> mscale 2.0 m0
(2><3)
 [ 0.0, 2.0,  4.0
 , 6.0, 8.0, 10.0 ]

-}

mscale :: Double -> Matrix R -> Matrix R
mscale s m = scale s m

{- |
mmul

>>> let v1 = vector [3.0,30.0,300.0]
>>> let m1 = fromLists [[1.0,10.0,100.0],[2.0,20.0,200.0]]
>>> mmul v1 m1
[30303.0,60606.0]

-}

mmul :: Vector R -> Matrix R -> Vector R
mmul v m = m #> v

{- |
mdot

>>> let m1 = fromLists [[1.0,2.0],[3.0,4.0]]
>>> let m2 = fromLists [[5.0,6.0],[7.0,8.0]]
>>> mdot m1 m2
70.0

-}

mdot :: Matrix R -> Matrix R -> Double
--mdot a b = sum $ concat $ toLists (a * b)
mdot a b = flatten a <.> flatten b

{- |
mavg

>>> let m1 = fromLists [[1.0,2.0],[3.0,4.0]]
>>> let m2 = fromLists [[4.0,3.0],[2.0,1.0]]
>>> mavg [m1,m2]
(2><2)
 [ 2.5, 2.5
 , 2.5, 2.5 ]
>>> let m1 = fromLists [[1.0,2.0,3.0],[4.0,5.0,6.0]]
>>> let m2 = fromLists [[6.0,5.0,4.0],[3.0,2.0,1.0]]
>>> let m3 = fromLists [[2.0,2.0,2.0],[2.0,2.0,2.0]]
>>> mavg [m1,m2,m3]
(2><3)
 [ 3.0, 3.0, 3.0
 , 3.0, 3.0, 3.0 ]

-}

mavg :: [Matrix R] -> Matrix R
mavg ms = mscale a ss
  where
    a = 1.0 / fromIntegral (length ms)
    ss = msum ms

{- |
transpose

>>> let m1 = fromLists [[1.0, 2.0], [3.0, 4.0]] :: Matrix R
>>> tr m1
(2><2)
 [ 1.0, 3.0
 , 2.0, 4.0 ]
>>> let m2 = fromLists [[1.0, 2.0, 3.0], [4.0, 5.0, 6.0], [7.0, 8.0, 9.0]] :: Matrix R
>>> tr m2
(3><3)
 [ 1.0, 4.0, 7.0
 , 2.0, 5.0, 8.0
 , 3.0, 6.0, 9.0 ]
>>> let m3 = fromLists [[1.0, 2.0, 3.0, 10.0], [4.0, 5.0, 6.0, 11.0], [7.0, 8.0, 9.0, 12.0]] :: Matrix R
>>> tr m3
(4><3)
 [  1.0,  4.0,  7.0
 ,  2.0,  5.0,  8.0
 ,  3.0,  6.0,  9.0
 , 10.0, 11.0, 12.0 ]

-}

transpose :: Eq a => [[a]] -> [[a]]
transpose fs
  | null l    = []
  | otherwise = l : transpose fs'
  where
    (ls, fs') = tr fs
    l = concat ls
    tr :: [[a]] -> ([[a]], [[a]])
    tr fs = unzip $ map (splitAt 1) fs

max_stars_repo_name: eijian/deeplearning

max_issues_repo_name: eijian/deeplearning

max_forks_repo_name: eijian/deeplearning


text:
module Vec where
import Linear.Epsilon (nearZero)
import qualified Data.Map.Strict as Map
import Data.Complex

-- http://blog.sigfpe.com/2007/03/monads-vector-spaces-and-quantum.html

data W b a = W { runW :: [(a,b)] } deriving (Eq,Show,Ord)
instance Semigroup (W b a) where
  (W x) <> (W y) = W (x <> y)
instance Monoid (W b a) where
  mempty = W mempty
mapW f (W l) = W $ map (\(a,b) -> (a,f b)) l
instance Functor (W b) where
    fmap f (W a) = W $ map (\(a,p) -> (f a,p)) a

instance Num b => Applicative (W b) where
  pure x = W [(x,1)]
  (W fs) <*> (W xs) = W [(f x, a * b) | (f, a) <- fs, (x, b) <- xs]
  -- the fs is a a vector of rows kind of.
    -- error "I don't feel like it"

instance Num b => Monad (W b) where
   return x = W [(x,1)]
   l >>= f = W $ concatMap (\(W d,p) -> map (\(x,q)->(x,p*q)) d) (runW $ fmap f l)

a .* b = mapW (a*) b

instance (Eq a,Show a,Num b) => Num (W b a) where
     W a + W b = W $ (a ++ b)
     a - b = a + (-1) .* b
     _ * _ = error "Num is annoying"
     abs _ = error "Num is annoying"
     signum _ = error "Num is annoying"
     fromInteger a = if a==0 then W [] else error "fromInteger can only take zero argument"

collect :: (Ord a, Num b) => W b a -> W b a
collect = W . Map.toList . Map.fromListWith (+) . runW

trimZero = W . filter (\(k,v) -> not $ nearZero v) . runW
simplify :: Ord a => Q a -> Q a
simplify = trimZero . collect
-- filter (not . nearZero . snd)

type P a = W Double a

type Q a = W (Complex Double) a

star :: Num b => W (Complex b) a -> W (Complex b) a
star = mapW conjugate

kron :: Num b => W b a -> W b c -> W b (a,c)
kron (W x) (W y) = W [((a,c), r1 * r2) | (a,r1) <- x , (c,r2) <- y ]


text:
{-# LANGUAGE BangPatterns #-}
{-# LANGUAGE RankNTypes  #-}

module Main where

import qualified Data.HashTable.Class                 as C
import           Data.HashTable.IO
import           Data.HashTable.Test.Common
import qualified Data.Vector.Unboxed as V
import qualified Data.Vector.Unboxed.Mutable as VM
import           Statistics.Quantile (continuousBy, cadpw)
import           Statistics.Sample
import           System.Environment
import           System.Random.MWC

overhead :: C.HashTable h =>
            FixedTableType h ->
            GenIO ->
            IO Double
overhead dummy rng = do
    size <- uniformR (1000,50000) rng
    !v <- replicateM' size $ uniform rng
    let _ = v :: [(Int,Int)]

    !ht <- fromList v
    forceType dummy ht

    x <- computeOverhead ht
    return x

  where
    replicateM' :: Int -> IO a -> IO [a]
    replicateM' !sz f = go sz []
      where
        go !i !l | i == 0 = return l
                 | otherwise = do
                     !x <- f
                     go (i-1) (x:l)

-- Returns mean / stddev
runTrials :: C.HashTable h =>
             FixedTableType h
          -> GenIO
          -> Int
          -> IO (Double, Double, Double, Double)
runTrials dummy rng ntrials = do
    sample <- rep ntrials $ overhead dummy rng

    let (m, v) = meanVarianceUnb sample
    return (m, sqrt v, p95 sample, pMax sample)
  where
    p95 sample = continuousBy cadpw 19 20 sample

    pMax sample = V.foldl' max (-1) sample

    rep !n !f = do
        mv <- VM.new n
        go mv

      where
        go !mv = go' 0
          where
            go' !i | i >= n = V.unsafeFreeze mv
                   | otherwise = do
                !d <- f
                VM.unsafeWrite mv i d
                go' $ i+1

main :: IO ()
main = do
    rng <- do
        args <- getArgs
        if null args
          then withSystemRandom (\x -> (return x) :: IO GenIO)
          else initialize $ V.fromList [read $ head args]

    runTrials dummyLinearTable rng nTrials >>= report "linear hash table"
    runTrials dummyBasicTable rng nTrials >>= report "basic hash table"
    runTrials dummyCuckooTable rng nTrials >>= report "cuckoo hash table"

  where
    nTrials = 200

    report name md = putStrLn msg
      where msg = concat [ "\n(Mean,StdDev,95%,Max) for overhead of "
                         , name
                         , " ("
                         , show nTrials
                         , " trials): "
                         , show md
                         , "\n" ]

    dummyBasicTable = dummyTable
                      :: forall k v . BasicHashTable k v

    dummyLinearTable = dummyTable
                       :: forall k v . LinearHashTable k v

    dummyCuckooTable = dummyTable
                       :: forall k v . CuckooHashTable k v

max_stars_repo_name: cornell-pl/HsAdapton

max_issues_repo_name: cornell-pl/HsAdapton

max_forks_repo_name: cornell-pl/HsAdapton


text:
-- ch03 example 3.14

import qualified Statistics.Distribution as D
import Statistics.Distribution.Binomial

main = do
  putStrLn $ "Unfortunately, we have no negative binomial distribution codes"

max_stars_repo_name: raydsameshima/Stat

max_issues_repo_name: raydsameshima/Stat

max_forks_repo_name: raydsameshima/Stat


text:
module Main where

import Numeric.LinearAlgebra
import LeastSquares

main :: IO ()
main = do
  let a = Mat $ (3><2) [1.0, 1.0, 1.0, 4.0, 2.0, -1.0]
  let b = Vec $ vector [13, 27, 1]
  let x = Var "x" 2
  print $ minimize (SumSquares (a * x - b)) []
  -- let a1 = Mat $ (3><1) [1.0, 1.0, 1.0]
  -- let a2 = Mat $ (3><1) [2.0, 4.0, -1.0]
  -- let x1 = Var "x1" 1
  -- let x2 = Var "x2" 1
  print $ minimize (SumSquares (a * x - b) + 2 * SumSquares x) []
  let c = Mat $ (1><2) [1.0, -1]
  let d = Vec $ vector [2]
  print $ minimize (SumSquares (a * x - b)) [c * x :==: d]

max_stars_repo_name: tridao/leastsquares

max_issues_repo_name: tridao/leastsquares

max_forks_repo_name: tridao/leastsquares