By the way, the variable "b" which "b{tv a1hL}" seems to be referring
to is on the last line of vecQ. Help is much appreciated!
Thanks,
Frederik
On Fri, Dec 16, 2005 at 09:54:28PM +0000, Frederik Eaton wrote:
> $ ghc --make Vector.hs -fth
> Chasing modules from: Vector.hs
> [1 of 2] Skipping Fu.Prepose ( Fu/Prepose.hs, Fu/Prepose.o )
> [2 of 2] Compiling Fu.Vector ( Vector.hs, Vector.o )
> ghc-6.5.20051208: panic! (the `impossible' happened, GHC version
> 6.5.20051208):
> Failed binder lookup: b{tv a1hL}
>
> Please report it as a compiler bug to [email protected],
> or http://sourceforge.net/projects/ghc/.
> {-# OPTIONS -fglasgow-exts -fallow-undecidable-instances #-}
> module Fu.Vector where
>
> import System.IO.Unsafe
> import Control.Exception
> import Data.Ix
> import Data.Array.IO
> import qualified Data.Array as Array
> import Data.Array hiding ((!))
> import Fu.Prepose hiding (__)
> import Prelude hiding (sum)
> import Language.Haskell.TH
> import Language.Haskell.TH.Syntax
> import qualified Data.List
>
> __ = __
>
> data L n = L Int deriving (Show, Eq, Ord)
> unL l@(L x) = check l x
>
> checkBounds :: ReflectNum n => L n -> a -> a
> checkBounds (l::L n) y =
> if l < (L 0) || l > maxBound then
> error ("Value "++show l++" out of bounds for L "++(show $ (maxBound
> :: L n)))
> else y
>
> -- need Ord too?
> class (Bounded a, Ix a, Eq a, Show a) => IxB a
>
> instance (Bounded a, Ix a, Eq a, Show a) => IxB a
>
> --instance (IxB a, IxB b) => IxB (a,b)
>
> --instance IxB ()
>
> instance (Bounded a, Bounded b) => Bounded (Either a b) where
> minBound = Left minBound
> maxBound = Right maxBound
>
> instance (IxB a, IxB b) => Ix (Either a b) where
> range (Left x, Left y) = map Left $ range (x,y)
> range (Right x, Right y) = map Right $ range (x,y)
> range (Left x, Right y) =
> (map Left $ range (x,maxBound))++
> (map Right $ range (minBound,y))
> range (Right _, Left _) = error "range endpoints out of order"
> index (Left x, Left y) (Left z) = index (x,y) z
> index (Left x, Left y) (Right z) = error "out of bounds"
> index (Right x, Right y) (Right z) = index (x,y) z
> index (Right x, Right y) (Left z) = error "out of bounds"
> index (Left x, Right y) (Left z) = index (x,maxBound) z
> index (Left x, Right y) (Right z) =
> rangeSize (x,maxBound) + index (minBound,y) z
> index (Right _, Left _) _ = error "range endpoints out of order"
> inRange (x,y) z = (x<=z) && (z<=y)
>
> --instance (IxB a, IxB b) => IxB (Either a b)
>
> domain :: IxB a => [a]
> domain = range (minBound, maxBound)
>
> instance ReflectNum n => Bounded (L n) where
> minBound = L 0
> maxBound = L $ reflectNum (__ :: n) - 1
>
> instance ReflectNum n => Check (L n) where
> check l = checkBounds l
>
> class Check a where
> check :: a -> (b -> b)
>
> instance ReflectNum n => Ix (L n) where
> index (a, b) (c) = index (unL a, unL b) (unL c)
> range (a, b) = map L $ range (unL a, unL b)
> inRange (a, b) c = inRange (unL a, unL b) (unL c)
> rangeSize (a, b) = rangeSize (unL a, unL b)
>
> -- needed?
> type One = Succ Zero
>
> -- needed?
> -- singleton domain
> type S = L One
>
> -- needed?
> -- empty domain
> type Z = L Zero
>
> -- like Array but uses IxB instead of Ix
> newtype Vector k i = Vector (Array.Array i k) deriving (Eq, Show)
>
> type Matrix k a b = Vector k (a,b)
>
> boundaries :: Bounded a => (a,a)
> boundaries = (minBound, maxBound)
>
> vector :: IxB a => (a -> k) -> Vector k a
> vector f = Vector
> (array boundaries (map (\i -> (i, f i)) domain))
>
> (!) :: IxB a => Vector k a -> a -> k
> (!) (Vector ax) i = (Array.!) ax i
>
> slice :: (IxB a, IxB b) => Vector k a -> (b -> a) -> Vector k b
> slice v f = vector ((v!) . f)
>
> updateArray ax i f = do
> e <- readArray ax i
> writeArray ax i (f e)
>
> margin :: (IxB a, IxB b, Num k) => Vector k a -> (a -> b) -> Vector k b
> margin (Vector v) f = unsafePerformIO $ (do
> (ax::IOArray a k) <- newArray (minBound, maxBound) 0
> mapM_ (\ (i,x) -> updateArray ax (f i) (+x)) (assocs v)
> ac <- getAssocs ax
> return $ Vector $ array boundaries ac
> )
>
> -- terminal map
> terminal :: a -> ()
> terminal x = ()
>
>
> sum :: (IxB a, Num k) => Vector k a -> k
> sum v = margin v (const ()) ! ()
>
> diag :: a -> (a,a)
> diag x = (x,x)
>
> trace m = sum $ slice m diag
>
> -- initial map
> initial :: Z -> a
> initial _ = undefined
>
> (*>) :: (Num k, IxB a, IxB b, IxB c) => Matrix k a b -> Matrix k b c ->
> Matrix k a c
> (*>) m n =
> vector (\ (i,k) -> Data.List.sum $ map (\j -> m ! (i,j) + n ! (j,k))
> domain)
>
> mapv :: (IxB a) => (k -> l) -> Vector k a -> Vector l a
> mapv f v = vector (\i -> f (v!i))
>
> mapv2 :: (IxB a) => (k -> l -> m) -> Vector k a -> Vector l a -> Vector m a
> mapv2 f u v = vector (\i -> f (u!i) (v!i))
>
> -- element-wise operations
> instance (Num k, IxB a) => Num (Vector k a) where
> (+) = mapv2 (+)
> (-) = mapv2 (-)
> (*) = mapv2 (*)
> negate = mapv negate
> signum = mapv signum
> abs = mapv abs
> fromInteger n = (fromInteger n) .* ones
>
> -- should probably be able to do this in constant time
> vec :: (Num k, IxB a, IxB b) => Matrix k a b -> Matrix k (a,b) S
> vec m = slice m (\ ((i,j),_) -> (i,j))
>
> -- better name?
> inv :: (RealFrac k, IxB a) => Matrix k a a -> Matrix k a a
> inv = undefined
>
> -- need least squares operator: />?
>
> det :: (Num k, IxB a) => Matrix k a a -> k
> det = undefined
>
> -- should reorder type args so this can be a Functor
> mapm :: (IxB a, IxB b) => (k -> l) -> Matrix k a b -> Matrix l a b
> -- in fact should it be a Monad? could be similar but types would need to
> change
> mapm = mapv
>
> vlen :: IxB a => Vector k a -> Int
> vlen (v::Vector k a) = rangeSize $ (boundaries :: (a,a))
>
> unif :: (Fractional k, IxB a) => Vector k a
> unif = v
> where
> p = 1/(fromIntegral $ vlen v)
> v = vector (const p)
>
> ones :: (Num k, IxB a) => Vector k a
> ones = vector (const 1)
>
> zeros :: (Num k, IxB a) => Vector k a
> zeros = vector (const 0)
>
> (.*) :: (Num k, IxB a) => k -> Vector k a -> Vector k a
> (.*) x v = vector (\i -> x*v!i)
>
> normalize :: Vector k a -> Vector k a
> normalize = undefined
>
> concat :: (IxB a, IxB b) => Vector k a -> Vector k b -> Vector k (Either a b)
> concat u v = vector $ \i ->
> case i of
> Left j -> u!j
> Right k -> v!k
>
> kronecker :: (IxB a, IxB b, Num k) => Vector k a -> Vector k b -> Vector k
> (a,b)
> kronecker u v = vector (\(i,j) -> u!i * v!j)
>
> -- range
> (|>) :: (IxB a, IxB b) => Vector k b -> Vector b a -> Vector k a
> (|>) u ix = vector (\i -> u ! (ix ! i))
>
> (||>) :: (IxB a, IxB b, IxB c, IxB d) => Matrix k c d -> (Vector c a, Vector
> d b) -> Matrix k a b
> (||>) m (ix0, ix1) = vector (\ (i,j) -> m ! (ix0 ! i, ix1 ! j))
>
> -- modify just the 'true' entries of the argument
> subsetModify :: Vector k b -> (b -> Bool) -> (forall a . Vector k a -> Vector
> k a) -> Vector k b
> subsetModify = undefined
>
> filterModify :: Vector k b -> (k -> b -> Bool) -> (forall a . Vector k a ->
> Vector k a) -> Vector k b
> filterModify = undefined
>
> -- need to inject values into a vector...
>
> -- permutations, size is n!
> data Perm n = Perm (Vector (L n) (L n))
>
> --det :: (Num k, IxB a) => Matrix a a k -> k
>
> listVec :: [k] -> (forall a . IxB a => Vector k a -> w) -> w
> listVec (l::[k]) f =
> let n = length l in
> reifyIntegral n (\ (_::rn) ->
> f (listToVec l::Vector k (L rn)))
>
> vecQ :: Lift k => [k] -> Q Exp
> vecQ (l::[k]) =
> let n = length l in
> reifyIntegral n (\ (x::rn) ->
> let r = T x
> rt = [| r |] in
> [| let (_::b) = $(rt) in listToVec l::Vector k (L b) |])
>
> listToVec :: IxB a => [k] -> Vector k a
> listToVec l::Vector k a = (v::Vector k a)
> where
> v = assert (length l == vlen v) $
> Vector $ listArray (boundaries :: (a,a)) l
>
> --listMat :: [[k]] -> (forall a b . (IxB a, IxB b) => Matrix k (a,b) -> w) ->
> w
>
> main = undefined
>
> {-# OPTIONS -fglasgow-exts -fallow-undecidable-instances #-}
>
> -- Big comment:
>
> -- Note that the terms "reify" and "reflect" may seem to be used in a way
> -- which is the reverse of what one might intuit! Here, the
> -- representation of a number in type-level binary is considered "real".
> -- Thus to turn a value into a type-level term is to "reify" it; to go
> -- the other way is to "reflect".
>
> module Fu.Prepose where
>
> import System.IO.Unsafe (unsafePerformIO)
> import Control.Exception (handle, handleJust, errorCalls)
> import Foreign.Marshal.Utils (with, new)
> import Foreign.Marshal.Array (peekArray, pokeArray)
> import Foreign.Marshal.Alloc (alloca)
> import Foreign.Ptr (Ptr, castPtr)
> import Foreign.Storable (Storable(sizeOf, peek))
> import Foreign.C.Types (CChar)
> import Foreign.StablePtr (StablePtr, newStablePtr,
> deRefStablePtr, freeStablePtr)
> import Data.Bits (Bits(..))
> import Prelude hiding (getLine)
>
> import Language.Haskell.TH hiding (reify)
> import Language.Haskell.TH.Syntax hiding (reify)
>
> import Debug.Trace
>
> newtype Modulus s a = Modulus a deriving (Eq, Show)
> newtype M s a = M a deriving (Eq, Show)
>
> unM :: M s a -> a
> unM (M a) = a
>
> __ = __
>
> class Modular s a | s -> a where modulus :: s -> a
>
> normalize :: (Modular s a, Integral a) => a -> M s a
> normalize a :: M s a = M (mod a (modulus (__ :: s)))
>
> instance (Modular s a, Integral a) => Num (M s a) where
> M a + M b = normalize (a + b)
> M a - M b = normalize (a - b)
> M a * M b = normalize (a * b)
> negate (M a) = normalize (negate a)
> fromInteger i = normalize (fromInteger i)
> signum = error "Modular numbers are not signed"
> abs = error "Modular numbers are not signed"
>
> aaa (M a :: M s1 a) (M b :: M s2 a) = normalize (a + b) :: M s1 a
> where _ = [__::s1, __ :: s2] -- make sure input moduli are equal
>
> withModulus :: a -> (forall s. Modular s a => s -> w) -> w
>
> data Zero; data Twice s; data Succ s; data Pred s
>
> data Positive s; data Negative s; data TwiceSucc s
>
> ----------------------------------------------------------------
>
> class GetType s where
> getType :: s -> Type
>
> instance GetType (Zero) where
> getType _ = ConT $ mkName "Zero"
>
> instance GetType s => GetType (Twice s) where
> getType _ = AppT (ConT $ mkName "Twice") (getType (__ :: s))
>
> instance GetType s => GetType (Succ s) where
> getType _ = AppT (ConT $ mkName "Succ") (getType (__ :: s))
>
> instance GetType s => GetType (Pred s) where
> getType _ = AppT (ConT $ mkName "Pred") (getType (__ :: s))
>
> instance GetType s => GetType (Positive s) where
> getType _ = AppT (ConT $ mkName "Positive") (getType (__ :: s))
>
> instance GetType s => GetType (Negative s) where
> getType _ = AppT (ConT $ mkName "Negative") (getType (__ :: s))
>
> instance GetType s => GetType (TwiceSucc s) where
> getType _ = AppT (ConT $ mkName "TwiceSucc") (getType (__ :: s))
>
> data T s = T s
>
> instance GetType s => Lift (T s) where
> lift _ = return $ SigE undE (getType (__ :: s))
>
> undE = VarE $ mkName "undefined"
>
> ----------------------------------------------------------------
>
> class GetType s => ReflectUnsigned s where reflectUnsigned :: Num a => s -> a
> instance ReflectUnsigned Zero where
> reflectUnsigned _ = 0
> instance ReflectUnsigned s => ReflectUnsigned (Twice s) where
> reflectUnsigned _ = reflectUnsigned (__ :: s) * 2
> instance ReflectUnsigned s => ReflectUnsigned (TwiceSucc s) where
> reflectUnsigned _ = reflectUnsigned (__ :: s) * 2 + 1
>
> instance ReflectUnsigned s => ReflectNum (Positive s) where
> reflectNum _ = reflectUnsigned (__ :: s)
> instance ReflectUnsigned s => ReflectNum (Negative s) where
> reflectNum _ = -1 - reflectUnsigned (__ :: s)
>
>
> class GetType s => ReflectNum s where reflectNum :: Num a => s -> a
> instance ReflectNum Zero where
> reflectNum _ = 0
> instance ReflectNum s => ReflectNum (Twice s) where
> reflectNum _ = reflectNum (__ :: s) * 2
> instance ReflectNum s => ReflectNum (Succ s) where
> reflectNum _ = reflectNum (__ :: s) + 1
> instance ReflectNum s => ReflectNum (Pred s) where
> reflectNum _ = reflectNum (__ :: s) - 1
>
>
> reifyIntegral :: Integral a =>
> a -> (forall s. ReflectNum s => s -> w) -> w
> reifyIntegral i k = --(trace $ "reifyIntegral "++show i) $
> case quotRem i 2 of
> (0, 0) -> k (__ :: Zero)
> (0, 1) -> k (__ :: Succ Zero)
> (0,- 1) -> k (__ :: Pred Zero)
> (j, 0) -> reifyIntegral j (\ (_ :: s) -> k (__ :: Twice s))
> (j, 1) -> reifyIntegral j (\ (_ :: s) -> k (__ :: Succ (Twice s)))
> (j,- 1) -> reifyIntegral j (\ (_ :: s) -> k (__ :: Pred (Twice s)))
>
>
> data ModulusNum s a
>
> instance (ReflectNum s, Num a) =>
> Modular (ModulusNum s a) a where
> modulus _ = reflectNum (__ :: s)
>
>
> withIntegralModulus :: Integral a =>
> a -> (forall s. Modular s a => s -> w) -> w
> withIntegralModulus i k =
> reifyIntegral i (\(_ :: s) -> k (__ :: ModulusNum s a))
>
>
> data Nil; data Cons s ss
>
> class ReflectNums ss where reflectNums :: Num a => ss -> [a]
> instance ReflectNums Nil where
> reflectNums _ = []
> instance (ReflectNum s, ReflectNums ss) =>
> ReflectNums (Cons s ss) where
> reflectNums _ = reflectNum (__ :: s) : reflectNums (__ :: ss)
>
> reifyIntegrals :: Integral a =>
> [a] -> (forall ss. ReflectNums ss => ss -> w) -> w
> reifyIntegrals [] k = k (__ :: Nil)
> reifyIntegrals (i:ii) k = reifyIntegral i (\(_ :: s) ->
> reifyIntegrals ii (\(_ :: ss) ->
> k (__ :: Cons s ss)))
>
>
> type Byte = CChar
>
> data Store s a
>
> class ReflectStorable s where
> reflectStorable :: Storable a => s a -> a
> instance ReflectNums s => ReflectStorable (Store s) where
> reflectStorable _ = unsafePerformIO
> $ alloca
> $ \p -> do pokeArray (castPtr p) bytes
> peek p
> where bytes = reflectNums (__ :: s) :: [Byte]
>
> reifyStorable :: Storable a =>
> a -> (forall s. ReflectStorable s => s a -> w) -> w
> reifyStorable a k =
> reifyIntegrals (bytes :: [Byte]) (\(_ :: s) -> k (__ :: Store s a))
> where bytes = unsafePerformIO
> $ with a (peekArray (sizeOf a) . castPtr)
>
>
> class Reflect s a | s -> a where reflect :: s -> a
>
> data Stable (s :: * -> *) a
>
> instance ReflectStorable s => Reflect (Stable s a) a where
> reflect = unsafePerformIO
> $ do a <- deRefStablePtr p
> return (const a)
> where p = reflectStorable (__ :: s p)
>
> reify :: a -> (forall s. Reflect s a => s -> w) -> w
> reify (a :: a) k = unsafePerformIO
> $ do p <- newStablePtr a
> reifyStorable p (\(_ :: s p) ->
> k' (__ :: Stable s a))
> where k' s = return (k s)
>
>
> data ModulusAny s
>
> instance Reflect s a => Modular (ModulusAny s) a where
> modulus _ = reflect (__ :: s)
>
> withModulus a k = reify a (\(_ :: s) -> k (__ :: ModulusAny s))
>
>
> withIntegralModulus' :: Integral a =>
> a -> (forall s. Modular s a => M s w) -> w
> withIntegralModulus' (i::a) k :: w =
> reifyIntegral i ( \(_ :: t) ->
> unM (k :: M (ModulusNum t a) w))
>
> test4' :: (Modular s a, Integral a) => M s a
> test4' = 3*3 + 5*5
>
> test4 = withIntegralModulus' 4 test4'
>
> data Even p q u v a = E a a deriving (Eq, Show)
>
> normalizeEven :: ( ReflectNum p, ReflectNum q, Integral a,
> Bits a) => a -> a -> Even p q u v a
> normalizeEven a b :: Even p q u v a =
> E (a .&. (shiftL 1 (reflectNum (__ :: p)) - 1)) -- $a \bmod 2^p$
> (mod b (reflectNum (__ :: q))) -- $b \bmod q$
>
> instance ( ReflectNum p, ReflectNum q,
> ReflectNum u, ReflectNum v,
> Integral a, Bits a) => Num (Even p q u v a) where
> E a1 b1 + E a2 b2 = normalizeEven (a1 + a2) (b1 + b2)
> {-"\raisebox{0pt}[2.5ex][0pt]{$\vdots$}"-}
>
>
> E a1 b1 - E a2 b2 = normalizeEven (a1 - a2) (b1 - b2)
> E a1 b1 * E a2 b2 = normalizeEven (a1 * a2) (b1 * b2)
> negate (E a b) = normalizeEven (-a) (-b)
> fromInteger i = normalizeEven (fromInteger i)
> (fromInteger i)
> signum = error "Modular numbers are not signed"
> abs = error "Modular numbers are not signed"
>
> gcd' x 0 = []
> gcd' x y = (x `div` y) : (gcd' y (x `mod` y))
>
> -- chinese remainder theorem
> crt a b =
> let (c,d) = snd $ foldl (\ ((a,b),(a',b')) p ->
> ((-(a*p+a'),b*p+b'),(-a,b)))
> ((1,0),(0,1)) (gcd' a b)
> in (-c*b,d*a)
>
> withIntegralModulus'' :: (Integral a, Bits a) =>
> a -> (forall s. Num (s a) => s a) -> a
> withIntegralModulus'' (i::a) k = case factor 0 i of
> (0, i) -> withIntegralModulus' i k -- odd case
> (p, q) -> let (u, v) = crt (2^p) q in -- even case: $i
> = 2^p q$
> -- trace ("(u,v)="++show(u,v)) $
> reifyIntegral p (\(_::p ) ->
> reifyIntegral q (\(_::q ) ->
> reifyIntegral u (\(_::u ) ->
> reifyIntegral v (\(_::v ) ->
> unEven (k :: Even p q u v a)))))
>
> factor :: (Num p, Integral q) => p -> q -> (p, q)
> factor p i = case quotRem i 2 of
> (0, 0) -> (0, 0) -- just zero
> (j, 0) -> factor (p+1) j -- accumulate powers of two
> _ -> (p, i) -- not even
>
> unEven ::( ReflectNum p, ReflectNum q, ReflectNum u,
> ReflectNum v, Integral a, Bits a) => Even p q u v a -> a
> unEven (E a b :: Even p q u v a) = -- trace "unEven" $
> mod (a * (reflectNum (__ :: u)) + b * (reflectNum (__ :: v)))
> (shiftL (reflectNum (__ :: q)) (reflectNum (__ :: p)))
>
>
> test5 :: Num (s a) => s a
> test5 = 1000*1000 + 513*513
>
> test5' = withIntegralModulus'' 1279 test5 :: Integer
> test5'' = withIntegralModulus'' 1280 test5 :: Integer
>
> _______________________________________________
> Glasgow-haskell-bugs mailing list
> [email protected]
> http://www.haskell.org/mailman/listinfo/glasgow-haskell-bugs
_______________________________________________
Glasgow-haskell-bugs mailing list
[email protected]
http://www.haskell.org/mailman/listinfo/glasgow-haskell-bugs