text:
-- Copyright 2016 TensorFlow authors.
--
-- Licensed under the Apache License, Version 2.0 (the "License");
-- you may not use this file except in compliance with the License.
-- You may obtain a copy of the License at
--
--     http://www.apache.org/licenses/LICENSE-2.0
--
-- Unless required by applicable law or agreed to in writing, software
-- distributed under the License is distributed on an "AS IS" BASIS,
-- WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
-- See the License for the specific language governing permissions and
-- limitations under the License.

{-# LANGUAGE FlexibleContexts #-}
{-# LANGUAGE OverloadedStrings #-}
{-# LANGUAGE TypeFamilies #-}
-- | Rendering of TensorFlow operations as Haskell functions.

module TensorFlow.OpGen
  ( OpGenFlags(..)
  , docOpList
  , flagParser)
  where

import Prelude hiding (head, tail)

import Control.Monad (guard)
import Data.Char (toLower, toUpper)
import Data.Foldable (toList)
import Data.Maybe (fromMaybe, maybeToList)
import Data.ProtoLens (def, showMessage)
import Data.List.NonEmpty (NonEmpty((:|)), head)
import qualified Data.List.NonEmpty as NE
import Lens.Family2 ((^.), (.~), (&), view)
import Options.Applicative (Parser, help, long, strOption, value)
import Proto.Tensorflow.Core.Framework.OpDef
  ( OpList
  , OpDef
  , OpDef'ArgDef
  , attr
  , description
  , inputArg
  , name
  , numberAttr
  , op
  , outputArg
  , summary
  , type'
  , typeAttr
  )
import Proto.Tensorflow.Core.Framework.Types (DataType(..))
import System.FilePath (takeBaseName)
import TensorFlow.OpGen.AttrVal
  (AttrDef
  , AttrCase(..)
  , AttrTemplate(..)
  , Template
  , attrDef
  , attrOriginal
  , attrTemplate
  , templateDefault
  , templateRestrictions
  )
import Text.PrettyPrint.Mainland
  ( Doc
  , (<>)
  , (<+>)
  , (</>)
  , (<+/>)
  , brackets
  , comma
  , commasep
  , dquotes
  , empty
  , enclose
  , flatten
  , folddoc
  , hang
  , indent
  , int
  , parens
  , sep
  , stack
  , strictText
  , tuple
  )
import qualified Data.Map.Strict as Map
import qualified Data.Set as Set
import qualified Data.Text as Text
import qualified Data.Semigroup as Semigroup
import Data.Text (Text)

data OpGenFlags = OpGenFlags
     { outputFile :: String
     , prefix :: String
     , excludeList :: String
     }

flagParser :: Parser OpGenFlags
flagParser = OpGenFlags
     <$> strOption (mconcat [ long "output"
                            , help "File to write."
                            ])
     <*> strOption (mconcat [ long "prefix"
                            , help "Haskell package prefix to use"
                            ])
     <*> strOption (mconcat [ long "exclude_list"
                            , value ""
                            , help "Comma separated Ops names to ignore"
                            ])

docOpList :: OpGenFlags -> OpList -> Doc
docOpList flags opList =
  stack [ "{-# LANGUAGE ConstraintKinds #-}"
        , "{-# LANGUAGE DataKinds #-}"
        , "{-# LANGUAGE FlexibleInstances #-}"
        , "{-# LANGUAGE OverloadedStrings #-}"
        , "{-# LANGUAGE RankNTypes #-}"
        , "{-# LANGUAGE ScopedTypeVariables #-}"
        , "module" <+> strictText moduleName <+> "where"
        , empty
        , imports
        , empty
        , folddoc (\x y -> x </> empty </> y)
                  (map renderDef $
                   filter (not . flip elem exclusions . view name) $
                   toList $ opList ^. op)
        ]
  where moduleName =
            Text.pack (prefix flags) <> "." <> camelCase
             -- Discards the optional trailing _op_lib
            (fromMaybe shortName (Text.stripSuffix "_op_lib" shortName))
        shortName = Text.pack (takeBaseName $ outputFile flags)
        exclusions = Text.splitOn "," $ Text.pack $ excludeList flags

camelCase s = Text.concat $ map upCase
                          $ filter (/= "ops")
                          $ Text.splitOn "_" s

-- | Upper-case the given text.
upCase :: Text -> Text
upCase = forceCase toUpper

-- | Lower-case the given name, and prevent it from overlapping with a reserved
-- Haskell name.
lowCase :: Text -> Text
lowCase = replaceReservedName . forceCase toLower

forceCase :: (Char -> Char) -> Text -> Text
forceCase convert s = maybe "" (\(c, cs) -> Text.cons (convert c) cs)
                      (Text.uncons s)

imports = stack [
      "import Data.ByteString (ByteString)"
    , "import Data.Complex (Complex)"
    , "import Data.Int (Int8, Int16, Int32, Int64)"
    , "import Data.Word (Word8, Word16)"
    , "import Lens.Family2 ((.~), (&))"
    , "import TensorFlow.Build"
    , "import TensorFlow.BuildOp"
    , "import TensorFlow.Tensor"
    , "import TensorFlow.Types"
      ]

renderDef :: OpDef -> Doc
renderDef d =
  stack [
      haddocks
    , n <+> "::" <+> hang 0 (typeSig d)
    , n <+> hang 0 args <+> "|" <+> funcGuard <+> "=" </>  -- args are indented
            -- the body needs to be indented wrt the name
            indent indentation functionBody
    , extras  -- just for debug
    ]
  where
    n = strictText $ fixOpName (d ^. name)
    args = sep $ [hsName | (_, hsName) <- mandatoryAttrs] ++ tensorArgs
    tensorArgs = [strictText $ lowCase (a ^. name) | a <- d ^. inputArg]
    fixOpName = lowCase
    funcGuard = "eqLengthGuard" <+> brackets (commasep entries)
      where
        entries =
            [ parens $ quotedText nAttr <> comma <+>
              brackets (commasep $ toList $
              NE.map renderTensorName tensorNames)
            | (nAttr, tensorNames) <- Map.toList $ numberAttrMap d
            ]
        renderTensorName x = parens $ quotedText x <> comma <+>
                             "length" <+> strictText x
    -- Uses hang 0 to align the argument vertically on multiple lines.
    functionBody = buildFunction <+> parens (hang 0 (stack buildOpParts))
                                 </> indent indentation (sep tensorArgs)
    buildFunction
        | null outputListsSizes = "buildOp"
        | otherwise = "buildListOp" <+> brackets (commasep outputListsSizes)
    outputListsSizes = [ strictText numberAttrName
                       | o <- d ^. outputArg
                       , let numberAttrName = o ^. numberAttr
                       , not (Text.null numberAttrName) &&
                         numberAttrName `Map.member` mandatoryAttrMap d
                       ]
    buildOpParts =
        "opDef" <+> quotedText (d ^. name) :
        -- Renders tensor arguments.
        [ "& opAttr" <+> quotedText tfName <+>
          ".~ tensorType (undefined ::" <+> strictText hsName <> ")"
        | (tfName, (hsName, _)) <- Map.toList typeMap
        ] ++
        -- Renders mandatory attributes as function parameters.
        [ "& opAttr" <+> dquotes tfName <+> ".~" <+> hsName
        | (tfName, hsName) <- mandatoryAttrs
        ] ++
        -- Renders sizes of tensor list types having number_attr.
        [ "& opAttr" <+> quotedText nAttr <+> ".~" <+>
          "(fromIntegral (length" <+> strictText (head tensorNames) <> ") :: Int64)"
        | (nAttr, tensorNames) <- Map.toList $ numberAttrMap d
        ]
    mandatoryAttrs = [(strictText tf, strictText hs)
                     | (tf, (hs, _, _)) <- Map.toList (mandatoryAttrMap d)
                     ]
    haddocks = "-- |" <+> multilineComment (d ^. summary) (d ^. description)
    extras = enclose "{-\n" "\n-}" $
             strictText $ Text.pack $
             showMessage ((def :: OpDef)
                          & inputArg .~ (d ^. inputArg)
                          & outputArg .~ (d ^. outputArg)
                          & attr .~ (d ^. attr))
    typeMap = opDefTypeMap d

-- | Makes a quoted string doc out of the given text value.
quotedText :: Text.Text -> Doc
quotedText = dquotes . strictText

-- | typeSig renders the type signature of the given OpDef.
typeSig :: OpDef -> Doc
typeSig d =
    foralls <+> constraints <+/>
    signatureFold (mandatoryAttrInputs ++ tensorInputs ++ [outputs])
  where
    foralls | Map.null typeMap = empty
            | otherwise =
              "forall"
              <+> sep (refTypes ++ map (strictText . fst) (Map.elems typeMap))
              <+> "."
    constraints | Map.null typeMap = empty
                | otherwise =
                  tuple (concatMap
                         (\(t, aDef) ->
                           "TensorType" <+> strictText t
                           : maybeToList (oneOfRestrictions aDef t))
                         (Map.elems typeMap)) <+> "=>"
    tensorInputs = zipWith tensorArg refTypes (d ^. inputArg)
    refTypes = map (\x -> "v" <> int x) [1..length (d ^. inputArg)]
    tensorArg refType arg = wrapArg refType arg <+>
                            hang 0 ("-- ^" <+> argComment arg)
    -- Argument type is a list of tensors if number_attr is set;
    -- otherwise it's a single Tensor.
    wrapArg refType arg =
        if Text.null (arg ^. numberAttr) then typ else brackets typ
      where typ = tensorType refType arg
    tensorType refType arg =
      "Tensor" <+> refType <+> maybe directType strictText indirectType
      where
        indirectType = fmap fst (Map.lookup (arg ^. typeAttr) typeMap)
        directType = dtTypeToDoc (arg ^. type')
    outputs =
      case d ^. outputArg of
        []  -> "ControlNode"
        [o] -> wrappedOutput o <+> "-- ^" <+> argComment o
        os  -> renderTupleResult os
    wrappedOutput = wrapArg "Value"
    -- Tuple result case is rendered differently to give
    -- individual elements their own comments.
    renderTupleResult os =
        stack $ [ tuple (map wrappedOutput os)
                , flatten commentSummary
                ] ++ map commentDetails os
      where
        commentSummary = "-- ^" <+> tuple [bold (o ^. name) | o <- os]
        commentDetails o =
          stack [ "--"
                , "-- *" <+> argComment o
                ]
    signatureFold = folddoc (\x y -> x </> "->" <+> y)
    mandatoryAttrInputs = [
      dtTypeToDoc dtType <+>
          hang 0 ("-- ^" <+> argComment' tfName descr)
      | (tfName, (_, dtType, descr)) <- Map.toList $ mandatoryAttrMap d]
    typeMap = opDefTypeMap d

-- | Returns the type restriction for the given tensor type if the
-- set of allowed types is not empty (i.e. restricted).
oneOfRestrictions :: AttrDef -> Text -> Maybe Doc
oneOfRestrictions aDef tName = do
    typs <- onAttrType (^. templateRestrictions) aDef
    guard $ not $ null typs
    let typeList = commasep $ map strictText $
                   Set.toList $ Set.fromList $
                   map dtTypeToHaskell typs
    return $ "OneOf" <+> "'" <> brackets typeList <+> strictText tName

-- | Identifies the attributes used as tensor cardinalities. In such
-- cases a list of tensors is supplied as an input_arg. The number of
-- such inputs is communicated as a separate opAttr.
-- The result key is TensorFlow attribute name and the value is the
-- tensor names which have number_attr set to the result key.
numberAttrMap :: OpDef -> Map.Map Text.Text (NonEmpty Text.Text)
numberAttrMap d = Map.fromListWith (Semigroup.<>) [
    (nAttr, replaceReservedName (inp ^. name) :| [])
    | inp <- d ^. inputArg
    , nAttr <- [inp ^. numberAttr]
    , not (Text.null nAttr)
    ]

argComment :: OpDef'ArgDef -> Doc
argComment arg = argComment' (arg ^. name) (arg ^. description)

argComment' :: Text.Text -> Text.Text -> Doc
argComment' argName argDesc =
    bold argName <> splitMultilineText (":" <+>) argDesc

bold :: Text.Text -> Doc
bold n = strictText ("__" <> n <> "__")

opDefTypeMap :: OpDef -> Map.Map Text.Text (Text.Text, AttrDef)
opDefTypeMap d =
    Map.fromList [(n, (lowCase n, a)) | (n, a) <- attrList d, isType a]

attrList :: OpDef -> [(Text.Text, AttrDef)]
attrList d = [(a ^. name, attrDef a) | a <- d ^. attr]

isType :: AttrDef -> Bool
isType = fromMaybe False . onAttrType (const True)

-- | Applies the given function to the data type. Is this a Prism?
onAttrType :: (Template DataType -> a) -> AttrDef -> Maybe a
onAttrType f x = case x ^. attrTemplate of
    AttrSingle (AttrType a) -> Just (f a)
    _ -> Nothing

-- | mandatoryAttrMap contains the attributes chosen by
-- isMandatoryAttr, excluding those which are derived from list of
-- tensor arguments. The key is the TF name of the attribute. The
-- value tuple is (haskell name, TF type, attribute comment).
mandatoryAttrMap :: OpDef -> Map.Map Text.Text (Text.Text, DataType, Text.Text)
mandatoryAttrMap d =
    Map.fromList [ (n, (lowCase n, dtType, a ^. attrOriginal.description))
                 | (n, a) <- attrList d
                 , Just dtType <- [isMandatoryAttr a]
                 -- Excludes the attributes rendered as list lengths.
                 , n `Map.notMember` numberAttrMap d
                 ]

-- | Inspects the attribute and if it is one of the implemented
-- non-tensor values lacking default, then returns Just the TF type.
isMandatoryAttr :: AttrDef -> Maybe DataType
isMandatoryAttr x =
   case x ^. attrTemplate of
     AttrSingle (AttrBool y)  -> noDefault DT_BOOL y
     AttrSingle (AttrInt64 y) -> noDefault DT_INT64 y
     AttrSingle (AttrFloat y) -> noDefault DT_FLOAT y
     _ -> Nothing
   where
     noDefault typ y = maybe (Just typ) (const Nothing) (y ^. templateDefault)

dtTypeToDoc :: DataType -> Doc
dtTypeToDoc = strictText . dtTypeToHaskell

-- NOTE: The cases of this function should be kept in sync with
-- TensorFlow.Types.AllTensorTypes.
dtTypeToHaskell :: DataType -> Text.Text
dtTypeToHaskell DT_BOOL = "Bool"
dtTypeToHaskell DT_BFLOAT16 = "Data.Word.Word16"
dtTypeToHaskell DT_COMPLEX128 = "(Data.Complex.Complex Double)"
dtTypeToHaskell DT_COMPLEX64 = "(Data.Complex.Complex Float)"
dtTypeToHaskell DT_DOUBLE = "Double"
dtTypeToHaskell DT_FLOAT = "Float"
dtTypeToHaskell DT_INT16 = "Data.Int.Int16"
dtTypeToHaskell DT_INT32 = "Data.Int.Int32"
dtTypeToHaskell DT_INT64 = "Data.Int.Int64"
dtTypeToHaskell DT_INT8 = "Data.Int.Int8"
dtTypeToHaskell DT_QINT32 = "Data.Int.Int32"  -- TODO(gnezdo): make unique
dtTypeToHaskell DT_QINT8 = "Data.Word.Word8"  -- TODO(gnezdo): make unique
dtTypeToHaskell DT_QINT16 = "Data.Int.Int16"  -- TODO(gnezdo): make unique
dtTypeToHaskell DT_QUINT16 = "Data.Word.Word16"  -- TODO(gnezdo): make unique
dtTypeToHaskell DT_QUINT8 = "Data.Word.Word8"  -- TODO(gnezdo): make unique
dtTypeToHaskell DT_STRING = "Data.ByteString.ByteString"
dtTypeToHaskell DT_UINT16 = "Data.Word.Word16"
dtTypeToHaskell DT_HALF = "Data.Word.Word16"  -- TODO(gnezdo): make unique
dtTypeToHaskell DT_UINT8 = "Data.Word.Word8"
dtTypeToHaskell x =
    Text.pack $ "Unsupported type in dtTypeToHaskell: " ++ show x

-- | haddockComment escapes TensorFlow doc strings into haddock.
-- TODO(gnezdo): deal with the markup.
haddockComment :: Text.Text -> Doc
haddockComment = strictText

multilineComment :: Text.Text -> Text.Text -> Doc
multilineComment summary' detail =
    haddockComment summary' </>
    splitMultilineText insertParagraphAndComment detail
  where insertParagraphAndComment x = "--" </> "--" <+> x

-- | Converts the given multi-line detail string into
-- a multi-line haddock. Applies the given lead to the
-- first line. Returns an empty document for empty detail.
splitMultilineText :: (Doc -> Doc) -> Text.Text -> Doc
splitMultilineText lead detail =
  case Text.lines detail of
    [] -> empty
    (l : ls) -> stack $ lead (haddockComment l)
                      : map (("--" <+>) . haddockComment) ls

replaceReservedName :: Text -> Text
replaceReservedName n
    | n `Set.member` reservedKeywords = n <> "'"
    | otherwise = n

indentation = 4

reservedKeywords :: Set.Set Text
reservedKeywords = Set.fromList $
    -- Haskell2010 keywords:
    -- https://www.haskell.org/onlinereport/haskell2010/haskellch2.html#x7-180002.4
    -- We don't include keywords that are allowed to be variable names,
    -- in particular: "as", "forall", and "hiding".
    [ "case"
    , "class"
    , "data"
    , "default"
    , "deriving"
    , "do"
    , "else"
    , "foreign"
    , "if"
    , "import"
    , "in"
    , "infix"
    , "infixl"
    , "infixr"
    , "instance"
    , "let"
    , "module"
    , "newtype"
    , "of"
    , "then"
    , "type"
    , "where"
    ]
    ++  -- Nonstandard extensions
    [ "mdo"   -- RecursiveDo
    , "rec"   -- Arrows, RecursiveDo
    , "proc"  -- Arrows
    ]

max_stars_repo_name: gleber/tensorflow-haskell

max_issues_repo_name: gleber/tensorflow-haskell

max_forks_repo_name: gleber/tensorflow-haskell


text:
-- | This module implements various functions that return a probabilistic result,
-- defined as unitary operators, and quantum computations.
module QIO.QIORandom where

import Data.Monoid as Monoid
import QIO.QioSyn
import QIO.Qdata
import QIO.Qio
import Data.Complex

-- | The exponentiated Pauli-X rotation
rX :: RR -> Rotation
rX r (x,y) = if x==y then (cos (r/2):+0) else (0:+ (-(sin (r/2))))

-- | The exponentiated Pauli-Y rotation
rY :: RR -> Rotation
rY r (x,y) = if x==y then (cos (r/2):+0) else (s * sin (r/2):+0) where s = if x then 1 else -1

-- | Applies a Hadamard rotation to each qubit in the given list of qubits
hadamards :: [Qbit] -> U
hadamards [] = mempty
hadamards (q:qs) = uhad q `mappend` hadamards qs

-- | returns the highest integer power of 2 that is less than or equal to \x\
pow2 :: Int -> Int
pow2 x = pow2' 0
  where pow2' y | 2^(y+1) > x = 2^y
            | otherwise = pow2' (y+1)

-- | A rotation that, given a qubit in state 0, leaves it in a super-position of
-- 0 and 1, such that the probability of measuring as state 0 is \ps\.
weightedU :: RR -> Qbit -> U
weightedU ps q | sqrt ps <= 1 = rot q (rX (2*(acos (sqrt ps))))
           | otherwise = error ("weightedU: Invalid Probability: " ++ show ps)
-- | A QIO computation that uses the "weightedU" unitary, to return a Bool that
-- has a probablity of \pf\ of being False.
weightedBool :: RR -> QIO Bool
weightedBool pf = do q <- mkQbit False
                     applyU (weightedU pf q)
                     measQ q

-- | removes any leading Falses from a list of booleans
rlf :: [Bool] -> [Bool]
rlf (False:bs) = rlf bs
rlf bs = bs

-- | removes any leading Falses from the (big-endian) bit-wise representation
-- of the given Int.
rlf_l :: Int -> [Bool]
rlf_l x = rlf (reverse (int2bits x))

-- | returns the number of bits left after calling the "flf_l" function
rlf_n :: Int -> Int
rlf_n x = length (rlf_l x)

-- | Given an Int \max\ that is the largest number required to be represented in
-- a quantum register, this function trims the front off the given register, to
-- leave the number of qubits required to represent \max\.
trim :: Int -> [Qbit] -> [Qbit]
trim max qbs = drop ((length qbs)-(rlf_n max)) qbs

-- | Given an Int \max\, and a quantum register in the state \max\, this function
-- defines a unitary operation that will leave the quantum register in state that
-- has equal probability of being measured in any of the states 0 to \max\.
randomU :: Int -> [Qbit] -> U
randomU max qbs = randomU' max (trim max qbs)
 where
  randomU' _ [] = mempty
  randomU' 0 _ = mempty
  randomU' max (q:qbs) = weightedU (fromIntegral ((max+1)-p)/fromIntegral (max+1)) q
                 `mappend`
                 condQ q (\x -> if x then (randomU (max-p) qbs)
                                else (hadamards qbs))
                  where p = pow2 max

-- | A quantum computation that will return a quantum integer in a state that
-- has equal probabilities of being measured in any of the state 0 to \max\.
randomQInt :: Int -> QIO QInt
randomQInt max = do
  qbs <- mkQ (reverse (int2bits max))
  applyU (randomU max qbs)
  return (QInt (reverse qbs))

-- | A quantum computation that will return a quantum integer in a state that
-- has equal probabilities of being measured in any of the state \min\ to \max\.
randomQIO :: (Int,Int) -> QIO Int
randomQIO (min,max) = do q <- randomInt (max-min)
                         return (q + min)

-- | A quantum computation that measures the outcome of "randomQInt"
randomInt :: Int -> QIO Int
randomInt max = do
  q <- randomQInt max
  measQ q

-- | A quantum computation that returns an integer that is equally likely to be
-- any number in the range 0 to \x\-1
random :: Int -> QIO Int
random x = randomInt (x-1)

-- | This function uses a Quantum computation to simulate the roll of a dice
dice :: IO Int
dice = do
  x <- run (randomInt 5)
  return (x+1)

-- | This function simulates the given number of repitions of dice rolls
dice_rolls :: Int -> IO [Int]
dice_rolls 0 = return []
dice_rolls y = do
  x <- dice
  xs <- dice_rolls (y-1)
  return (x:xs)

-- | Returns the number of occurences of 1 through 6 in the given list of Ints
occs :: [Int] -> (Int,Int,Int,Int,Int,Int)
occs rs = (rs' 1,rs' 2,rs' 3,rs' 4,rs' 5,rs' 6)
  where
-- "unterminated quasiquotation" -> space needed
   rs' x = length [ y | y <-rs, y == x]

-- | Returns the number of occurences of 1 through 6 in the given number of
-- rolls of the dice.
probs' :: Int -> IO (Int,Int,Int,Int,Int,Int)
probs' x = do
  xs <- dice_rolls x
  return (occs xs)

-- | Returns the percentage of occurences of 1 through 6, after the given number
-- of rolls of the dice.
probs :: Int -> IO (RR,RR,RR,RR,RR,RR)
probs x = do
  (a,b,c,d,e,f) <- probs' x
  return (fromIntegral a/x',fromIntegral b/x',fromIntegral c/x',fromIntegral d/x',fromIntegral e/x',fromIntegral f/x')
    where x' = fromIntegral x

ext: hs

max_stars_repo_name: probprob/qio-haskell

max_issues_repo_name: probprob/qio-haskell

max_forks_repo_name: probprob/qio-haskell


text:
{-# language TypeFamilies          #-}
{-# language FlexibleInstances     #-}
{-# language FlexibleContexts      #-}
{-# language MultiParamTypeClasses #-}
{-# language UndecidableInstances  #-}
{-# language Rank2Types            #-}
{-# language ScopedTypeVariables   #-}
{-# language TypeOperators #-}

module Feldspar.Software.Frontend where

import Feldspar.Sugar
import Feldspar.Representation
import Feldspar.Frontend
import Feldspar.Array.Vector hiding (reverse, (++))
import Feldspar.Array.Buffered (ArraysSwap(..))
import Feldspar.Software.Primitive
import Feldspar.Software.Primitive.Backend ()
import Feldspar.Software.Expression
import Feldspar.Software.Representation
import Data.Struct

import Data.Bits (Bits, FiniteBits)
import Data.Complex
import Data.Constraint hiding (Sub)
import Data.Proxy
import Data.List (genericLength)
import Data.Word hiding (Word)

-- syntactic.
import Language.Syntactic (Syntactic(..))
import Language.Syntactic.Functional
import qualified Language.Syntactic as Syn

-- operational-higher.
import qualified Control.Monad.Operational.Higher as Oper

-- imperative-edsl.
import Language.Embedded.Imperative.Frontend.General hiding (Ref, Arr, IArr)
import qualified Language.Embedded.Imperative     as Imp
import qualified Language.Embedded.Imperative.CMD as Imp

-- hardware-edsl.
import qualified Language.Embedded.Hardware.Command.CMD as Hard

-- hmm!
import Feldspar.Hardware.Primitive  (HardwarePrimType(..), HardwarePrimTypeRep(..))
import Feldspar.Hardware.Expression (HType')
import Feldspar.Hardware.Frontend   (HSig, withHType')

import Prelude hiding (length, Word, (<=), (<), (>=), (>))
import qualified Prelude as P

--------------------------------------------------------------------------------
-- * Expressions.
--------------------------------------------------------------------------------

instance Value SExp
  where
    value = sugarSymSoftware . Lit

instance Share SExp
  where
    share = sugarSymSoftware (Let "")

instance Iterate SExp
  where
    loop = sugarSymSoftware ForLoop

instance Cond SExp
  where
    cond = sugarSymSoftware Cond

instance Equality SExp
  where
    (==) = sugarSymPrimSoftware Eq

instance Ordered SExp
  where
    (<)  = sugarSymPrimSoftware Lt
    (<=) = sugarSymPrimSoftware Lte
    (>)  = sugarSymPrimSoftware Gt
    (>=) = sugarSymPrimSoftware Gte

instance Logical SExp
  where
    not  = sugarSymPrimSoftware Not
    (&&) = sugarSymPrimSoftware And
    (||) = sugarSymPrimSoftware Or

instance Multiplicative SExp
  where
    mult = sugarSymPrimSoftware Mul
    div  = sugarSymPrimSoftware Div
    mod  = sugarSymPrimSoftware Mod

instance Bitwise SExp
  where
    complement = sugarSymPrimSoftware BitCompl
    (.&.) = sugarSymPrimSoftware BitAnd
    (.|.) = sugarSymPrimSoftware BitOr
    xor   = sugarSymPrimSoftware BitXor
    sll   = sugarSymPrimSoftware ShiftL
    srl   = sugarSymPrimSoftware ShiftR
    rol   = sugarSymPrimSoftware RotateL
    ror   = sugarSymPrimSoftware RotateR

instance Casting SExp
  where
    i2n = sugarSymPrimSoftware I2N
    i2b = sugarSymPrimSoftware I2B
    b2i = sugarSymPrimSoftware B2I

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

instance (Bounded a, SType a) => Bounded (SExp a)
  where
    minBound = value minBound
    maxBound = value maxBound

instance (Num a, SType' a) => Num (SExp a)
  where
    fromInteger = value . fromInteger
    (+)         = sugarSymPrimSoftware Add
    (-)         = sugarSymPrimSoftware Sub
    (*)         = sugarSymPrimSoftware Mul
    negate      = sugarSymPrimSoftware Neg
    abs         = error "todo: abs not implemeted for `SExp`"
    signum      = error "todo: signum not implemented for `SExp`"

instance (Fractional a, SType' a) => Fractional (SExp a)
  where
    fromRational = value . fromRational
    (/)          = sugarSymPrimSoftware FDiv

instance (Floating a, SType' a) => Floating (SExp a)
  where
    pi    = sugarSymPrimSoftware Pi
    exp   = sugarSymPrimSoftware Exp
    log   = sugarSymPrimSoftware Log
    sqrt  = sugarSymPrimSoftware Sqrt
    (**)  = sugarSymPrimSoftware Pow
    sin   = sugarSymPrimSoftware Sin
    cos   = sugarSymPrimSoftware Cos
    tan   = sugarSymPrimSoftware Tan
    asin  = sugarSymPrimSoftware Asin
    acos  = sugarSymPrimSoftware Acos
    atan  = sugarSymPrimSoftware Atan
    sinh  = sugarSymPrimSoftware Sinh
    cosh  = sugarSymPrimSoftware Cosh
    tanh  = sugarSymPrimSoftware Tanh
    asinh = sugarSymPrimSoftware Asinh
    acosh = sugarSymPrimSoftware Acosh
    atanh = sugarSymPrimSoftware Atanh

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

complex :: (Num a, SType' a, SType' (Complex a)) =>
  SExp a -> -- ^ Real
  SExp a -> -- ^ Imaginary
  SExp (Complex a)
complex = sugarSymPrimSoftware Complex

polar :: (Floating a, SType' a, SType' (Complex a)) =>
  SExp a -> -- ^ Magnitude
  SExp a -> -- ^ Phase
  SExp (Complex a)
polar = sugarSymPrimSoftware Polar

real :: (SType' a, SType' (Complex a)) => SExp (Complex a) -> SExp a
real = sugarSymPrimSoftware Real

imaginary :: (SType' a, SType' (Complex a)) => SExp (Complex a) -> SExp a
imaginary = sugarSymPrimSoftware Imag

magnitude :: (RealFloat a, SType' a, SType' (Complex a)) => SExp (Complex a) -> SExp a
magnitude = sugarSymPrimSoftware Magnitude

phase :: (RealFloat a, SType' a, SType' (Complex a)) => SExp (Complex a) -> SExp a
phase = sugarSymPrimSoftware Phase

conjugate :: (RealFloat a, SType' a, SType' (Complex a)) => SExp (Complex a) -> SExp (Complex a)
conjugate = sugarSymPrimSoftware Conjugate

ilog2 :: (FiniteBits a, Integral a, SType' a) => SExp a -> SExp a
ilog2 a = snd $ P.foldr (\ffi vr -> share vr (step ffi)) (a,0) ffis
  where
    step (ff,i) (v,r) = share (b2i (v > fromInteger ff) .<<. value i) $ \shift ->
      (v .>>. i2n shift, r .|. shift)

    -- [(0x1, 0), (0x3, 1), (0xF, 2), (0xFF, 3), (0xFFFF, 4), ...]
    ffis = (`P.zip` [0..])
         $ P.takeWhile (P.<= (2 P.^ (bitSize a `P.div` 2) - 1 :: Integer))
         $ P.map ((subtract 1) . (2 P.^) . (2 P.^))
         $ [(0::Integer)..]
    -- Based on: http://graphics.stanford.edu/~seander/bithacks.html#IntegerLog

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

foreignImport
  :: ( Syn.Signature sig
     , fi  ~ Syn.SmartFun dom sig
     , sig ~ Syn.SmartSig fi
     , dom ~ Syn.SmartSym fi
     , dom ~ SoftwareDomain
     , Syn.SyntacticN f fi
     , Type SoftwarePrimType (Syn.DenResult sig)
     )
  => String -> Denotation sig -> f
foreignImport str f = sugarSymSoftware (Construct str f)

--------------------------------------------------------------------------------
-- * Instructions.
--------------------------------------------------------------------------------

desugar :: (Syntactic a, Domain a ~ SoftwareDomain) => a -> SExp (Internal a)
desugar = SExp . Syn.desugar

sugar   :: (Syntactic a, Domain a ~ SoftwareDomain) => SExp (Internal a) -> a
sugar   = Syn.sugar . unSExp

resugar
  :: ( Syntactic a
     , Syntactic b
     , Internal a ~ Internal b
     , Domain a   ~ SoftwareDomain
     , Domain b   ~ SoftwareDomain
     )
  => a -> b
resugar = Syn.resugar

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

instance References Software
  where
    type Reference Software = Ref
    initRef = Software . fmap Ref . mapStructA (Imp.initRef) . resugar
    newRef  = Software . fmap Ref . mapStructA (const Imp.newRef) $ typeRep
    getRef  = Software . fmap resugar . mapStructA getRef' . unRef
    setRef ref
      = Software
      . sequence_
      . zipListStruct setRef' (unRef ref)
      . resugar
    unsafeFreezeRef
      = Software
      . fmap resugar
      . mapStructA freezeRef'
      . unRef

-- Imp.getRef specialized to software.
getRef' :: forall b . SoftwarePrimType b => Imp.Ref b -> Oper.Program SoftwareCMD (Oper.Param2 SExp SoftwarePrimType) (SExp b)
getRef' = withSType (Proxy :: Proxy b) Imp.getRef

-- Imp.setRef specialized to software.
setRef' :: forall b . SoftwarePrimType b => Imp.Ref b -> SExp b -> Oper.Program SoftwareCMD (Oper.Param2 SExp SoftwarePrimType) ()
setRef' = withSType (Proxy :: Proxy b) Imp.setRef

-- 'Imp.unsafeFreezeRef' specialized to software.
freezeRef' :: forall b . SoftwarePrimType b => Imp.Ref b -> Oper.Program SoftwareCMD (Oper.Param2 SExp SoftwarePrimType) (SExp b)
freezeRef' = withSType (Proxy :: Proxy b) Imp.unsafeFreezeRef

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

instance Slicable SExp (Arr a)
  where
    slice from len (Arr o l arr) = Arr (o+from) len arr

instance Finite SExp (Arr a)
  where
    length = arrLength

instance Arrays Software
  where
    type Array Software = Arr
    newArr len
      = Software
      $ fmap (Arr 0 len)
      $ mapStructA (const (Imp.newArr len))
      $ typeRep
    initArr elems
      = Software
      $ fmap (Arr 0 len . Node)
      $ Imp.constArr elems
      where len = value $ genericLength elems
    getArr arr ix
      = Software
      $ fmap resugar
      $ mapStructA (flip getArr' (ix + arrOffset arr))
      $ unArr arr
    setArr arr ix a
      = Software
      $ sequence_
      $ zipListStruct
         (\a' arr' -> setArr' arr' (ix + arrOffset arr) a')
         (resugar a)
      $ unArr arr
    copyArr arr brr
      = Software
      $ sequence_
      $ zipListStruct (\a b ->
          Imp.copyArr (a, arrOffset arr) (b, arrOffset brr) (length brr))
        (unArr arr)
        (unArr brr)

-- 'Imp.getArr' specialized to software.
getArr' :: forall b . SoftwarePrimType b
  => Imp.Arr Index b -> SExp Index
  -> Oper.Program SoftwareCMD (Oper.Param2 SExp SoftwarePrimType) (SExp b)
getArr' = withSType (Proxy :: Proxy b) Imp.getArr

-- 'Imp.setArr' specialized to software.
setArr' :: forall b . SoftwarePrimType b
  => Imp.Arr Index b -> SExp Index -> SExp b
  -> Oper.Program SoftwareCMD (Oper.Param2 SExp SoftwarePrimType) ()
setArr' = withSType (Proxy :: Proxy b) Imp.setArr

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

instance Syntax SExp a => Indexed SExp (IArr a)
  where
    type ArrElem (IArr a) = a
    (!) (IArr off len a) ix = resugar $ mapStruct index a
      where
        index :: forall b . SoftwarePrimType b => Imp.IArr Index b -> SExp b
        index arr = sugarSymPrimSoftware
          (GuardVal InternalAssertion "arrIndex: index out of bounds.")
          (ix < len)
          (sugarSymPrimSoftware (ArrIx arr) (ix + off) :: SExp b)

instance Slicable SExp (IArr a)
  where
    slice from len (IArr o l arr) = IArr (o+from) len arr

instance Finite SExp (IArr a)
  where
    length = iarrLength

instance IArrays Software
  where
    type IArray Software = IArr
    unsafeFreezeArr arr
      = Software
      $ fmap (IArr (arrOffset arr) (length arr))
      $ mapStructA (Imp.unsafeFreezeArr)
      $ unArr arr
    unsafeThawArr iarr
      = Software
      $ fmap (Arr (iarrOffset iarr) (length iarr))
      $ mapStructA (Imp.unsafeThawArr)
      $ unIArr iarr

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

-- | Short-hand for software pull vectors.
type SPull a = Pull SExp a

-- | Short-hand for software push vectors.
type SPush a = Push Software a

-- | Short-hand for software manifest vectors.
type SManifest a = Manifest Software a

instance Syntax SExp (SExp a) => Pushy Software (IArr (SExp a)) (SExp a)
  where
    toPush iarr = toPush (M iarr :: Manifest Software (SExp a))

instance ViewManifest Software (IArr (SExp a)) (SExp a)
  where
    viewManifest = Just . M

instance Manifestable Software (IArr (SExp a)) (SExp a)

instance ArraysSwap Software
  where
    unsafeArrSwap arr brr = Software $ sequence_ $ zipListStruct Imp.unsafeSwapArr
      (unArr arr)
      (unArr brr)

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

instance Control Software
  where
    iff c t f
      = Software
      $ Imp.iff (resugar c)
          (unSoftware t)
          (unSoftware f)

instance Loop Software
  where
    while c body
      = Software
      $ Imp.while
          (fmap resugar $ unSoftware c)
          (unSoftware body)
    for lower step upper body
      = Software
      $ Imp.for
          (resugar lower, step, Imp.Incl $ resugar upper)
          (unSoftware . body . resugar)

instance Assert Software
  where
    assert = assertLabel $ UserAssertion ""

assertLabel :: AssertionLabel -> SExp Bool -> String -> Software ()
assertLabel lbl cond msg = Software $ Oper.singleInj $ Assert lbl cond msg

--------------------------------------------------------------------------------
-- ** Software instructions.
--------------------------------------------------------------------------------

--------------------------------------------------------------------------------
-- *** Assertions.

guard :: Syntax SExp a => SExp Bool -> String -> a -> a
guard = guardLabel $ UserAssertion ""

guardLabel :: Syntax SExp a => AssertionLabel -> SExp Bool -> String -> a -> a
guardLabel lbl cond msg = sugarSymSoftware (GuardVal lbl msg) cond

hint :: (Syntax SExp a, Syntax SExp b, Primitive SExp (Internal a))
  => a -- ^ Value to be used in invariant.
  -> b -- ^ Result value.
  -> b
hint x y = sugarSymSoftware HintVal x y

--------------------------------------------------------------------------------
-- *** File handling.

-- | Open a file.
fopen :: FilePath -> IOMode -> Software Handle
fopen file = Software . Imp.fopen file

-- | Close a file.
fclose :: Handle -> Software ()
fclose = Software . Imp.fclose

-- | Check for end of file.
feof :: Handle -> Software (SExp Bool)
feof = Software . Imp.feof

-- | Put a primitive value to a handle.
fput :: (Formattable a, SType' a)
    => Handle
    -> String  -- ^ Prefix.
    -> SExp a  -- ^ Expression to print.
    -> String  -- ^ Suffix.
    -> Software ()
fput h pre e post = Software $ Imp.fput h pre e post

-- | Get a primitive value from a handle.
fget :: (Formattable a, SType' a) => Handle -> Software (SExp a)
fget = Software . Imp.fget

-- | Handle to \stdin\.
stdin :: Handle
stdin = Imp.stdin

-- | Handle to \stdout\.
stdout :: Handle
stdout = Imp.stdout

--------------------------------------------------------------------------------
-- *** Printing.

class PrintfType r
  where
    fprf :: Handle -> String -> [Imp.PrintfArg SExp SoftwarePrimType] -> r

instance (a ~ ()) => PrintfType (Software a)
  where
    fprf h form = Software . Oper.singleInj . Imp.FPrintf h form . reverse

instance (Formattable a, SType' a, PrintfType r) => PrintfType (SExp a -> r)
  where
    fprf h form as = \a -> fprf h form (Imp.PrintfArg a : as)

-- | Print to a handle. Accepts a variable number of arguments.
fprintf :: PrintfType r => Handle -> String -> r
fprintf h format = fprf h format []

-- | Print to @stdout@. Accepts a variable number of arguments.
printf :: PrintfType r => String -> r
printf = fprintf Imp.stdout

--------------------------------------------------------------------------------
-- *** Memory.

-- | Software argument specialized to software primitives.
type SArg = Argument SoftwarePrimType

-- | Establish a memory-mapping to a hardware signature.
mmap :: String -> HSig a -> Software (Address a)
mmap address sig =
  do pointer <- Software $ Oper.singleInj $ MMap address sig
     return $ Address pointer sig

-- | Call a memory-mapped component.
call :: Address a -> SArg (Soften a) -> Software ()
call address arg = Software $ Oper.singleInj $ Call address arg

-- | ...
nil :: SArg ()
nil = Nil

-- | ...
(>:) :: forall a b . (SType' a, HType' a, Integral a)
  => Ref (SExp a) -> SArg b -> SArg (Ref (SExp a) -> b)
(>:) = withHType' (Proxy :: Proxy a) ARef

(>.) :: forall a b . (SType' a, HType' a, Integral a)
  => SExp a -> SArg b -> SArg (Ref (SExp a) -> b)
(>.) v = undefined

-- | ...
(>>:) :: forall a b . (SType' a, HType' a, Integral a)
  => Arr (SExp a) -> SArg b -> SArg (Arr (SExp a) -> b)
(>>:) = withHType' (Proxy :: Proxy a) AArr

(>>.) :: forall a b . (SType' a, HType' a, Integral a)
  => IArr (SExp a) -> SArg b -> SArg (Arr (SExp a) -> b)
(>>.) = undefined

infixr 1 >:, >>:

--------------------------------------------------------------------------------
-- *** C specific.

-- | Create a null pointer
newPtr :: SType' a => Software (Ptr a)
newPtr = newNamedPtr "p"

-- | Create a named null pointer
--
-- The provided base name may be appended with a unique identifier to avoid name
-- collisions.
newNamedPtr :: SType' a => String -> Software (Ptr a)
newNamedPtr = Software . Imp.newNamedPtr

-- | Cast a pointer to an array
ptrToArr :: SType' a => Ptr a -> SExp Length -> Software (Arr (SExp a))
ptrToArr ptr len = fmap (Arr 0 len . Node) $ Software $ Imp.ptrToArr ptr

-- | Create a pointer to an abstract object. The only thing one can do with such
-- objects is to pass them to 'callFun' or 'callProc'.
newObject
    :: String  -- ^ Object type
    -> Bool    -- ^ Pointed?
    -> Software Object
newObject = newNamedObject "obj"

-- | Create a pointer to an abstract object. The only thing one can do with such
-- objects is to pass them to 'callFun' or 'callProc'.
--
-- The provided base name may be appended with a unique identifier to avoid name
-- collisions.
newNamedObject
    :: String  -- ^ Base name
    -> String  -- ^ Object type
    -> Bool    -- ^ Pointed?
    -> Software Object
newNamedObject base t p = Software $ Imp.newNamedObject base t p

-- | Add an @#include@ statement to the generated code
addInclude :: String -> Software ()
addInclude = Software . Imp.addInclude

-- | Add a global definition to the generated code
--
-- Can be used conveniently as follows:
--
-- > {-# LANGUAGE QuasiQuotes #-}
-- >
-- > import Feldspar.IO
-- >
-- > prog = do
-- >     ...
-- >     addDefinition myCFunction
-- >     ...
-- >   where
-- >     myCFunction = [cedecl|
-- >       void my_C_function( ... )
-- >       {
-- >           // C code
-- >           // goes here
-- >       }
-- >       |]
addDefinition :: Imp.Definition -> Software ()
addDefinition = Software . Imp.addDefinition

-- | Declare an external function
addExternFun :: SType' res
    => String    -- ^ Function name
    -> proxy res -- ^ Proxy for expression and result type
    -> [FunArg SExp SoftwarePrimType]
                 -- ^ Arguments (only used to determine types)
    -> Software ()
addExternFun fun res args = Software $ Imp.addExternFun fun res args

-- | Declare an external procedure
addExternProc
    :: String -- ^ Procedure name
    -> [FunArg SExp SoftwarePrimType]
              -- ^ Arguments (only used to determine types)
    -> Software ()
addExternProc proc args = Software $ Imp.addExternProc proc args

-- | Call a function
callFun :: SType' a
    => String -- ^ Function name
    -> [FunArg SExp SoftwarePrimType]
              -- ^ Arguments
    -> Software (SExp a)
callFun fun as = Software $ Imp.callFun fun as

-- | Call a procedure
callProc
    :: String -- ^ Function name
    -> [FunArg SExp SoftwarePrimType]
              -- ^ Arguments
    -> Software ()
callProc fun as = Software $ Imp.callProc fun as

-- | Call a procedure and assign its result
callProcAssign :: Assignable obj
    => obj    -- ^ Object to which the result should be assigned
    -> String -- ^ Procedure name
    -> [FunArg SExp SoftwarePrimType]
              -- ^ Arguments
    -> Software ()
callProcAssign obj fun as = Software $ Imp.callProcAssign obj fun as

-- | Declare and call an external function
externFun :: SType' res
    => String -- ^ Procedure name
    -> [FunArg SExp SoftwarePrimType]
              -- ^ Arguments
    -> Software (SExp res)
externFun fun args = Software $ Imp.externFun fun args

-- | Declare and call an external procedure
externProc
    :: String -- ^ Procedure name
    -> [FunArg SExp SoftwarePrimType]
              -- ^ Arguments
    -> Software ()
externProc proc args = Software $ Imp.externProc proc args

-- | Generate code into another translation unit
inModule :: String -> Software () -> Software ()
inModule mod = Software . Imp.inModule mod . unSoftware

-- | Get current time as number of seconds passed today
getTime :: Software (SExp Double)
getTime = Software Imp.getTime

-- | Constant string argument
strArg :: String -> FunArg SExp SoftwarePrimType
strArg = Imp.strArg

-- | Value argument
valArg :: SoftwarePrimType a => SExp a -> FunArg SExp SoftwarePrimType
valArg = Imp.valArg

-- | Reference argument
refArg :: SoftwarePrimType (Internal a) => Ref a -> FunArg SExp SoftwarePrimType
refArg (Ref r) = Imp.refArg (extractNode r)

-- | Mutable array argument
arrArg :: SoftwarePrimType (Internal a) => Arr a -> FunArg SExp SoftwarePrimType
arrArg (Arr o _ a) = Imp.offset (Imp.arrArg (extractNode a)) o

-- | Immutable array argument
iarrArg :: SoftwarePrimType (Internal a) => IArr a -> FunArg SExp SoftwarePrimType
iarrArg (IArr o _ a) = Imp.offset (Imp.iarrArg (extractNode a)) o

-- | Abstract object argument
objArg :: Object -> FunArg SExp SoftwarePrimType
objArg = Imp.objArg

-- | Named constant argument
constArg
    :: String  -- ^ Type
    -> String  -- ^ Named constant
    -> FunArg SExp SoftwarePrimType
constArg = Imp.constArg

-- | Modifier that takes the address of another argument
addr :: FunArg SExp SoftwarePrimType -> FunArg SExp SoftwarePrimType
addr = Imp.addr

-- | Modifier that dereferences another argument
deref :: FunArg SExp SoftwarePrimType -> FunArg SExp SoftwarePrimType
deref = Imp.deref

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

-- Swap an `Imp.FreePred` constraint with a `SoftwarePrimType` one.
withSType :: forall a b . Proxy a
  -> (Imp.FreePred SExp a => b)
  -> (SoftwarePrimType  a => b)
withSType _ f = case predicateDict (softwareRep :: SoftwarePrimTypeRep a) of
  Dict -> f

-- Proves that a type from `SoftwarePrimTypeRep` satisfies `Imp.FreePred`.
predicateDict :: SoftwarePrimTypeRep a -> Dict (Imp.FreePred SExp a)
predicateDict rep = case rep of
  BoolST   -> Dict
  Int8ST   -> Dict
  Int16ST  -> Dict
  Int32ST  -> Dict
  Int64ST  -> Dict
  Word8ST  -> Dict
  Word16ST -> Dict
  Word32ST -> Dict
  Word64ST -> Dict
  FloatST  -> Dict
  ComplexFloatST  -> Dict
  ComplexDoubleST -> Dict

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

max_stars_repo_name: markus-git/co-feldspar

max_issues_repo_name: markus-git/co-feldspar

max_forks_repo_name: markus-git/co-feldspar


text:
module Sim.Sim
  ( stepOde,
  )
where

import Linear hiding (cross)
import Math.Frames
import Math.Spatial
import Math.Vectorize (Vectorize (..))
import Models.Aircraft
import Numeric.GSL.ODE (odeSolve)
import Numeric.LinearAlgebra (Vector (..), fromList, toList, toRows)
import Physics.RigidBody

runOde ::
  RigidBodyConstants ->
  V3T Body Double ->
  V3T Body Double ->
  AircraftState Double ->
  AircraftState Double
runOde constants forces moments state =
  AircraftState $ stepRigidBodyOde constants forces moments (dsRigidBody state)

-- | Steps the simulation forward by one timestep us RK45 ode solver.
--
-- This is pretty ugly because we have to convert back and forth from a list of
-- doubles in order to interact with GSL's ode library.
stepOde ::
  RigidBodyConstants ->
  V3T Body Double ->
  V3T Body Double ->
  AircraftState Double ->
  AircraftState Double
stepOde constants forces moments state = devectorize $ toList solution
  where
    solution :: Vector Double
    (_ : solution : _) =
      toRows $
        odeSolve
          wrapOde
          (vectorize state)
          (fromList [0.0, 0.01])

    wrapOde :: Double -> [Double] -> [Double]
    wrapOde t0 x0 =
      vectorize $
        runOde
          constants
          forces
          moments
          (devectorize x0)


text:
{-# LANGUAGE ForeignFunctionInterface #-}
module Data.Matrix.Static.Internal
    ( c_dd_mul
    , c_ds_mul
    , c_sd_mul
    , c_ss_mul
    , c_ss_cmul
    , c_sd_plus
    , c_ss_plus
    , c_inverse
    , c_cholesky
    , c_eig
    , c_eigs
    , c_eigsh
    , c_seigs
    , c_seigsh
    , c_geigsh
    , c_bdcsvd

    , computationInfo
    ) where

import Data.Complex (Complex)
import Foreign
import Foreign.C.Types
import Foreign.C.String

-------------------------------------------------------------------------------
-- Arithmetic
-------------------------------------------------------------------------------
foreign import ccall "eigen_dd_mul"
    c_dd_mul :: CInt
          -> Ptr a -> CInt -> CInt
          -> Ptr a -> CInt -> CInt
          -> Ptr a -> CInt -> CInt
          -> IO CString

foreign import ccall "eigen_ds_mul"
    c_ds_mul :: CInt
          -> Ptr a -> CInt -> CInt
          -> Ptr a -> CInt -> CInt
          -> Ptr a -> Ptr CInt -> Ptr CInt -> CInt -> CInt -> CInt
          -> IO CString

foreign import ccall "eigen_sd_mul"
    c_sd_mul :: CInt
          -> Ptr a -> CInt -> CInt
          -> Ptr a -> Ptr CInt -> Ptr CInt -> CInt -> CInt -> CInt
          -> Ptr a -> CInt -> CInt
          -> IO CString

foreign import ccall "eigen_sd_plus"
    c_sd_plus :: CInt
          -> Ptr a -> CInt -> CInt
          -> Ptr a -> Ptr CInt -> Ptr CInt -> CInt -> CInt -> CInt
          -> Ptr a -> CInt -> CInt
          -> IO CString

foreign import ccall "eigen_ss_mul"
    c_ss_mul :: CInt
          -> Ptr (Ptr a) -> Ptr CInt -> Ptr (Ptr CInt) -> CInt -> CInt -> Ptr CInt
          -> Ptr a -> Ptr CInt -> Ptr CInt -> CInt -> CInt -> CInt
          -> Ptr a -> Ptr CInt -> Ptr CInt -> CInt -> CInt -> CInt
          -> IO CString

foreign import ccall "eigen_ss_cmul"
    c_ss_cmul :: CInt
          -> Ptr (Ptr a) -> Ptr CInt -> Ptr (Ptr CInt) -> CInt -> CInt -> Ptr CInt
          -> Ptr a -> Ptr CInt -> Ptr CInt -> CInt -> CInt -> CInt
          -> Ptr a -> Ptr CInt -> Ptr CInt -> CInt -> CInt -> CInt
          -> IO CString

foreign import ccall "eigen_ss_plus"
    c_ss_plus :: CInt
          -> Ptr (Ptr a) -> Ptr CInt -> Ptr (Ptr CInt) -> CInt -> CInt -> Ptr CInt
          -> Ptr a -> Ptr CInt -> Ptr CInt -> CInt -> CInt -> CInt
          -> Ptr a -> Ptr CInt -> Ptr CInt -> CInt -> CInt -> CInt
          -> IO CString

foreign import ccall "eigen_inverse"
    c_inverse :: CInt
              -> Ptr a -> CInt -> CInt
              -> Ptr a -> CInt -> CInt
              -> IO CString

foreign import ccall "eigen_cholesky"
    c_cholesky :: CInt
               -> Ptr a -> Ptr a ->  CInt -> IO CString

foreign import ccall "eigen_eig"
    c_eig :: Ptr (Complex Double) -> Ptr (Complex Double)
          -> Ptr Double -> CInt -> IO CString

foreign import ccall "spectral_eigs"
    c_eigs :: CInt      -- ^ The number of eigenvectors to search
           -> Ptr (Complex Double)   -- ^
           -> Ptr (Complex Double)   -- ^
           -> Ptr Double   -- ^ Input matrix
           -> CInt         -- ^ Matrix size
           -> CInt   -- ^ ncv
           -> CInt   -- ^ max iterations
           -> Double  -- ^ tolerance
           -> CInt
           -> Double
           -> CInt
           -> IO CInt

foreign import ccall "spectral_eigsh"
    c_eigsh :: CInt
            -> Ptr Double -> Ptr Double
            -> Ptr Double -> CInt
            -> CInt -> CInt -> Double -> CInt -> Double -> CInt
            -> IO CInt

foreign import ccall "spectral_seigs"
    c_seigs :: CInt
            -> Ptr (Complex Double) -> Ptr (Complex Double)
            -> Ptr Double -> Ptr CInt -> Ptr CInt -> CInt -> CInt
            -> CInt -> CInt -> Double -> CInt -> Double -> CInt
            -> IO CInt

foreign import ccall "spectral_seigsh"
    c_seigsh :: CInt
             -> Ptr Double -> Ptr Double
             -> Ptr Double -> Ptr CInt -> Ptr CInt -> CInt -> CInt
             -> CInt -> CInt -> Double -> CInt -> Double -> CInt
             -> IO CInt

foreign import ccall "spectral_geigsh"
    c_geigsh :: CInt
             -> Ptr Double -> Ptr Double
             -> Ptr Double -> CInt
             -> Ptr Double -> Ptr CInt -> Ptr CInt -> CInt
             -> CInt -> CInt -> Double -> CInt
             -> IO CInt

foreign import ccall "eigen_bdcsvd"
    c_bdcsvd :: CInt -> Ptr a -> Ptr b -> Ptr a
             -> Ptr a -> CInt -> CInt -> IO CString

computationInfo :: CInt -> Maybe String
computationInfo 0 = Nothing
computationInfo 1 = Just "NOT_COMPUTED"
computationInfo 2 = Just "NOT_CONVERGING"
computationInfo 3 = Just "NUMERICAL_ISSUE"
computationInfo _ = Just "UNKNOWN ERROR"
{-# INLINE computationInfo #-}

max_stars_repo_name: kaizhang/matrix-sized

max_issues_repo_name: kaizhang/matrix-sized

max_forks_repo_name: kaizhang/matrix-sized


text:
{-# OPTIONS_GHC -Wall -Wno-orphans #-}
{-# LANGUAGE ScopedTypeVariables #-}
{-# LANGUAGE TemplateHaskell #-}
module TestZDD (zddTestGroup) where

import Control.DeepSeq
import Control.Monad
import Control.Monad.ST
import Data.IntMap (IntMap)
import qualified Data.IntMap as IntMap
import Data.IntSet (IntSet)
import qualified Data.IntSet as IntSet
import Data.IORef
import Data.List
import qualified Data.Map.Strict as Map
import Data.Proxy
import Data.Set (Set)
import qualified Data.Set as Set
import Data.Vector (Vector)
import Data.Word
import qualified GHC.Exts as Exts
import Statistics.Distribution
import Statistics.Distribution.ChiSquared (chiSquared)
import System.IO.Unsafe
import qualified System.Random.MWC as Rand
import Test.QuickCheck.Function (apply)
import Test.QuickCheck.Instances.Vector ()
import Test.Tasty
import Test.Tasty.HUnit
import Test.Tasty.QuickCheck
import Test.Tasty.TH

import Data.DecisionDiagram.BDD.Internal.ItemOrder (OrderedItem (..))
import Data.DecisionDiagram.ZDD (ZDD (..), ItemOrder (..))
import qualified Data.DecisionDiagram.ZDD as ZDD

import Utils

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

instance ZDD.ItemOrder a => Arbitrary (ZDD a) where
  arbitrary = do
    vars <- liftM (sortBy (ZDD.compareItem (Proxy :: Proxy a)) . IntSet.toList) arbitrary
    let f vs n = oneof $
          [ return ZDD.empty
          , return ZDD.base
          ]
          ++
          [ do v <- elements vs
               let vs' = dropWhile (\v' -> compareItem (Proxy :: Proxy a) v' v  /= GT) vs
               p0 <- f vs' (n `div` 2)
               p1 <- f vs' (n `div` 2) `suchThat` (/= ZDD.empty)
               return (ZDD.Branch v p0 p1)
          | n > 0, not (null vs)
          ]
    sized (f vars)

  shrink (ZDD.Empty) = []
  shrink (ZDD.Base) = []
  shrink (ZDD.Branch x p0 p1) =
    [p0, p1]
    ++
    [ ZDD.Branch x p0' p1'
    | (p0', p1') <- shrink (p0, p1), p1' /= ZDD.empty
    ]

arbitraryMember :: ZDD.ItemOrder a => ZDD a -> Gen IntSet
arbitraryMember zdd = do
  (seed :: Vector Word32) <- arbitrary
  return $ runST $ do
    gen <- Rand.initialize seed
    ZDD.uniformM zdd gen

-- ------------------------------------------------------------------------
-- Union
-- ------------------------------------------------------------------------

prop_union_unitL :: Property
prop_union_unitL =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      ZDD.empty `ZDD.union` a === a

prop_union_unitR :: Property
prop_union_unitR =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      a `ZDD.union` ZDD.empty === a

prop_union_comm :: Property
prop_union_comm =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b) ->
      a `ZDD.union` b === b `ZDD.union` a

prop_union_assoc :: Property
prop_union_assoc =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b, c) ->
      a `ZDD.union` (b `ZDD.union` c) === (a `ZDD.union` b) `ZDD.union` c

prop_union_idempotent :: Property
prop_union_idempotent =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      a `ZDD.union` a === a

-- ------------------------------------------------------------------------
-- Intersection
-- ------------------------------------------------------------------------

prop_intersection_comm :: Property
prop_intersection_comm =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b) ->
      a `ZDD.intersection` b === b `ZDD.intersection` a

prop_intersection_assoc :: Property
prop_intersection_assoc =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b, c) ->
      a `ZDD.intersection` (b `ZDD.intersection` c) === (a `ZDD.intersection` b) `ZDD.intersection` c

prop_intersection_idempotent :: Property
prop_intersection_idempotent =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      a `ZDD.intersection` a === a

prop_intersection_emptyL :: Property
prop_intersection_emptyL =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      ZDD.empty `ZDD.intersection` a === ZDD.empty

prop_intersection_emptyR :: Property
prop_intersection_emptyR =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      a `ZDD.intersection` ZDD.empty === ZDD.empty

prop_intersection_distL :: Property
prop_intersection_distL =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b, c) ->
      a `ZDD.intersection` (b `ZDD.union` c) === (a `ZDD.intersection` b) `ZDD.union` (a `ZDD.intersection` c)

prop_intersection_distR :: Property
prop_intersection_distR =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b, c) ->
      (a `ZDD.union` b) `ZDD.intersection` c === (a `ZDD.intersection` c) `ZDD.union` (b `ZDD.intersection` c)

-- ------------------------------------------------------------------------
-- Difference
-- ------------------------------------------------------------------------

prop_difference_cancel :: Property
prop_difference_cancel =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      a ZDD.\\ a === ZDD.empty

prop_difference_unit :: Property
prop_difference_unit =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      a ZDD.\\ ZDD.empty === a

prop_union_difference :: Property
prop_union_difference =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b, c) ->
      (a `ZDD.union` b) ZDD.\\ c === (a ZDD.\\ c) `ZDD.union` (b ZDD.\\ c)

prop_difference_union :: Property
prop_difference_union =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b, c) ->
      a ZDD.\\ (b `ZDD.union` c) === (a ZDD.\\ b) ZDD.\\ c

-- ------------------------------------------------------------------------
-- Non-superset
-- ------------------------------------------------------------------------

prop_nonSuperset :: Property
prop_nonSuperset =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b) ->
      let a' = ZDD.toSetOfIntSets a
          b' = ZDD.toSetOfIntSets b
          p xs = and [not (IntSet.isSubsetOf ys xs) | ys <- Set.toList b']
       in ZDD.toSetOfIntSets (a `ZDD.nonSuperset` b) === Set.filter p a'

prop_nonSuperset_cancel :: Property
prop_nonSuperset_cancel =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      a `ZDD.nonSuperset` a === ZDD.empty

prop_nonSuperset_unit :: Property
prop_nonSuperset_unit =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      a `ZDD.nonSuperset` ZDD.empty === a

prop_union_nonSuperset :: Property
prop_union_nonSuperset =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b, c) ->
      (a `ZDD.union` b) `ZDD.nonSuperset` c === (a `ZDD.nonSuperset` c) `ZDD.union` (b `ZDD.nonSuperset` c)

prop_nonSuperset_union :: Property
prop_nonSuperset_union =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b, c) ->
      a `ZDD.nonSuperset` (b `ZDD.union` c) === (a `ZDD.nonSuperset` b) `ZDD.nonSuperset` c

prop_nonSuperset_difference :: Property
prop_nonSuperset_difference =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b) ->
      let c = a `ZDD.nonSuperset` b
          d = a ZDD.\\ b
       in counterexample (show (c, d)) $ c `ZDD.isSubsetOf` d

-- ------------------------------------------------------------------------
-- Minimal hitting sets
-- ------------------------------------------------------------------------

isHittingSetOf :: IntSet -> Set IntSet -> Bool
isHittingSetOf s g = all (\e -> not (IntSet.null (s `IntSet.intersection` e))) g

prop_minimalHittingSetsKnuth_isHittingSet :: Property
prop_minimalHittingSetsKnuth_isHittingSet =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      let b = ZDD.minimalHittingSetsKnuth a
          a' = ZDD.toSetOfIntSets a
          b' = ZDD.toSetOfIntSets b
       in all (`isHittingSetOf` a') b'

prop_minimalHittingSetsImai_isHittingSet :: Property
prop_minimalHittingSetsImai_isHittingSet =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      let b = ZDD.minimalHittingSetsImai a
          a' = ZDD.toSetOfIntSets a
          b' = ZDD.toSetOfIntSets b
       in all (`isHittingSetOf` a') b'

prop_minimalHittingSetsToda_isHittingSet :: Property
prop_minimalHittingSetsToda_isHittingSet =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      let b = ZDD.minimalHittingSetsToda a
          a' = ZDD.toSetOfIntSets a
          b' = ZDD.toSetOfIntSets b
       in all (`isHittingSetOf` a') b'

prop_minimalHittingSetsKnuth_duality :: Property
prop_minimalHittingSetsKnuth_duality =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      let b = ZDD.minimalHittingSetsKnuth a
       in ZDD.minimalHittingSetsKnuth (ZDD.minimalHittingSetsKnuth b) === b

prop_minimalHittingSetsImai_duality :: Property
prop_minimalHittingSetsImai_duality =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      let b = ZDD.minimalHittingSetsImai a
       in ZDD.minimalHittingSetsImai (ZDD.minimalHittingSetsImai b) === b

prop_minimalHittingSetsToda_duality :: Property
prop_minimalHittingSetsToda_duality =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      let b = ZDD.minimalHittingSetsToda a
       in ZDD.minimalHittingSetsToda (ZDD.minimalHittingSetsToda b) === b

prop_minimalHittingSets_Imai_equal_Knuth :: Property
prop_minimalHittingSets_Imai_equal_Knuth =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      ZDD.minimalHittingSetsImai a === ZDD.minimalHittingSetsKnuth a

prop_minimalHittingSets_Toda_equal_Knuth :: Property
prop_minimalHittingSets_Toda_equal_Knuth =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      ZDD.minimalHittingSetsToda a === ZDD.minimalHittingSetsKnuth a

-- ------------------------------------------------------------------------
-- Misc
-- ------------------------------------------------------------------------

prop_empty :: Property
prop_empty =
  forAllItemOrder $ \(_ :: Proxy o) ->
    ZDD.toSetOfIntSets (ZDD.empty :: ZDD o) === Set.empty

prop_base :: Property
prop_base =
  forAllItemOrder $ \(_ :: Proxy o) ->
    ZDD.toSetOfIntSets (ZDD.base :: ZDD o) === Set.singleton IntSet.empty

prop_singleton :: Property
prop_singleton =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \xs ->
      let a :: ZDD o
          a = ZDD.singleton xs
       in counterexample (show a) $ ZDD.toSetOfIntSets a === Set.singleton xs

prop_subsets_member :: Property
prop_subsets_member =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \xs ->
      let a :: ZDD o
          a = ZDD.subsets xs
       in counterexample (show a) $ forAll (subsetOf xs) $ \ys ->
            ys `ZDD.member` a

prop_subsets_member_empty :: Property
prop_subsets_member_empty =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \xs ->
      let a :: ZDD o
          a = ZDD.subsets xs
       in counterexample (show a) $ IntSet.empty `ZDD.member` a

prop_subsets_member_itself :: Property
prop_subsets_member_itself =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \xs ->
      let a :: ZDD o
          a = ZDD.subsets xs
       in counterexample (show a) $ xs `ZDD.member` a

prop_subsets_size :: Property
prop_subsets_size =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \xs ->
      let a :: ZDD o
          a = ZDD.subsets xs
       in counterexample (show a) $ ZDD.size a === (2 :: Integer) ^ (IntSet.size xs)

prop_subsetsAtLeast :: Property
prop_subsetsAtLeast =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrarySmallIntMap $ \(xs :: IntMap Integer) ->
      forAll arbitrary $ \k ->
        let a :: ZDD o
            a = ZDD.subsetsAtLeast xs k
         in counterexample (show a) $
              if ZDD.null a then
                property (k > sum [max 0 w | (_,w) <- IntMap.toList xs])
              else
                forAll (arbitraryMember a) $ \ys ->
                  (ys `IntSet.isSubsetOf` IntMap.keysSet xs)
                  .&&.
                  sum [xs IntMap.! y | y <- IntSet.toList ys] >= k

prop_subsetsAtMost :: Property
prop_subsetsAtMost =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrarySmallIntMap $ \(xs :: IntMap Integer) ->
      forAll arbitrary $ \k ->
        let a :: ZDD o
            a = ZDD.subsetsAtMost xs k
         in counterexample (show a) $
              if ZDD.null a then
                property (k < sum [min 0 w | (_,w) <- IntMap.toList xs])
              else
                forAll (arbitraryMember a) $ \ys ->
                  (ys `IntSet.isSubsetOf` IntMap.keysSet xs)
                  .&&.
                  sum [xs IntMap.! y | y <- IntSet.toList ys] <= k

prop_subsetsExactly :: Property
prop_subsetsExactly =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrarySmallIntMap $ \(xs :: IntMap Integer) ->
      forAll arbitrary $ \k ->
        let a :: ZDD o
            a = ZDD.subsetsExactly xs k
         in counterexample (show a) $
              if ZDD.null a then
                property True
              else
                forAll (arbitraryMember a) $ \ys ->
                  (ys `IntSet.isSubsetOf` IntMap.keysSet xs)
                  .&&.
                  sum [xs IntMap.! y | y <- IntSet.toList ys] === k

prop_subsetsExactly_2 :: Property
prop_subsetsExactly_2 =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrarySmallIntMap $ \(xs :: IntMap Integer) ->
      forAll (gen xs) $ \(ys, k) ->
        let a :: ZDD o
            a = ZDD.subsetsExactly xs k
         in counterexample (show a) $ ys `ZDD.member` a
  where
    gen xs = do
      ys <- sublistOf (IntMap.toList xs)
      return (IntSet.fromList [y | (y,_) <- ys], sum [w | (_,w) <- ys])

prop_subsetsExactlyIntegral :: Property
prop_subsetsExactlyIntegral =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrarySmallIntMap $ \(xs :: IntMap Integer) ->
      forAll arbitrary $ \k ->
        (ZDD.subsetsExactlyIntegral xs k :: ZDD o) === ZDD.subsetsExactly xs k

prop_combinations_are_combinations :: Property
prop_combinations_are_combinations =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \xs ->
      forAll arbitrary $ \(NonNegative k) ->
        let a :: ZDD o
            a = ZDD.combinations xs k
         in counterexample (show a) $
              not (ZDD.null a)
              ==>
              (forAll (arbitraryMember a) $ \ys -> (ys `IntSet.isSubsetOf` xs) .&&. (IntSet.size ys === k))

prop_combinations_size :: Property
prop_combinations_size =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \xs ->
      forAll arbitrary $ \(NonNegative k) ->
        let a :: ZDD o
            a = ZDD.combinations xs k
            n = toInteger $ IntSet.size xs
         in counterexample (show a) $ ZDD.size a === (product [(n - toInteger k + 1)..n] `div` (product [1..toInteger k]))

case_toList_lazyness :: Assertion
case_toList_lazyness = do
  let xss :: ZDD ZDD.AscOrder
      xss = ZDD.subsets (IntSet.fromList [1..128])
  deepseq (take 100 (Exts.toList xss)) $ return ()

prop_toList_sorted :: Property
prop_toList_sorted =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(xss :: ZDD o) ->
      let yss :: [[OrderedItem o]]
          yss = map (sort . map OrderedItem . IntSet.toList) $ take 100 $ Exts.toList xss
       in yss ===  sort yss

prop_toList_fromList :: Property
prop_toList_fromList =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \xss ->
      let a :: ZDD o
          a = Exts.fromList xss
          f = Set.fromList
       in counterexample (show a) $ f (Exts.toList a) === f xss

prop_fromList_toList :: Property
prop_fromList_toList =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      let xss = Exts.toList a
       in counterexample (show xss) $ Exts.fromList xss === a

prop_toSetOfIntSets_fromSetOfIntSets :: Property
prop_toSetOfIntSets_fromSetOfIntSets =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \xss ->
      let a :: ZDD o
          a = ZDD.fromSetOfIntSets xss
       in counterexample (show a) $ ZDD.toSetOfIntSets a === xss

prop_fromSetOfIntSets_toSetOfIntSets :: Property
prop_fromSetOfIntSets_toSetOfIntSets =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      let xss = ZDD.toSetOfIntSets a
       in counterexample (show xss) $ ZDD.fromSetOfIntSets xss === a

prop_insert :: Property
prop_insert =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      forAll arbitrary $ \xs ->
        ZDD.toSetOfIntSets (ZDD.insert xs a) === Set.insert xs (ZDD.toSetOfIntSets a)

prop_insert_idempotent :: Property
prop_insert_idempotent =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      forAll arbitrary $ \xs ->
        let b = ZDD.insert xs a
         in counterexample (show b) $ ZDD.insert xs b === b

prop_delete :: Property
prop_delete =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      forAll (subsetOf (ZDD.flatten a)) $ \xs ->
        ZDD.toSetOfIntSets (ZDD.delete xs a) === Set.delete xs (ZDD.toSetOfIntSets a)

prop_delete_idempotent :: Property
prop_delete_idempotent =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      forAll (subsetOf (ZDD.flatten a)) $ \xs ->
        let b = ZDD.delete xs a
         in counterexample (show b) $ ZDD.delete xs b === b

prop_mapInsert :: Property
prop_mapInsert =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      forAll arbitrary $ \x ->
        ZDD.toSetOfIntSets (ZDD.mapInsert x a) === Set.map (IntSet.insert x) (ZDD.toSetOfIntSets a)

prop_mapInsert_idempotent :: Property
prop_mapInsert_idempotent =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      forAll arbitrary $ \x ->
        let b = ZDD.mapInsert x a
         in counterexample (show b) $ ZDD.mapInsert x b === b

prop_mapDelete :: Property
prop_mapDelete =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      forAll arbitrary $ \x ->
        ZDD.toSetOfIntSets (ZDD.mapDelete x a) === Set.map (IntSet.delete x) (ZDD.toSetOfIntSets a)

prop_mapDelete_idempotent :: Property
prop_mapDelete_idempotent =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      forAll arbitrary $ \x ->
        let b = ZDD.mapDelete x a
         in counterexample (show b) $ ZDD.mapDelete x b === b

prop_change :: Property
prop_change =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
    forAll arbitrary $ \x ->
      let f xs
            | IntSet.member x xs = IntSet.delete x xs
            | otherwise = IntSet.insert x xs
       in ZDD.toSetOfIntSets (ZDD.change x a) === Set.map f (ZDD.toSetOfIntSets a)

prop_change_involution :: Property
prop_change_involution =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      forAll arbitrary $ \x ->
        ZDD.change x (ZDD.change x a) === a

prop_subset1 :: Property
prop_subset1 =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
    forAll arbitrary $ \x ->
      ZDD.toSetOfIntSets (ZDD.subset1 x a) === Set.map (IntSet.delete x) (Set.filter (IntSet.member x) (ZDD.toSetOfIntSets a))

prop_subset0 :: Property
prop_subset0 =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
    forAll arbitrary $ \x ->
      ZDD.toSetOfIntSets (ZDD.subset0 x a) === Set.filter (IntSet.notMember x) (ZDD.toSetOfIntSets a)

prop_member_1 :: Property
prop_member_1 =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      conjoin [counterexample (show xs) (ZDD.member xs a) | xs <- ZDD.toListOfIntSets a]

prop_member_2 :: Property
prop_member_2 =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      forAll (subsetOf (ZDD.flatten a)) $ \s2 ->
        (s2 `ZDD.member` a) === (s2 `Set.member` ZDD.toSetOfIntSets a)

prop_size :: Property
prop_size =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      ZDD.size a === Set.size (ZDD.toSetOfIntSets a)

prop_null_size :: Property
prop_null_size =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      ZDD.null a === (ZDD.size a == (0 :: Int))

prop_isSubsetOf_refl :: Property
prop_isSubsetOf_refl =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      a `ZDD.isSubsetOf` a

prop_isSubsetOf_empty :: Property
prop_isSubsetOf_empty =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      ZDD.empty `ZDD.isSubsetOf` a

prop_isSubsetOf_and_isProperSubsetOf :: Property
prop_isSubsetOf_and_isProperSubsetOf =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b) ->
      (a `ZDD.isSubsetOf` b) === (a `ZDD.isProperSubsetOf` b || a == b)

prop_isProperSubsetOf_not_refl :: Property
prop_isProperSubsetOf_not_refl =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      not (ZDD.isProperSubsetOf a a)

prop_disjoint :: Property
prop_disjoint =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o, b) ->
      ZDD.disjoint a b === ZDD.null (a `ZDD.intersection` b)

case_numNodes :: Assertion
case_numNodes = do
  let bdd, bdd1 :: ZDD ZDD.AscOrder
      bdd = ZDD.Branch 0 (ZDD.Branch 1 ZDD.empty bdd1) (ZDD.Branch 2 ZDD.empty bdd1)
      bdd1 = ZDD.Branch 3 ZDD.empty ZDD.base
  ZDD.numNodes bdd @?= 6

prop_flatten :: Property
prop_flatten =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      ZDD.flatten a === IntSet.unions (ZDD.toListOfIntSets a)

prop_uniformM :: Property
prop_uniformM =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll (arbitrary `suchThat` ((>= (2::Integer)) . ZDD.size)) $ \(a :: ZDD o) ->
      forAll arbitrary $ \(seed :: Vector Word32) ->
        let m :: Integer
            m = ZDD.size a
            n = 1000
            samples = runST $ do
              gen <- Rand.initialize seed
              replicateM n $ ZDD.uniformM a gen
            hist_actual = Map.fromListWith (+) [(s, 1) | s <- samples]
            hist_expected = [(s, fromIntegral n / fromIntegral m) | s <- ZDD.toListOfIntSets a]
            chi_sq = sum [(Map.findWithDefault 0 s hist_actual - cnt) ** 2 / cnt | (s, cnt) <- hist_expected]
            threshold = complQuantile (chiSquared (fromIntegral m - 1)) 0.0001
         in counterexample (show hist_actual ++ " /= " ++ show (Map.fromList hist_expected)) $
              and [xs `ZDD.member` a | xs <- Map.keys hist_actual] .&&.
              counterexample ("χ² = " ++ show chi_sq ++ " >= " ++ show threshold) (chi_sq < threshold)

prop_findMinSum :: Property
prop_findMinSum =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll (arbitrary `suchThat` (not . ZDD.null)) $ \(a :: ZDD o) ->
      forAll arbitrary $ \(weight' :: Fun Int Integer) ->
        let weight = apply weight'
            setWeight s = sum [weight x | x <- IntSet.toList s]
            (obj0, s0) = ZDD.findMinSum weight a
         in counterexample (show [(x, weight x) | x <- IntSet.toList (ZDD.flatten a)]) $
            counterexample (show (obj0, s0)) $
              setWeight s0 === obj0
              .&&.
              conjoin [counterexample (show (s, setWeight s)) $ obj0 <= setWeight s | s <- ZDD.toListOfIntSets a]

prop_findMaxSum :: Property
prop_findMaxSum =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll (arbitrary `suchThat` (not . ZDD.null)) $ \(a :: ZDD o) ->
      forAll arbitrary $ \(weight' :: Fun Int Integer) ->
        let weight = apply weight'
            setWeight s = sum [weight x | x <- IntSet.toList s]
            (obj0, s0) = ZDD.findMaxSum weight a
         in counterexample (show [(x, weight x) | x <- IntSet.toList (ZDD.flatten a)]) $
            counterexample (show (obj0, s0)) $
              setWeight s0 === obj0
              .&&.
              conjoin [counterexample (show (s, setWeight s)) $ obj0 >= setWeight s | s <- ZDD.toListOfIntSets a]

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

case_fold_laziness :: Assertion
case_fold_laziness = do
  let bdd :: ZDD ZDD.AscOrder
      bdd = ZDD.Branch 0 (ZDD.Branch 1 ZDD.Empty ZDD.Base) (ZDD.Branch 2 ZDD.Empty ZDD.Base)
      f x lo _hi =
        if x == 2 then
          error "unused value should not be evaluated"
        else
          lo
  seq (ZDD.fold f id bdd) $ return ()

case_fold'_strictness :: Assertion
case_fold'_strictness = do
  ref <- newIORef False
  let bdd :: ZDD ZDD.AscOrder
      bdd = ZDD.Branch 0 (ZDD.Branch 1 ZDD.Empty ZDD.Base) (ZDD.Branch 2 ZDD.Empty ZDD.Base)
      f x lo _hi = unsafePerformIO $ do
        when (x==2) $ writeIORef ref True
        return lo
  seq (ZDD.fold' f id bdd) $ do
    flag <- readIORef ref
    assertBool "unused value should be evaluated" flag

prop_fold_inSig :: Property
prop_fold_inSig =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(zdd :: ZDD o) ->
      ZDD.fold (\x lo hi -> ZDD.inSig (ZDD.SBranch x lo hi)) (ZDD.inSig . ZDD.SLeaf) zdd
      ===
      zdd

prop_fold'_inSig :: Property
prop_fold'_inSig =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(zdd :: ZDD o) ->
      ZDD.fold' (\x lo hi -> ZDD.inSig (ZDD.SBranch x lo hi)) (ZDD.inSig . ZDD.SLeaf) zdd
      ===
      zdd

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

prop_unfoldHashable_outSig :: Property
prop_unfoldHashable_outSig =
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(zdd :: ZDD o) ->
      ZDD.unfoldHashable ZDD.outSig zdd === zdd

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

prop_toGraph_fromGraph :: Property
prop_toGraph_fromGraph = do
  forAllItemOrder $ \(_ :: Proxy o) ->
    forAll arbitrary $ \(a :: ZDD o) ->
      ZDD.fromGraph (ZDD.toGraph a) === a

-- ------------------------------------------------------------------------
-- Test cases from JDD library
-- https://bitbucket.org/vahidi/jdd/src/21e103689c697fa40022294a829cab04add8ff79/src/jdd/zdd/ZDD.java

case_jdd_test_1 :: Assertion
case_jdd_test_1 = do
  let x1 = 1
      x2 = 2

      a, b, c, d, e, f, g :: ZDD ZDD.DescOrder
      a = ZDD.empty          -- {}
      b = ZDD.base           -- {{}}
      c = ZDD.change x1 b    -- {{x1}}
      d = ZDD.change x2 b    -- {{x2}}
      e = ZDD.union c d      -- {{x1}, {x2}}
      f = ZDD.union b e      -- {{}, {x1}, {x2}}
      g = ZDD.difference f c -- {{}, {x2}}

  -- directly from minatos paper, figure 9
  -- [until we find a better way to test isomorphism...]
  a @?= ZDD.Empty
  b @?= ZDD.Base
  c @?= ZDD.Branch x1 ZDD.empty ZDD.base
  d @?= ZDD.Branch x2 ZDD.empty ZDD.base
  e @?= ZDD.Branch x2 c ZDD.base
  f @?= ZDD.Branch x2 (ZDD.Branch x1 ZDD.base ZDD.base) ZDD.base
  g @?= ZDD.Branch x2 ZDD.base ZDD.base

  -- intersect
  ZDD.intersection a b @?= a
  ZDD.intersection a ZDD.base @?= a
  ZDD.intersection b b @?= b
  ZDD.intersection c e @?= c
  ZDD.intersection e f @?= e
  ZDD.intersection e g @?= d

  -- union
  ZDD.union a a @?= a
  ZDD.union b b @?= b
  ZDD.union a b @?= b
  ZDD.union g e @?= f

  -- diff
  ZDD.difference a a @?= a
  ZDD.difference b b @?= a
  ZDD.difference d c @?= d
  ZDD.difference c d @?= c
  ZDD.difference e c @?= d
  ZDD.difference e d @?= c
  ZDD.difference g b @?= d

  ZDD.difference g d @?= b
  ZDD.difference f g @?= c
  ZDD.difference f e @?= b

  -- subset0
  ZDD.subset0 x1 b @?= b
  ZDD.subset0 x2 b @?= b
  ZDD.subset0 x2 d @?= a
  ZDD.subset0 x2 e @?= c

  -- subset1
  ZDD.subset1 x1 b @?= ZDD.empty
  ZDD.subset1 x2 b @?= ZDD.empty
  ZDD.subset1 x2 d @?= b
  ZDD.subset1 x2 g @?= b
  ZDD.subset1 x1 g @?= a
  ZDD.subset1 x2 e @?= b

case_jdd_test_2 :: Assertion
case_jdd_test_2 = do
  let [a, b, c, d] = [4, 3, 2, 1]

  -- this is the exact construction sequence of Fig.14 in "Zero-suppressed BDDs and their application" by Minato
  let fig14 :: ZDD ZDD.DescOrder
      fig14 = ZDD.union z00__ z0100
        where
          z0000 = ZDD.base
          z000_ = ZDD.union z0000 (ZDD.change d z0000)
          z00__ = ZDD.union z000_ (ZDD.change c z000_)
          z0100 = ZDD.change b z0000

  -- this is the exact construction sequence of Fig.13 in "Zero-suppressed BDDs and their application" by Minato
  let fig13 :: ZDD ZDD.DescOrder
      fig13 = ZDD.intersection z0___ tmp
        where
          z___0 = ZDD.subsets (IntSet.fromList [a, b, c])
          z__0_ = ZDD.subsets (IntSet.fromList [a, b, d])
          z_0__ = ZDD.subsets (IntSet.fromList [a, c, d])
          z0___ = ZDD.subsets (IntSet.fromList [b, c, d])
          z__00 = ZDD.intersection z___0 z__0_
          tmp = ZDD.union z_0__ z__00

  fig14 @?= fig13

case_jdd_test_3 :: Assertion
case_jdd_test_3 = do
  let tmp :: ZDD ZDD.DescOrder
      tmp  = ZDD.fromListOfIntSets $ map IntSet.fromList [[2], [0,1], [1]] -- "100 011 010"
      tmp2 = ZDD.union (ZDD.singleton (IntSet.fromList [0, 1])) ZDD.base   -- union("11", base)
  -- 1. INTERSECT
  ZDD.intersection tmp tmp2 @?= ZDD.singleton (IntSet.fromList [0, 1])
  -- 2. UNION
  ZDD.union tmp tmp2 @?= ZDD.union tmp ZDD.base
  -- 3. DIFF
  ZDD.difference tmp tmp2 @?= ZDD.fromListOfIntSets (map IntSet.fromList [[1], [2]])

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

subsetOf :: IntSet -> Gen IntSet
subsetOf = liftM IntSet.fromList . sublistOf . IntSet.toList

arbitrarySmallIntMap :: Arbitrary a => Gen (IntMap a)
arbitrarySmallIntMap = do
  n <- choose (0, 12)
  liftM IntMap.fromList $ replicateM n $ do
    k <- arbitrary
    v <- arbitrary
    return (k, v)

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

zddTestGroup :: TestTree
zddTestGroup = $(testGroupGenerator)

max_stars_repo_name: msakai/haskell-decision-diagrams

max_issues_repo_name: msakai/haskell-decision-diagrams

max_forks_repo_name: msakai/haskell-decision-diagrams


text:
module TensorUtils where

import Prelude hiding ((.), id)
import System.Random
import Control.Monad

import Numeric.LinearAlgebra.Array
import Numeric.LinearAlgebra.Array.Util

import CategoricDefinitions
import Autodiff.GAD
import Autodiff.D
import Autodiff.Dual
import Autodiff.Additive

import OnesLike

{-
All tensor stuff is pretty much ad-hoc, a complete rewrite is eventually needed
-}

class ArrShow a where
    arrShow :: a -> String

instance {-# OVERLAPPABLE #-} Show a => ArrShow a where
    arrShow = show

instance  ArrShow Tensor where
    arrShow = formatFixed 2

instance (ArrShow a, ArrShow b) => ArrShow (a, b) where
    arrShow (a, b) = let ls = replicate 10 '-' ++ "\n"
                     in ls ++ arrShow a ++ ", \n" ++ arrShow b ++ "\n" ++ ls

sh :: ArrShow a => a -> IO ()
sh a = putStr $ arrShow a ++ "\n"

infixl 8 −|
--(−|) :: Name → [Array Double ] → Array Double
(−|) = index

axes :: Tensor -> String -> (String, [Int])
axes t axesNames = let toSum = filter (\idx -> head (iName idx) `elem` axesNames) (dims t)
                   in (map (head . iName) toSum, map iDim toSum)

class TensorContractable k where
    sumAxes :: String -> Tensor `k` Tensor

instance TensorContractable (->) where
    sumAxes axesNames = \t -> let (names, lengths) = axes t axesNames
                                  u = lengths # repeat 1 ! names
                              in t*u

instance (Scalable k Tensor, TensorContractable k)
    => TensorContractable (GADType k) where
    sumAxes axesNames = GAD $ \t -> let (names, lengths) = axes t axesNames
                                        u = lengths # repeat 1 ! names
                                    in (t*u, scale u)

instance (TensorContractable k) => TensorContractable (DualType k) where
    sumAxes axesNames = Dual (sumAxes axesNames)

instance TensorContractable DType where
    sumAxes axesNames = D (sumAxes axesNames)

instance TensorContractable (->+) where
    sumAxes axesNames = AddFun (sumAxes axesNames)

--meanAxes :: _ => String -> Tensor -> Tensor
--meanAxes axesNames t = let (names, lengths) = axes t axesNames
--                           n = product lengths
--                           v = scalar . fromIntegral $ n
--                       in divide (sumAxes axesNames t, v)

-- rand array with given shape and axis names
randomTensor :: [Int] -> String -> IO Tensor
randomTensor xs cs = do
    rs <- replicateM (product xs) randomIO
    return $ xs # rs ! cs


text:
module Main where

--import Data.Complex
import Graphics.Gloss
import Data.Complex
import Data.Ratio
import Fractals

window :: Display
window = InWindow "Mandelbrot set" (1000, 600) (10, 10)

background :: Color
background = black

drawing :: Complex Double -> Picture
--drawing = bmp 1000 600 (((-0.1011):+0.9563)::Complex Double) 0.0000001
drawing z0 = bmp 1000 600 z0 0.0001

main :: IO ()
--main = print $ order (mandelbrot ((-1.6) :+ 0.0)) (2.0::Double) 256 (0.0 :+ 0.0)
main = print point >> (display window background $ drawing point)
    where point = (tangentPoint (9 % 17))::Complex Double

max_stars_repo_name: arekfu/fractals-haskell

max_issues_repo_name: arekfu/fractals-haskell

max_forks_repo_name: arekfu/fractals-haskell


text:
{-# LANGUAGE CPP #-}

{- |
Module      : Language.Scheme.Compiler.Types
Copyright   : Justin Ethier
Licence     : MIT (see LICENSE in the distribution)

Maintainer  : github.com/justinethier
Stability   : experimental
Portability : portable

This module contains data types used by the compiler.
-}

module Language.Scheme.Compiler.Types
    (
    -- * Data types
      CompOpts (..)
    , CompLibOpts (..)
    , defaultCompileOptions
    , HaskAST (..)
    -- * Utility functions
    , ast2Str
    , asts2Str
    , createAstFunc
    , createAstCont
    , joinL
    , moduleRuntimeVar
    , showValAST
    -- * Headers appended to output file
    , header
    , headerComment
    , headerModule
    , headerImports
    )
where
import qualified Language.Scheme.Core as LSC (version)
import Language.Scheme.Types
import qualified Language.Scheme.Util (escapeBackslashes)
import qualified Data.Array
import qualified Data.ByteString as BS
import qualified Data.Complex as DC
import qualified Data.List
import qualified Data.Map
import qualified Data.Ratio as DR

-- |A type to store options passed to compile.
--  Eventually all of this might be able to be
--  integrated into a Compile monad.
data CompOpts = CompileOptions {
    coptsThisFunc :: String,
    -- ^Immediate name to use when creating a compiled function.
    --  Presumably there is other code that is expecting
    --  to call into it.

    coptsThisFuncUseValue :: Bool,
    -- ^Whether to include the /value/ parameter in the current function

    coptsThisFuncUseArgs :: Bool,
    -- ^Whether to include the /args/ parameter in the current function

    coptsNextFunc :: Maybe String
    -- ^The name to use for the next function after the current
    --  compiler recursion is finished. For example, after compiling
    --  a block of code, the control flow would be expected to go
    --  to this function.
    }

-- |The default compiler options
defaultCompileOptions :: String -> CompOpts
defaultCompileOptions thisFunc = CompileOptions thisFunc False False Nothing

-- |Options passed to the compiler library module
data CompLibOpts = CompileLibraryOptions {
    compBlock :: String -> Maybe String -> Env
              -> [HaskAST] -> [LispVal] -> IOThrowsError [HaskAST],
    compLisp :: Env -> String -> String -> Maybe String
              -> IOThrowsError [HaskAST]
    }

-- |Runtime reference to module data structure
moduleRuntimeVar :: String
moduleRuntimeVar = " modules "

-- |Create code for a function
createAstFunc
  :: CompOpts  -- ^ Compilation options
  -> [HaskAST] -- ^ Body of the function
  -> HaskAST -- ^ Complete function code
createAstFunc (CompileOptions thisFunc useVal useArgs _) funcBody = do
  let val = if useVal then "value" else "_"
      args = if useArgs then "(Just args)" else "_"
  AstFunction thisFunc (" env cont " ++ val ++ " " ++ args ++ " ") funcBody

-- |Create code for a continutation
createAstCont
  :: CompOpts -- ^ Compilation options
  -> String -- ^ Value to send to the continuation
  -> String -- ^ Extra leading indentation (or blank string if none)
  -> HaskAST -- ^ Generated code
createAstCont (CompileOptions _ _ _ (Just nextFunc)) var indentation = do
  AstValue $ indentation ++ "  " ++ nextFunc ++ " env cont " ++ var ++ " (Just [])"
createAstCont (CompileOptions _ _ _ Nothing) var indentation = do
  AstValue $ indentation ++ "  continueEval env cont " ++ var ++ " Nothing"

--  FUTURE: is this even necessary? Would just a string be good enough?

-- |A very basic type to store a Haskell AST.
data HaskAST = AstAssignM String HaskAST
  | AstFunction {astfName :: String,
--                 astfType :: String,
                 astfArgs :: String,
                 astfCode :: [HaskAST]
                }
 | AstValue String
 | AstRef String
 | AstContinuation {astcNext :: String,
                    astcArgs :: String
                   }

-- |Generate code based on the given Haskell AST
showValAST :: HaskAST -> String
showValAST (AstAssignM var val) = "  " ++ var ++ " <- " ++ show val
showValAST (AstFunction name args code) = do
  let typeSig = "\n" ++ name ++ " :: Env -> LispVal -> LispVal -> Maybe [LispVal] -> IOThrowsError LispVal "
  let fheader = "\n" ++ name ++ args ++ " = do "
  let fbody = unwords . map (\x -> '\n' : x ) $ map showValAST code
#ifdef UseDebug
  let appendArg arg = do
        if Data.List.isInfixOf arg args
           then " ++ \" \" ++ " ++ (show arg) ++
                " ++ \" [\" ++ (show " ++ arg ++ ")" ++
                " ++ \"] \""
           else ""
  let fdebug = "\n  _ <- liftIO $ (trace (\"" ++
               name ++ "\"" ++
               (appendArg "value") ++
               (appendArg "args") ++
               ") getCPUTime)"
  typeSig ++ fheader ++ fdebug ++ fbody
#else
  typeSig ++ fheader ++ fbody
#endif
showValAST (AstValue v) = v
showValAST (AstRef v) = v
showValAST (AstContinuation nextFunc args) =
    "  continueEval env (makeCPSWArgs env cont " ++
       nextFunc ++ " " ++ args ++ ") (Nil \"\") Nothing "

instance Show HaskAST where show = showValAST

-- |A utility function to join list members together
joinL
  :: forall a. [[a]] -- ^ Original list-of-lists
  -> [a] -- ^ Separator
  -> [a] -- ^ Joined list
joinL ls sep = Data.List.intercalate sep ls

-- |Convert abstract syntax tree to a string
ast2Str :: LispVal -> String
ast2Str (String s) = "String " ++ show s
ast2Str (Char c) = "Char " ++ show c
ast2Str (Atom a) = "Atom " ++ show a
ast2Str (Number n) = "Number (" ++ show n ++ ")"
ast2Str (Complex c) = "Complex $ (" ++ (show $ DC.realPart c) ++ ") :+ (" ++ (show $ DC.imagPart c) ++ ")"
ast2Str (Rational r) = "Rational $ (" ++ (show $ DR.numerator r) ++ ") % (" ++ (show $ DR.denominator r) ++ ")"
ast2Str (Float f) = "Float (" ++ show f ++ ")"
ast2Str (Bool True) = "Bool True"
ast2Str (Bool False) = "Bool False"
ast2Str (HashTable ht) = do
 let ls = Data.Map.toList ht
     conv (a, b) = "(" ++ ast2Str a ++ "," ++ ast2Str b ++ ")"
 "HashTable $ Data.Map.fromList $ [" ++ joinL (map conv ls) "," ++ "]"
ast2Str (Vector v) = do
  let ls = Data.Array.elems v
      size = (length ls) - 1
  "Vector (listArray (0, " ++ show size ++ ")" ++ "[" ++ joinL (map ast2Str ls) "," ++ "])"
ast2Str (ByteVector bv) = do
  let ls = BS.unpack bv
  "ByteVector ( BS.pack " ++ "[" ++ joinL (map show ls) "," ++ "])"
ast2Str (List ls) = "List [" ++ joinL (map ast2Str ls) "," ++ "]"
ast2Str (DottedList ls l) =
  "DottedList [" ++ joinL (map ast2Str ls) "," ++ "] $ " ++ ast2Str l
ast2Str l = show l -- Error?

-- |Convert a list of abstract syntax trees to a list of strings
asts2Str :: [LispVal] -> String
asts2Str ls = do
    "[" ++ (joinL (map ast2Str ls) ",") ++ "]"

-- |Header comment used at the top of a Haskell program generated
--  by the compiler
headerComment:: [String]
headerComment = [
   "--"
 , "-- This file was automatically generated by the husk scheme compiler (huskc)"
 , "--"
 , "--  http://justinethier.github.io/husk-scheme "
 , "--  (c) 2010 Justin Ethier "
 , "--  Version " ++ LSC.version
 , "--"]

-- |Main module used in a compiled Haskell program
headerModule :: [String]
headerModule = ["module Main where "]

-- |Imports used for a compiled program
headerImports :: [String]
headerImports = [
   "Language.Scheme.Core "
 , "Language.Scheme.Numerical "
 , "Language.Scheme.Macro "
 , "Language.Scheme.Primitives "
 , "Language.Scheme.Types     -- Scheme data types "
 , "Language.Scheme.Variables -- Scheme variable operations "
 , "Control.Monad.Error "
 , "Data.Array "
 , " qualified Data.ByteString as BS "
 , "Data.Complex "
 , " qualified Data.Map "
 , "Data.Ratio "
 , "Data.Word "
 , "System.IO "
#ifdef UseDebug
 , "System.CPUTime "
 , "Debug.Trace "
#endif
 ]

-- |Block of code used in the header of a Haskell program
--  generated by the compiler.
header :: String -> Bool -> String -> [String]
header filepath useCompiledLibs langRev = do
  let env = if useCompiledLibs
            then "primitiveBindings"
            else case langRev of
                   "7" -> "r7rsEnv"
                   _ -> "r5rsEnv"
      initSrfi55 =
        case langRev of
          "7" -> []
          _ -> [ "exec55_3 env cont _ _ = do "
               , "  liftIO $ registerExtensions env getDataFileName' "
               , "  continueEval env (makeCPSWArgs env cont exec []) (Nil \"\") Nothing"]
  [ " "
    , " "
    , "-- |Get variable at runtime "
    , "getRTVar env var = do "
    , "  v <- getVar env var "
    , "  return $ case v of "
    , "    List _ -> Pointer var env "
    , "    DottedList _ _ -> Pointer var env "
    , "    String _ -> Pointer var env "
    , "    Vector _ -> Pointer var env "
    , "    ByteVector _ -> Pointer var env "
    , "    HashTable _ -> Pointer var env "
    , "    _ -> v "
    , " "
-- TODO:  this is just a temporary function until calls to continueEval can be purged from the compiler
    , "continueEval' env cont value = continueEval env cont value Nothing "
    , " "
    , "applyWrapper env cont (Nil _) (Just (a:as))  = do "
    , "  apply cont a as "
    , " "
    , "applyWrapper env cont value (Just (a:as))  = do "
    , "  apply cont a $ as ++ [value] "
    , " "
    , "getDataFileName' :: FilePath -> IO FilePath "
    , "getDataFileName' name = return $ \"" ++ (Language.Scheme.Util.escapeBackslashes filepath) ++ "\" ++ name "
    , " "]
    ++ initSrfi55 ++
    [ " "
    , "main :: IO () "
    , "main = do "
    , "  env <- " ++ env ++ " "
    , "  result <- (runIOThrows $ liftM show $ hsInit env (makeNullContinuation env) (Nil \"\") Nothing) "
    , "  case result of "
    , "    Just errMsg -> putStrLn errMsg "
    , "    _ -> return () "
    , " "
    , "hsInit env cont _ _ = do "
    , "  _ <- defineVar env \"" ++ moduleRuntimeVar ++ "\" $ HashTable $ Data.Map.fromList [] "
    , "  run env cont (Nil \"\") (Just [])"
    , " "]

max_stars_repo_name: svenkeidel/husk-scheme

max_issues_repo_name: haskell-lisp/husk-scheme

max_forks_repo_name: haskell-lisp/husk-scheme


text:
module HScheme.Data where

import Control.Monad.Except
import Data.Complex ( Complex )
import Text.ParserCombinators.Parsec ( ParseError )
import qualified Data.Vector as V

data Value = Atom String
    | List [Value]
    | DottedList [Value] Value
    | Vector (V.Vector Value)
    | Number Integer
    | Float Double
    | Ratio Rational
    | Complex (Complex Double)
    | String String
    | Character Char
    | Bool Bool
    deriving (Eq)

instance Show Value where
    show (String s) = "\"" ++ s ++ "\""
    show (Character c) = "#\\" ++ [c]
    show (Atom name) = name
    show (Number n) = show n
    show (Float n) = show n
    show (Ratio n) = show n
    show (Complex n) = show n
    show (Bool True) = "#t"
    show (Bool False) = "#f"
    show (List xs) = "(" ++ unwordsList xs ++ ")"
    show (DottedList head tail) = "(" ++ unwordsList head ++ " . " ++ show tail ++ ")"
    show (Vector vs) = "#(" ++ unwordsList (V.toList vs) ++ ")"

unwordsList :: [Value] -> String
unwordsList = unwords . map show

{-
    -- Error handling --
-}

data Error = NumArgs Integer [Value]
    | TypeMismatch String Value
    | Parser ParseError
    | BadSpecialForm String Value
    | NotFunction String String
    | UnboundVar String String
    | Default String

instance Show Error where
    show (NumArgs expected found)       = "Expected " ++ show expected
                                        ++ " args; found values " ++ unwordsList found
    show (TypeMismatch expected found)  = "Invalid Type: expected " ++ expected
                                        ++ ", found " ++ show found
    show (Parser err)                   = "Parse error at " ++ show err
    show (BadSpecialForm message form)  = message ++ ": " ++ show form
    show (NotFunction message func)     = message ++ ": " ++ show func
    show (UnboundVar message varname)   = message ++ ": " ++ varname
    show (Default message)              = message

type ThrowsErr = Either Error

trapErr :: (MonadError a m, Show a) => m String -> m String
trapErr action = catchError action (return . show)

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

max_stars_repo_name: dvdvgt/hScheme

max_issues_repo_name: dvdvgt/hScheme

max_forks_repo_name: dvdvgt/hScheme


text:
{-# LANGUAGE CPP #-}
{-# OPTIONS_GHC -fno-warn-unused-imports #-}

module CGX.Prelude
  ( module X
#if __GLASGOW_HASKELL__ <= 800
  , fromLeft
  , fromRight
#endif
  , mk
  , ps
  , psu
  ) where

import qualified Data.Text as Text
import qualified Prelude

import Prelude as X
  ( Bounded(..)
  , Enum(..)
  , iterate
  , maxBound
  )

import Control.Applicative as X
  ( Applicative(..)
  , (<$)
  , (<$>)
  , (<*)
  , (*>)
  , pure
  )

import Control.Monad as X
  ( Functor
  , Monad
  , forever
  , fmap
  , liftM
  , return
  , when
  , (>>=)
  , (=<<)
  , (>>)
  )

import Control.Monad.State as X
  ( State
  , get
  , put
  , runState
  )

import Control.Exception.Base as X
  ( SomeException(..)
  , handle
  )

import Control.Concurrent as X

import Debug.Trace as X
  ( trace
  , traceIO
  , traceShow
  )

import Safe as X
  ( headMay
  , abort
  )

import Data.Eq as X
import Data.Ord as X
import Data.Monoid as X
import Data.Traversable as X
import Data.Foldable as X hiding
  ( foldr1
  , foldl1
  , maximum
  , maximumBy
  , minimum
  , minimumBy
  )

import Data.Int as X
import Data.Bits as X
import Data.Word as X
import Data.Bool as X hiding (bool)
import Data.Char as X (Char)
import Data.Maybe as X hiding (fromJust)
import Data.Either as X
--import Data.Complex as X

import Data.Function as X
  ( id
  , const
  , (.)
  , ($)
  , flip
  , fix
  , on
  )

import Data.Tuple as X
import Data.List as X
  ( filter
  , intersect
  , iterate
  , reverse
  , take
  , takeWhile
  , (\\)
  , (++)
  )

import Data.Text as X
  ( Text
  , pack
  )
import qualified Data.Text.Lazy
import qualified Data.Text.IO
import Data.String.Conv as X
  ( toS
  )
import Data.IntMap as X (IntMap)

import Data.Tagged as X

import GHC.IO as X (IO)
import GHC.Num as X
import GHC.Real as X
import GHC.Float as X
import GHC.Generics as X
import GHC.Show as X
import GHC.Exts as X
  ( Constraint
  , Ptr
  , FunPtr
  , the
  )

#if __GLASGOW_HASKELL__ <= 800
fromLeft :: a -> Either a b -> a
fromLeft _ (Left a) = a
fromLeft a _        = a

fromRight :: b -> Either a b -> b
fromRight _ (Right b) = b
fromRight b _         = b
#endif

mk :: a -> Tagged b a
mk = X.Tagged

psu :: Show a => Tagged b a -> Text
psu = Text.pack . show . X.untag

ps :: Show a => a -> Text
ps = Text.pack . show

max_stars_repo_name: crappygraphix/junks

max_issues_repo_name: crappygraphix/junks

max_forks_repo_name: crappygraphix/junks