$ 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 glasgow-haskell-bugs@haskell.org, 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 Glasgow-haskell-bugs@haskell.org http://www.haskell.org/mailman/listinfo/glasgow-haskell-bugs