lol-cpp (empty) → 0.0.0.1
raw patch · 22 files changed
+3058/−0 lines, 22 filesdep +DRBGdep +MonadRandomdep +arithmoisetup-changed
Dependencies added: DRBG, MonadRandom, arithmoi, base, bytestring, constraints, containers, crypto-api, data-default, deepseq, lol, lol-benches, lol-cpp, lol-tests, monadcryptorandom, mtl, numeric-prelude, protocol-buffers, protocol-buffers-descriptor, random, reflection, repa, singletons, tagged-transformer, template-haskell, th-desugar, transformers, vector, vector-th-unbox
Files
- CHANGES.md +6/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP.hs +454/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/Backend.hs +337/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/Extension.hs +164/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/Instances.hs +30/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/common.cpp +27/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/common.h +20/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/crt.cpp +598/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/g.cpp +273/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/l.cpp +98/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/mul.cpp +35/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/norm.cpp +80/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/random.cpp +64/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/tensor.h +67/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/types.h +170/−0
- Crypto/Lol/Cyclotomic/Tensor/CPP/zq.cpp +50/−0
- LICENSE +339/−0
- README +5/−0
- Setup.hs +2/−0
- benchmarks/BenchCPPMain.hs +83/−0
- lol-cpp.cabal +135/−0
- tests/TestCPPMain.hs +21/−0
+ CHANGES.md view
@@ -0,0 +1,6 @@+Changelog for lol-cpp project+================================++0.0.0.1+----+ * Initial split from lol package. Now requires GHC >= 8.0.
+ Crypto/Lol/Cyclotomic/Tensor/CPP.hs view
@@ -0,0 +1,454 @@+{-|+Module : Crypto.Lol.Cyclotomic.Tensor.CPP+Description : Wrapper for a C++ implementation of the 'Tensor' interface.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX++Wrapper for a C++ implementation of the 'Tensor' interface.+-}++{-# LANGUAGE ConstraintKinds #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE GADTs #-}+{-# LANGUAGE GeneralizedNewtypeDeriving #-}+{-# LANGUAGE InstanceSigs #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE PolyKinds #-}+{-# LANGUAGE RankNTypes #-}+{-# LANGUAGE RebindableSyntax #-}+{-# LANGUAGE RoleAnnotations #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE StandaloneDeriving #-}+{-# LANGUAGE TupleSections #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeOperators #-}+{-# LANGUAGE UndecidableInstances #-}++module Crypto.Lol.Cyclotomic.Tensor.CPP (CT) where++import Algebra.Additive as Additive (C)+import Algebra.Module as Module (C)+import Algebra.ZeroTestable as ZeroTestable (C)++import Control.Applicative hiding ((*>))+import Control.Arrow ((***))+import Control.DeepSeq+import Control.Monad.Except+import Control.Monad.Identity (Identity (..), runIdentity)+import Control.Monad.Random+import Control.Monad.Trans as T (lift)++import Data.Coerce+import Data.Constraint hiding ((***))+import Data.Int+import Data.Maybe+import Data.Traversable as T+import Data.Vector.Generic as V (fromList, toList, unzip)+import Data.Vector.Storable as SV (Vector, convert, foldl',+ fromList, generate,+ length, map, replicate,+ replicateM, thaw, thaw, toList,+ unsafeFreeze,+ unsafeWith, zipWith, (!))+import Data.Vector.Storable.Mutable as SM hiding (replicate)++import Foreign.Marshal.Utils (with)+import Foreign.Ptr++import Crypto.Lol.CRTrans+import Crypto.Lol.Cyclotomic.Tensor+import Crypto.Lol.Cyclotomic.Tensor.CPP.Backend+import Crypto.Lol.Cyclotomic.Tensor.CPP.Extension+import Crypto.Lol.Cyclotomic.Tensor.CPP.Instances ()+import Crypto.Lol.GaussRandom+import Crypto.Lol.Prelude as LP hiding+ (replicate,+ unzip, zip)+import Crypto.Lol.Types.FiniteField+import Crypto.Lol.Types.IZipVector+import Crypto.Lol.Types.Proto+import Crypto.Lol.Utils.ShowType++import Data.Foldable as F++import System.IO.Unsafe (unsafePerformIO)++-- | Newtype wrapper around a Vector.+newtype CT' (m :: Factored) r = CT' { unCT :: Vector r }+ deriving (Show, Eq, NFData)++-- the first argument, though phantom, affects representation+type role CT' representational nominal++-- GADT wrapper that distinguishes between Unbox and unrestricted+-- element types++-- | An implementation of 'Tensor' backed by C++ code.+data CT (m :: Factored) r where+ CT :: Storable r => CT' m r -> CT m r+ ZV :: IZipVector m r -> CT m r++deriving instance Show r => Show (CT m r)++instance Show (ArgType CT) where+ show _ = "CT"++instance Eq r => Eq (CT m r) where+ (ZV x) == (ZV y) = x == y+ (CT x) == (CT y) = x == y+ x@(CT _) == y = x == toCT y+ y == x@(CT _) = x == toCT y++instance (Protoable (IZipVector m r), Fact m, Storable r) => Protoable (CT m r) where+ type ProtoType (CT m r) = ProtoType (IZipVector m r)++ toProto x@(CT _) = toProto $ toZV x+ toProto (ZV x) = toProto x++ fromProto x = toCT <$> ZV <$> fromProto x++toCT :: (Storable r) => CT m r -> CT m r+toCT v@(CT _) = v+toCT (ZV v) = CT $ zvToCT' v++toZV :: (Fact m) => CT m r -> CT m r+toZV (CT (CT' v)) = ZV $ fromMaybe (error "toZV: internal error") $+ iZipVector $ convert v+toZV v@(ZV _) = v++zvToCT' :: forall m r . (Storable r) => IZipVector m r -> CT' m r+zvToCT' v = coerce (convert $ unIZipVector v :: Vector r)++wrap :: (Storable s, Storable r) => (CT' l s -> CT' m r) -> (CT l s -> CT m r)+{-# INLINABLE wrap #-}+wrap f (CT v) = CT $ f v+wrap f (ZV v) = CT $ f $ zvToCT' v++wrapM :: (Storable s, Storable r, Monad mon) => (CT' l s -> mon (CT' m r))+ -> (CT l s -> mon (CT m r))+{-# INLINABLE wrapM #-}+wrapM f (CT v) = CT <$> f v+wrapM f (ZV v) = CT <$> f (zvToCT' v)++-- convert an CT' *twace* signature to Tagged one+type family Tw (r :: *) :: * where+ Tw (CT' m' r -> CT' m r) = Tagged '(m,m') (Vector r -> Vector r)+ Tw (Maybe (CT' m' r -> CT' m r)) = TaggedT '(m,m') Maybe (Vector r -> Vector r)++type family Em r where+ Em (CT' m r -> CT' m' r) = Tagged '(m,m') (Vector r -> Vector r)+ Em (Maybe (CT' m r -> CT' m' r)) = TaggedT '(m,m') Maybe (Vector r -> Vector r)+++---------- NUMERIC PRELUDE INSTANCES ----------++-- CJP: Additive, Ring are not necessary when we use zipWithT+-- EAC: This has performance implications for the CT backend,+-- which used a (very fast) C function for (*) and (+)+instance (Additive r, Storable r, Fact m)+ => Additive.C (CT m r) where+ (CT (CT' a)) + (CT (CT' b)) = CT $ CT' $ SV.zipWith (+) a b+ a + b = toCT a + toCT b+ negate (CT (CT' a)) = CT $ CT' $ SV.map negate a -- EAC: This probably should be converted to C code+ negate a = negate $ toCT a++ zero = CT $ repl zero++{-+instance (Fact m, Ring r, Storable r, Dispatch r)+ => Ring.C (CT m r) where+ (CT a@(CT' _)) * (CT b@(CT' _)) = CT $ (untag $ cZipDispatch dmul) a b++ fromInteger = CT . repl . fromInteger+-}++instance (ZeroTestable r, Storable r)+ => ZeroTestable.C (CT m r) where+ --{-# INLINABLE isZero #-}+ isZero (CT (CT' a)) = SV.foldl' (\ b x -> b && isZero x) True a+ isZero (ZV v) = isZero v++instance (GFCtx fp d, Fact m, Additive (CT m fp))+ => Module.C (GF fp d) (CT m fp) where++ r *> v = case v of+ CT (CT' arr) -> CT $ CT' $ SV.fromList $ unCoeffs $ r *> Coeffs $ SV.toList arr+ ZV zv -> ZV $ fromJust $ iZipVector $ V.fromList $ unCoeffs $ r *> Coeffs $ V.toList $ unIZipVector zv++---------- Category-theoretic instances ----------++instance Fact m => Functor (CT m) where+ -- Functor instance is implied by Applicative laws+ fmap f x = pure f <*> x++instance Fact m => Applicative (CT m) where+ pure = ZV . pure++ (ZV f) <*> (ZV a) = ZV (f <*> a)+ f@(ZV _) <*> v@(CT _) = f <*> toZV v++instance Fact m => Foldable (CT m) where+ -- Foldable instance is implied by Traversable+ foldMap = foldMapDefault++instance Fact m => Traversable (CT m) where+ traverse f r@(CT _) = T.traverse f $ toZV r+ traverse f (ZV v) = ZV <$> T.traverse f v++instance Tensor CT where++ type TElt CT r = (Storable r, Dispatch r)++ entailIndexT = tag $ Sub Dict+ entailEqT = tag $ Sub Dict+ entailZTT = tag $ Sub Dict+ -- entailRingT = tag $ Sub Dict+ entailNFDataT = tag $ Sub Dict+ entailRandomT = tag $ Sub Dict+ entailShowT = tag $ Sub Dict+ entailModuleT = tag $ Sub Dict++ scalarPow = CT . scalarPow' -- Vector code++ l = wrap $ basicDispatch dl+ lInv = wrap $ basicDispatch dlinv++ mulGPow = wrap $ basicDispatch dmulgpow+ mulGDec = wrap $ basicDispatch dmulgdec++ divGPow = wrapM $ dispatchGInv dginvpow+ divGDec = wrapM $ dispatchGInv dginvdec++ crtFuncs = (,,,,) <$>+ return (CT . repl) <*>+ (wrap . untag (cZipDispatch dmul) <$> gCRT) <*>+ (wrap . untag (cZipDispatch dmul) <$> gInvCRT) <*>+ (wrap <$> untagT ctCRT) <*>+ (wrap <$> untagT ctCRTInv)++ twacePowDec = wrap $ runIdentity $ coerceTw twacePowDec'+ embedPow = wrap $ runIdentity $ coerceEm embedPow'+ embedDec = wrap $ runIdentity $ coerceEm embedDec'++ tGaussianDec v = CT <$> cDispatchGaussian v+ --tGaussianDec v = CT <$> coerceT' (gaussianDec v)++ -- we do not wrap this function because (currently) it can only be called on lifted types+ gSqNormDec (CT v) = untag gSqNormDec' v+ gSqNormDec (ZV v) = gSqNormDec (CT $ zvToCT' v)++ crtExtFuncs = (,) <$> (wrap <$> coerceTw twaceCRT')+ <*> (wrap <$> coerceEm embedCRT')++ coeffs = wrapM $ coerceCoeffs coeffs'++ powBasisPow = (CT <$>) <$> coerceBasis powBasisPow'++ crtSetDec = (CT <$>) <$> coerceBasis crtSetDec'++ fmapT f = wrap $ coerce (SV.map f)++ zipWithT f v1' v2' =+ let (CT (CT' v1)) = toCT v1'+ (CT (CT' v2)) = toCT v2'+ in CT $ CT' $ SV.zipWith f v1 v2++ unzipT v =+ let (CT (CT' x)) = toCT v+ in (CT . CT') *** (CT . CT') $ unzip x++ {-# INLINABLE entailIndexT #-}+ {-# INLINABLE entailEqT #-}+ {-# INLINABLE entailZTT #-}+ {-# INLINABLE entailNFDataT #-}+ {-# INLINABLE entailRandomT #-}+ {-# INLINABLE entailShowT #-}+ {-# INLINABLE scalarPow #-}+ {-# INLINABLE l #-}+ {-# INLINABLE lInv #-}+ {-# INLINABLE mulGPow #-}+ {-# INLINABLE mulGDec #-}+ {-# INLINABLE divGPow #-}+ {-# INLINABLE divGDec #-}+ {-# INLINABLE crtFuncs #-}+ {-# INLINABLE twacePowDec #-}+ {-# INLINABLE embedPow #-}+ {-# INLINABLE embedDec #-}+ {-# INLINABLE tGaussianDec #-}+ {-# INLINABLE gSqNormDec #-}+ {-# INLINE crtExtFuncs #-}+ {-# INLINABLE coeffs #-}+ {-# INLINABLE powBasisPow #-}+ {-# INLINABLE crtSetDec #-}+ {-# INLINABLE fmapT #-}+ {-# INLINE zipWithT #-}+ {-# INLINE unzipT #-}++coerceTw :: (Functor mon) => TaggedT '(m, m') mon (Vector r -> Vector r) -> mon (CT' m' r -> CT' m r)+coerceTw = (coerce <$>) . untagT++coerceEm :: (Functor mon) => TaggedT '(m, m') mon (Vector r -> Vector r) -> mon (CT' m r -> CT' m' r)+coerceEm = (coerce <$>) . untagT++-- | Useful coersion for defining @coeffs@ in the @Tensor@+-- interface. Using 'coerce' alone is insufficient for type inference.+coerceCoeffs :: Tagged '(m,m') (Vector r -> [Vector r]) -> CT' m' r -> [CT' m r]+coerceCoeffs = coerce++-- | Useful coersion for defining @powBasisPow@ and @crtSetDec@ in the @Tensor@+-- interface. Using 'coerce' alone is insufficient for type inference.+coerceBasis :: Tagged '(m,m') [Vector r] -> Tagged m [CT' m' r]+coerceBasis = coerce++dispatchGInv :: forall m r . (Storable r, Fact m)+ => (Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO Int16)+ -> CT' m r -> Maybe (CT' m r)+dispatchGInv f =+ let factors = proxy (marshalFactors <$> ppsFact) (Proxy::Proxy m)+ totm = proxy (fromIntegral <$> totientFact) (Proxy::Proxy m)+ numFacts = fromIntegral $ SV.length factors+ in \(CT' x) -> unsafePerformIO $ do+ yout <- SV.thaw x+ ret <- SM.unsafeWith yout (\pout ->+ SV.unsafeWith factors (\pfac ->+ f pout totm pfac numFacts))+ if ret /= 0+ then Just . CT' <$> unsafeFreeze yout+ else return Nothing++withBasicArgs :: forall m r . (Fact m, Storable r)+ => (Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO ())+ -> CT' m r -> IO (CT' m r)+withBasicArgs f =+ let factors = proxy (marshalFactors <$> ppsFact) (Proxy::Proxy m)+ totm = proxy (fromIntegral <$> totientFact) (Proxy::Proxy m)+ numFacts = fromIntegral $ SV.length factors+ in \(CT' x) -> do+ yout <- SV.thaw x+ SM.unsafeWith yout (\pout ->+ SV.unsafeWith factors (\pfac ->+ f pout totm pfac numFacts))+ CT' <$> unsafeFreeze yout++basicDispatch :: (Storable r, Fact m)+ => (Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO ())+ -> CT' m r -> CT' m r+basicDispatch f = unsafePerformIO . withBasicArgs f++gSqNormDec' :: (Storable r, Fact m, Dispatch r)+ => Tagged m (CT' m r -> r)+gSqNormDec' = return $ (!0) . unCT . unsafePerformIO . withBasicArgs dnorm++ctCRT :: (Storable r, CRTrans mon r, Dispatch r, Fact m)+ => TaggedT m mon (CT' m r -> CT' m r)+ctCRT = do+ ru' <- ru+ return $ \x -> unsafePerformIO $+ withPtrArray ru' (flip withBasicArgs x . dcrt)++-- CTensor CRT^(-1) functions take inverse rus+ctCRTInv :: (Storable r, CRTrans mon r, Dispatch r, Fact m)+ => TaggedT m mon (CT' m r -> CT' m r)+ctCRTInv = do+ mhatInv <- snd <$> crtInfo+ ruinv' <- ruInv+ return $ \x -> unsafePerformIO $+ withPtrArray ruinv' (\ruptr -> with mhatInv (flip withBasicArgs x . dcrtinv ruptr))++cZipDispatch :: (Storable r, Fact m)+ => (Ptr r -> Ptr r -> Int64 -> IO ())+ -> Tagged m (CT' m r -> CT' m r -> CT' m r)+cZipDispatch f = do -- in Tagged m+ totm <- fromIntegral <$> totientFact+ return $ coerce $ \a b -> unsafePerformIO $ do+ yout <- SV.thaw a+ SM.unsafeWith yout (\pout ->+ SV.unsafeWith b (\pin ->+ f pout pin totm))+ unsafeFreeze yout++cDispatchGaussian :: forall m r var rnd .+ (Storable r, Transcendental r, Dispatch r, Ord r,+ Fact m, ToRational var, Random r, MonadRandom rnd)+ => var -> rnd (CT' m r)+cDispatchGaussian var = flip proxyT (Proxy::Proxy m) $ do -- in TaggedT m rnd+ -- get rus for (Complex r)+ -- takes ru (not ruInv) to match RT+ ruinv' <- mapTaggedT (return . fromMaybe (error "complexGaussianRoots")) ru+ totm <- pureT totientFact+ mval <- pureT valueFact+ rad <- pureT radicalFact+ yin <- T.lift $ realGaussians (var * fromIntegral (mval `div` rad)) totm+ return $ unsafePerformIO $+ withPtrArray ruinv' (\ruptr -> withBasicArgs (dgaussdec ruptr) (CT' yin))++instance (Storable r, Random r, Fact m) => Random (CT' m r) where+ --{-# INLINABLE random #-}+ random = runRand $ replM (liftRand random)++ randomR = error "randomR nonsensical for CT'"++instance (Storable r, Random (CT' m r)) => Random (CT m r) where+ --{-# INLINABLE random #-}+ random = runRand $ CT <$> liftRand random++ randomR = error "randomR nonsensical for CT"++instance (NFData r) => NFData (CT m r) where+ rnf (CT v) = rnf v+ rnf (ZV v) = rnf v++repl :: forall m r . (Fact m, Storable r) => r -> CT' m r+repl = let n = proxy totientFact (Proxy::Proxy m)+ in coerce . SV.replicate n++replM :: forall m r mon . (Fact m, Storable r, Monad mon)+ => mon r -> mon (CT' m r)+replM = let n = proxy totientFact (Proxy::Proxy m)+ in fmap coerce . SV.replicateM n++scalarPow' :: forall m r . (Fact m, Additive r, Storable r) => r -> CT' m r+-- constant-term coefficient is first entry wrt powerful basis+scalarPow' =+ let n = proxy totientFact (Proxy::Proxy m)+ in \r -> CT' $ generate n (\i -> if i == 0 then r else zero)++ru, ruInv :: (CRTrans mon r, Fact m, Storable r)+ => TaggedT m mon [Vector r]+ru = do+ mval <- pureT valueFact+ wPow <- fst <$> crtInfo+ LP.map+ (\(p,e) -> do+ let pp = p^e+ pow = mval `div` pp+ generate pp (wPow . (*pow))) <$>+ pureT ppsFact++ruInv = do+ mval <- pureT valueFact+ wPow <- fst <$> crtInfo+ LP.map+ (\(p,e) -> do+ let pp = p^e+ pow = mval `div` pp+ generate pp (\i -> wPow $ -i*pow)) <$>+ pureT ppsFact++wrapVector :: forall mon m r . (Monad mon, Fact m, Ring r, Storable r)+ => TaggedT m mon (Kron r) -> mon (CT' m r)+wrapVector v = do+ vmat <- proxyT v (Proxy::Proxy m)+ let n = proxy totientFact (Proxy::Proxy m)+ return $ CT' $ generate n (flip (indexK vmat) 0)++gCRT, gInvCRT :: (Storable r, CRTrans mon r, Fact m)+ => mon (CT' m r)+gCRT = wrapVector gCRTK+gInvCRT = wrapVector gInvCRTK
+ Crypto/Lol/Cyclotomic/Tensor/CPP/Backend.hs view
@@ -0,0 +1,337 @@+{-|+Module : Crypto.Lol.Cyclotomic.Tensor.CPP.Backend+Description : Transforms Haskell types into C counterparts.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX++This module contains the functions to transform Haskell types into their+C counterpart, and to transform polymorphic Haskell functions into C funtion+calls in a type-safe way.+-}++{-# LANGUAGE ConstraintKinds #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE FlexibleInstances #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE PolyKinds #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeOperators #-}+{-# LANGUAGE UndecidableInstances #-}++{-# OPTIONS_GHC -fno-warn-unticked-promoted-constructors #-}++module Crypto.Lol.Cyclotomic.Tensor.CPP.Backend+( Dispatch+, dcrt, dcrtinv+, dgaussdec+, dl, dlinv+, dnorm+, dmulgpow, dmulgdec+, dginvpow, dginvdec+, dmul+, marshalFactors+, CPP+, withArray, withPtrArray+) where++import Crypto.Lol.Prelude as LP (Complex, PP, Proxy (..), Tagged,+ map, mapM_, proxy, tag)+import Crypto.Lol.Reflects+import Crypto.Lol.Types.Unsafe.RRq+import Crypto.Lol.Types.Unsafe.ZqBasic++import Data.Int+import Data.Vector.Storable as SV (Vector, fromList,+ unsafeToForeignPtr0)+import Data.Vector.Storable.Internal (getPtr)++import Foreign.ForeignPtr (touchForeignPtr)+import Foreign.Marshal.Array (withArray)+import Foreign.Marshal.Utils (with)+import Foreign.Ptr (Ptr, castPtr, plusPtr)+import Foreign.Storable (Storable (..))++import GHC.TypeLits -- for error message++-- | Convert a list of prime powers to a suitable C representation.+marshalFactors :: [PP] -> Vector CPP+marshalFactors = SV.fromList . LP.map (\(p,e) -> (fromIntegral p, fromIntegral e))++-- http://stackoverflow.com/questions/6517387/vector-vector-foo-ptr-ptr-foo-io-a-io-a+-- | Evaluates a C function that takes an "a** ptr" on a list of Vectors.+withPtrArray :: (Storable a) => [Vector a] -> (Ptr (Ptr a) -> IO b) -> IO b+withPtrArray v f = do+ let vs = LP.map SV.unsafeToForeignPtr0 v+ ptrV = LP.map (\(fp,_) -> getPtr fp) vs+ res <- withArray ptrV f+ LP.mapM_ (\(fp,_) -> touchForeignPtr fp) vs+ return res++-- Note: These types need to be the same, otherwise something goes wrong on the C end...+-- | C representation of a prime power.+type CPP = (Int16, Int16)++instance (Storable a, Storable b)+ => Storable (a,b) where+ sizeOf _ = sizeOf (undefined :: a) + sizeOf (undefined :: b)+ alignment _ = max (alignment (undefined :: a)) (alignment (undefined :: b))+ peek p = do+ a <- peek (castPtr p :: Ptr a)+ b <- peek (castPtr (plusPtr p (sizeOf a)) :: Ptr b)+ return (a,b)+ poke p (a,b) = do+ poke (castPtr p :: Ptr a) a+ poke (castPtr (plusPtr p (sizeOf a)) :: Ptr b) b++data ZqB64D -- for type safety purposes+data ComplexD+data DoubleD+data Int64D+data RRqD++type family CTypeOf x where+ CTypeOf (a,b) = EqCType a b (CTypeOf a) (CTypeOf b)+ CTypeOf (ZqBasic (q :: k) Int64) = ZqB64D+ CTypeOf Double = DoubleD+ CTypeOf Int64 = Int64D+ CTypeOf (Complex Double) = ComplexD+ CTypeOf (RRq (q :: k) Double) = RRqD++ -- EAC: See #12237 and #11990+ CTypeOf (ZqBasic (q :: k) i) = TypeError (Text "Unsupported C type: " :<>: ShowType (ZqBasic q i) :$$: Text "Use Int64 as the base ring")+ CTypeOf (Complex i) = TypeError (Text "Unsupported C type: " :<>: ShowType (Complex i) :$$: Text "Use Double as the base ring")+ CTypeOf (RRq (q :: k) i) = TypeError (Text "Unsupported C type: " :<>: ShowType (RRq q i) :$$: Text "Use Double as the base ring")+ CTypeOf a = TypeError (Text "Unsupported C type: " :<>: ShowType a)++type family EqCType a b c d where+ EqCType a b ZqB64D ZqB64D = ZqB64D+ EqCType a b RRqD RRqD = RRqD+ EqCType a b ComplexD ComplexD = ComplexD+ EqCType a b c c = TypeError (Text "Cannot call C code on a tuple of type " :<>: ShowType a)+ EqCType a b c d = TypeError (Text "You are trying to use CTensor on a tuple," :<>:+ Text " but the tuple contains two different C types: " :$$:+ ShowType a :<>: Text " and " :<>: ShowType b)++-- returns the modulus as a nested list of moduli+class (Tuple a) => ZqTuple a where+ type ModPairs a+ getModuli :: Tagged a (ModPairs a)++instance (Reflects q Int64) => ZqTuple (ZqBasic q Int64) where+ type ModPairs (ZqBasic q Int64) = Int64+ getModuli = tag $ proxy value (Proxy::Proxy q)++instance (Reflects q r, RealFrac r) => ZqTuple (RRq q r) where+ type ModPairs (RRq q r) = Int64+ getModuli = tag $ round (proxy value (Proxy::Proxy q) :: r)++instance (ZqTuple a, ZqTuple b) => ZqTuple (a, b) where+ type ModPairs (a,b) = (ModPairs a, ModPairs b)+ getModuli =+ let as = proxy getModuli (Proxy::Proxy a)+ bs = proxy getModuli (Proxy :: Proxy b)+ in tag (as,bs)++-- counts components in a nested tuple+class Tuple a where+ numComponents :: Tagged a Int16++instance {-# Overlappable #-} Tuple a where+ numComponents = tag 1++instance (Tuple a, Tuple b) => Tuple (a,b) where+ numComponents = tag $ proxy numComponents (Proxy::Proxy a) + proxy numComponents (Proxy::Proxy b)++-- | Single-argument synonym for @Dispatch'@.+type Dispatch r = (Dispatch' (CTypeOf r) r)++-- | Class to safely match Haskell types with the appropriate C function.+class (repr ~ CTypeOf r) => Dispatch' repr r where+ -- | Equivalent to 'Tensor's @crt@.+ dcrt :: Ptr (Ptr r) -> Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO ()+ -- | Equivalent to 'Tensor's @crtInv@.+ dcrtinv :: Ptr (Ptr r) -> Ptr r -> Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO ()+ -- | Equivalent to 'Tensor's @tGaussianDec@.+ dgaussdec :: Ptr (Ptr (Complex r)) -> Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO ()+ -- | Equivalent to 'Tensor's @l@.+ dl :: Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO ()+ -- | Equivalent to 'Tensor's @lInv@.+ dlinv :: Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO ()+ -- | Equivalent to 'Tensor's @gSqNormDec@.+ dnorm :: Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO ()+ -- | Equivalent to 'Tensor's @mulGPow@.+ dmulgpow :: Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO ()+ -- | Equivalent to 'Tensor's @mulGDec@.+ dmulgdec :: Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO ()+ -- | Equivalent to 'Tensor's @divGPow@.+ dginvpow :: Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO Int16+ -- | Equivalent to 'Tensor's @divGDec@.+ dginvdec :: Ptr r -> Int64 -> Ptr CPP -> Int16 -> IO Int16+ -- | Equivalent to @zipWith (*)@+ dmul :: Ptr r -> Ptr r -> Int64 -> IO ()++instance (ZqTuple r, Storable (ModPairs r), CTypeOf r ~ RRqD)+ => Dispatch' RRqD r where+ dcrt = error "cannot call CT CRT on type RRq"+ dcrtinv = error "cannot call CT CRTInv on type RRq"+ dl = error "cannot call CT L on type RRq (though you probably should be able to)"+ dlinv = error "cannot call CT LInv on type RRq (though you probably should be able to)"+ dnorm = error "cannto call CT normSq on type RRq"+ dmulgpow = error "cannot call CT mulGPow on type RRq"+ dmulgdec = error "cannot call CT mulGDec on type RRq"+ dginvpow = error "cannot call CT divGPow on type RRq"+ dginvdec = error "cannot call CT divGDec on type RRq"+ dmul = error "cannot call CT mul on type RRq"+ dgaussdec = error "cannot call CT gaussianDec on type RRq"++instance (ZqTuple r, Storable (ModPairs r), CTypeOf r ~ ZqB64D)+ => Dispatch' ZqB64D r where+ dcrt ruptr pout totm pfac numFacts =+ let qs = proxy getModuli (Proxy::Proxy r)+ numPairs = proxy numComponents (Proxy::Proxy r)+ in with qs $ \qsptr ->+ tensorCRTRq numPairs (castPtr pout) totm pfac numFacts (castPtr ruptr) (castPtr qsptr)+ dcrtinv ruptr minv pout totm pfac numFacts =+ let qs = proxy getModuli (Proxy::Proxy r)+ numPairs = proxy numComponents (Proxy::Proxy r)+ in with qs $ \qsptr ->+ tensorCRTInvRq numPairs (castPtr pout) totm pfac numFacts (castPtr ruptr) (castPtr minv) (castPtr qsptr)+ dl pout totm pfac numFacts =+ let qs = proxy getModuli (Proxy::Proxy r)+ numPairs = proxy numComponents (Proxy::Proxy r)+ in with qs $ \qsptr ->+ tensorLRq numPairs (castPtr pout) totm pfac numFacts (castPtr qsptr)+ dlinv pout totm pfac numFacts =+ let qs = proxy getModuli (Proxy::Proxy r)+ numPairs = proxy numComponents (Proxy::Proxy r)+ in with qs $ \qsptr ->+ tensorLInvRq numPairs (castPtr pout) totm pfac numFacts (castPtr qsptr)+ dnorm = error "cannot call CT normSq on type ZqBasic"+ dmulgpow pout totm pfac numFacts =+ let qs = proxy getModuli (Proxy::Proxy r)+ numPairs = proxy numComponents (Proxy::Proxy r)+ in with qs $ \qsptr ->+ tensorGPowRq numPairs (castPtr pout) totm pfac numFacts (castPtr qsptr)+ dmulgdec pout totm pfac numFacts =+ let qs = proxy getModuli (Proxy::Proxy r)+ numPairs = proxy numComponents (Proxy::Proxy r)+ in with qs $ \qsptr ->+ tensorGDecRq numPairs (castPtr pout) totm pfac numFacts (castPtr qsptr)+ dginvpow pout totm pfac numFacts =+ let qs = proxy getModuli (Proxy::Proxy r)+ numPairs = proxy numComponents (Proxy::Proxy r)+ in with qs $ \qsptr ->+ tensorGInvPowRq numPairs (castPtr pout) totm pfac numFacts (castPtr qsptr)+ dginvdec pout totm pfac numFacts =+ let qs = proxy getModuli (Proxy::Proxy r)+ numPairs = proxy numComponents (Proxy::Proxy r)+ in with qs $ \qsptr ->+ tensorGInvDecRq numPairs (castPtr pout) totm pfac numFacts (castPtr qsptr)+ dmul aout bout totm =+ let qs = proxy getModuli (Proxy::Proxy r)+ numPairs = proxy numComponents (Proxy::Proxy r)+ in with qs $ \qsptr ->+ mulRq numPairs (castPtr aout) (castPtr bout) totm (castPtr qsptr)+ dgaussdec = error "cannot call CT gaussianDec on type ZqBasic"++-- products of Complex correspond to CRTExt of a Zq product+instance (Tuple r, CTypeOf r ~ ComplexD) => Dispatch' ComplexD r where+ dcrt ruptr pout totm pfac numFacts =+ tensorCRTC (proxy numComponents (Proxy::Proxy r)) (castPtr pout) totm pfac numFacts (castPtr ruptr)+ dcrtinv ruptr minv pout totm pfac numFacts =+ tensorCRTInvC (proxy numComponents (Proxy::Proxy r)) (castPtr pout) totm pfac numFacts (castPtr ruptr) (castPtr minv)+ dl pout =+ tensorLC (proxy numComponents (Proxy::Proxy r)) (castPtr pout)+ dlinv pout =+ tensorLInvC (proxy numComponents (Proxy::Proxy r)) (castPtr pout)+ dnorm = error "cannot call CT normSq on type Complex Double"+ dmulgpow pout =+ tensorGPowC (proxy numComponents (Proxy::Proxy r)) (castPtr pout)+ dmulgdec pout =+ tensorGDecC (proxy numComponents (Proxy::Proxy r)) (castPtr pout)+ dginvpow pout =+ tensorGInvPowC (proxy numComponents (Proxy::Proxy r)) (castPtr pout)+ dginvdec pout =+ tensorGInvDecC (proxy numComponents (Proxy::Proxy r)) (castPtr pout)+ dmul aout bout =+ mulC (proxy numComponents (Proxy::Proxy r)) (castPtr aout) (castPtr bout)+ dgaussdec = error "cannot call CT gaussianDec on type Comple Double"++-- no support for products of Double+instance Dispatch' DoubleD Double where+ dcrt = error "cannot call CT Crt on type Double"+ dcrtinv = error "cannot call CT CrtInv on type Double"+ dl pout =+ tensorLDouble 1 (castPtr pout)+ dlinv pout =+ tensorLInvDouble 1 (castPtr pout)+ dnorm pout = tensorNormSqD 1 (castPtr pout)+ dmulgpow = error "cannot call CT mulGPow on type Double"+ dmulgdec = error "cannot call CT mulGDec on type Double"+ dginvpow = error "cannot call CT divGPow on type Double"+ dginvdec = error "cannot call CT divGDec on type Double"+ dmul = error "cannot call CT (*) on type Double"+ dgaussdec ruptr pout totm pfac numFacts =+ tensorGaussianDec 1 (castPtr pout) totm pfac numFacts (castPtr ruptr)++-- no support for products of Z+instance Dispatch' Int64D Int64 where+ dcrt = error "cannot call CT Crt on type Int64"+ dcrtinv = error "cannot call CT CrtInv on type Int64"+ dl pout =+ tensorLR 1 (castPtr pout)+ dlinv pout =+ tensorLInvR 1 (castPtr pout)+ dnorm pout =+ tensorNormSqR 1 (castPtr pout)+ dmulgpow pout =+ tensorGPowR 1 (castPtr pout)+ dmulgdec pout =+ tensorGDecR 1 (castPtr pout)+ dginvpow pout =+ tensorGInvPowR 1 (castPtr pout)+ dginvdec pout =+ tensorGInvDecR 1 (castPtr pout)+ dmul = error "cannot call CT (*) on type Int64"+ dgaussdec = error "cannot call CT gaussianDec on type Int64"++foreign import ccall unsafe "tensorLR" tensorLR :: Int16 -> Ptr Int64 -> Int64 -> Ptr CPP -> Int16 -> IO ()+foreign import ccall unsafe "tensorLInvR" tensorLInvR :: Int16 -> Ptr Int64 -> Int64 -> Ptr CPP -> Int16 -> IO ()+foreign import ccall unsafe "tensorLRq" tensorLRq :: Int16 -> Ptr (ZqBasic q Int64) -> Int64 -> Ptr CPP -> Int16 -> Ptr Int64 -> IO ()+foreign import ccall unsafe "tensorLInvRq" tensorLInvRq :: Int16 -> Ptr (ZqBasic q Int64) -> Int64 -> Ptr CPP -> Int16 -> Ptr Int64 -> IO ()+foreign import ccall unsafe "tensorLDouble" tensorLDouble :: Int16 -> Ptr Double -> Int64 -> Ptr CPP -> Int16 -> IO ()+foreign import ccall unsafe "tensorLInvDouble" tensorLInvDouble :: Int16 -> Ptr Double -> Int64 -> Ptr CPP -> Int16 -> IO ()+foreign import ccall unsafe "tensorLC" tensorLC :: Int16 -> Ptr (Complex Double) -> Int64 -> Ptr CPP -> Int16 -> IO ()+foreign import ccall unsafe "tensorLInvC" tensorLInvC :: Int16 -> Ptr (Complex Double) -> Int64 -> Ptr CPP -> Int16 -> IO ()++foreign import ccall unsafe "tensorNormSqR" tensorNormSqR :: Int16 -> Ptr Int64 -> Int64 -> Ptr CPP -> Int16 -> IO ()+foreign import ccall unsafe "tensorNormSqD" tensorNormSqD :: Int16 -> Ptr Double -> Int64 -> Ptr CPP -> Int16 -> IO ()++foreign import ccall unsafe "tensorGPowR" tensorGPowR :: Int16 -> Ptr Int64 -> Int64 -> Ptr CPP -> Int16 -> IO ()+foreign import ccall unsafe "tensorGPowRq" tensorGPowRq :: Int16 -> Ptr (ZqBasic q Int64) -> Int64 -> Ptr CPP -> Int16 -> Ptr Int64 -> IO ()+foreign import ccall unsafe "tensorGPowC" tensorGPowC :: Int16 -> Ptr (Complex Double) -> Int64 -> Ptr CPP -> Int16 -> IO ()+foreign import ccall unsafe "tensorGDecR" tensorGDecR :: Int16 -> Ptr Int64 -> Int64 -> Ptr CPP -> Int16 -> IO ()+foreign import ccall unsafe "tensorGDecRq" tensorGDecRq :: Int16 -> Ptr (ZqBasic q Int64) -> Int64 -> Ptr CPP -> Int16 -> Ptr Int64 -> IO ()+foreign import ccall unsafe "tensorGDecC" tensorGDecC :: Int16 -> Ptr (Complex Double) -> Int64 -> Ptr CPP -> Int16 -> IO ()+foreign import ccall unsafe "tensorGInvPowR" tensorGInvPowR :: Int16 -> Ptr Int64 -> Int64 -> Ptr CPP -> Int16 -> IO Int16+foreign import ccall unsafe "tensorGInvPowRq" tensorGInvPowRq :: Int16 -> Ptr (ZqBasic q Int64) -> Int64 -> Ptr CPP -> Int16 -> Ptr Int64 -> IO Int16+foreign import ccall unsafe "tensorGInvPowC" tensorGInvPowC :: Int16 -> Ptr (Complex Double) -> Int64 -> Ptr CPP -> Int16 -> IO Int16+foreign import ccall unsafe "tensorGInvDecR" tensorGInvDecR :: Int16 -> Ptr Int64 -> Int64 -> Ptr CPP -> Int16 -> IO Int16+foreign import ccall unsafe "tensorGInvDecRq" tensorGInvDecRq :: Int16 -> Ptr (ZqBasic q Int64) -> Int64 -> Ptr CPP -> Int16 -> Ptr Int64 -> IO Int16+foreign import ccall unsafe "tensorGInvDecC" tensorGInvDecC :: Int16 -> Ptr (Complex Double) -> Int64 -> Ptr CPP -> Int16 -> IO Int16++foreign import ccall unsafe "tensorCRTRq" tensorCRTRq :: Int16 -> Ptr (ZqBasic q Int64) -> Int64 -> Ptr CPP -> Int16 -> Ptr (Ptr (ZqBasic q Int64)) -> Ptr Int64 -> IO ()+foreign import ccall unsafe "tensorCRTC" tensorCRTC :: Int16 -> Ptr (Complex Double) -> Int64 -> Ptr CPP -> Int16 -> Ptr (Ptr (Complex Double)) -> IO ()+foreign import ccall unsafe "tensorCRTInvRq" tensorCRTInvRq :: Int16 -> Ptr (ZqBasic q Int64) -> Int64 -> Ptr CPP -> Int16 -> Ptr (Ptr (ZqBasic q Int64)) -> Ptr (ZqBasic q Int64) -> Ptr Int64 -> IO ()+foreign import ccall unsafe "tensorCRTInvC" tensorCRTInvC :: Int16 -> Ptr (Complex Double) -> Int64 -> Ptr CPP -> Int16 -> Ptr (Ptr (Complex Double)) -> Ptr (Complex Double) -> IO ()++foreign import ccall unsafe "tensorGaussianDec" tensorGaussianDec :: Int16 -> Ptr Double -> Int64 -> Ptr CPP -> Int16 -> Ptr (Ptr (Complex Double)) -> IO ()++foreign import ccall unsafe "mulRq" mulRq :: Int16 -> Ptr (ZqBasic q Int64) -> Ptr (ZqBasic q Int64) -> Int64 -> Ptr Int64 -> IO ()+foreign import ccall unsafe "mulC" mulC :: Int16 -> Ptr (Complex Double) -> Ptr (Complex Double) -> Int64 -> IO ()
+ Crypto/Lol/Cyclotomic/Tensor/CPP/Extension.hs view
@@ -0,0 +1,164 @@+{-|+Module : Crypto.Lol.Cyclotomic.Tensor.CPP.Extension+Description : Embedding/twacing in various bases for CPP.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX++CPP Tensor-specific functions for embedding/twacing in various bases.+-}++{-# LANGUAGE ConstraintKinds #-}+{-# LANGUAGE DataKinds #-}+{-# LANGUAGE FlexibleContexts #-}+{-# LANGUAGE MultiParamTypeClasses #-}+{-# LANGUAGE PolyKinds #-}+{-# LANGUAGE RebindableSyntax #-}+{-# LANGUAGE ScopedTypeVariables #-}+{-# LANGUAGE TypeFamilies #-}+{-# LANGUAGE TypeOperators #-}++module Crypto.Lol.Cyclotomic.Tensor.CPP.Extension+( embedPow', embedDec', embedCRT'+, twacePowDec', twaceCRT'+, coeffs', powBasisPow'+, crtSetDec'+, backpermute'+) where++import Crypto.Lol.CRTrans+import Crypto.Lol.Cyclotomic.Tensor as T+import Crypto.Lol.Prelude as LP hiding (lift, null)+import Crypto.Lol.Types.FiniteField+import Crypto.Lol.Types.ZmStar+++import Control.Applicative hiding (empty)+import Control.Monad.Trans (lift)++import Data.Maybe+import Data.Reflection (reify)+import qualified Data.Vector as V+import Data.Vector.Storable as SV+import qualified Data.Vector.Unboxed as U+++-- | /O(n)/ Yield the vector obtained by replacing each element @i@ of the+-- index vector by @xs'!'i@. This is equivalent to @'map' (xs'!') is@ but is+-- often much more efficient.+--+-- > backpermute <a,b,c,d> <0,3,2,3,1,0> = <a,d,c,d,b,a>+backpermute' :: (Storable a) =>+ U.Vector Int -- ^ @is@ index vector (of length @n@)+ -> Vector a -- ^ @xs@ value vector+ -> Vector a+{-# INLINABLE backpermute' #-}+backpermute' is v = generate (U.length is) (\i -> v ! (is U.! i))++embedPow', embedDec' :: (Additive r, Storable r, m `Divides` m')+ => Tagged '(m, m') (Vector r -> Vector r)+{-# INLINABLE embedPow' #-}+{-# INLINABLE embedDec' #-}+-- | Embeds an vector in the powerful basis of the the mth cyclotomic ring+-- to an vector in the powerful basis of the m'th cyclotomic ring when @m | m'@+embedPow' = (\indices arr -> generate (U.length indices) $ \idx ->+ let (j0,j1) = indices U.! idx+ in if j0 == 0+ then arr ! j1+ else zero) <$> baseIndicesPow+-- | Embeds an vector in the decoding basis of the the mth cyclotomic ring+-- to an vector in the decoding basis of the m'th cyclotomic ring when @m | m'@+embedDec' = (\indices arr -> generate (U.length indices)+ (\idx -> maybe LP.zero+ (\(sh,b) -> if b then negate (arr ! sh) else arr ! sh)+ (indices U.! idx))) <$> baseIndicesDec++-- | Embeds an vector in the CRT basis of the the mth cyclotomic ring+-- to an vector in the CRT basis of the m'th cyclotomic ring when @m | m'@+embedCRT' :: forall mon m m' r . (CRTrans mon r, Storable r, m `Divides` m')+ => TaggedT '(m, m') mon (Vector r -> Vector r)+embedCRT' =+ (lift (proxyT crtInfo (Proxy::Proxy m') :: mon (CRTInfo r))) >>+ (pureT $ backpermute' <$> baseIndicesCRT)++-- | maps a vector in the powerful/decoding basis, representing an+-- O_m' element, to a vector of arrays representing O_m elements in+-- the same type of basis+coeffs' :: (Storable r, m `Divides` m')+ => Tagged '(m, m') (Vector r -> [Vector r])+coeffs' = flip (\x -> V.toList . V.map (`backpermute'` x))+ <$> extIndicesCoeffs++-- | The "tweaked trace" function in either the powerful or decoding+-- basis of the m'th cyclotomic ring to the mth cyclotomic ring when+-- @m | m'@.+twacePowDec' :: forall m m' r . (Storable r, m `Divides` m')+ => Tagged '(m, m') (Vector r -> Vector r)+{-# INLINABLE twacePowDec' #-}+twacePowDec' = backpermute' <$> extIndicesPowDec++kronToVec :: forall mon m r . (Monad mon, Fact m, Ring r, Storable r)+ => TaggedT m mon (Kron r) -> TaggedT m mon (Vector r)+kronToVec v = do+ vmat <- v+ let n = proxy totientFact (Proxy::Proxy m)+ return $ generate n (flip (indexK vmat) 0)++twaceCRT' :: forall mon m m' r .+ (Storable r, CRTrans mon r, m `Divides` m')+ => TaggedT '(m, m') mon (Vector r -> Vector r)+{-# INLINE twaceCRT' #-}+twaceCRT' = tagT $ do+ g' <- proxyT (kronToVec gCRTK) (Proxy::Proxy m')+ gInv <- proxyT (kronToVec gInvCRTK) (Proxy::Proxy m)+ embed <- proxyT embedCRT' (Proxy::Proxy '(m,m'))+ indices <- pure $ proxy extIndicesCRT (Proxy::Proxy '(m,m'))+ (_, m'hatinv) <- proxyT crtInfo (Proxy::Proxy m')+ let phi = proxy totientFact (Proxy::Proxy m)+ phi' = proxy totientFact (Proxy::Proxy m')+ mhat = fromIntegral $ proxy valueHatFact (Proxy::Proxy m)+ hatRatioInv = m'hatinv * mhat+ reltot = phi' `div` phi+ -- tweak = mhat * g' / (m'hat * g)+ tweak = SV.map (* hatRatioInv) $ SV.zipWith (*) (embed gInv) g'+ return $ \ arr -> -- take true trace after mul-by-tweak+ let v = backpermute' indices (SV.zipWith (*) tweak arr)+ in generate phi $ \i -> foldl1' (+) $ SV.unsafeSlice (i*reltot) reltot v++-- | The powerful extension basis, wrt the powerful basis.+-- Outputs a list of vectors in O_m' that are an O_m basis for O_m'+powBasisPow' :: forall m m' r . (m `Divides` m', Ring r, SV.Storable r)+ => Tagged '(m, m') [SV.Vector r]+powBasisPow' = do+ (_, phi, phi', _) <- indexInfo+ idxs <- baseIndicesPow+ return $ LP.map (\k -> generate phi' $ \j ->+ let (j0,j1) = idxs U.! j+ in if j0==k && j1==0 then one else zero)+ [0..phi' `div` phi - 1]++-- | A list of vectors representing the mod-p CRT set of the+-- extension O_m'/O_m+crtSetDec' :: forall m m' fp .+ (m `Divides` m', PrimeField fp, Coprime (PToF (CharOf fp)) m', SV.Storable fp)+ => Tagged '(m, m') [SV.Vector fp]+crtSetDec' =+ let m'p = Proxy :: Proxy m'+ p = proxy valuePrime (Proxy::Proxy (CharOf fp))+ phi = proxy totientFact m'p+ d = proxy (order p) m'p+ h :: Int = proxy valueHatFact m'p+ hinv = recip $ fromIntegral h+ in reify d $ \(_::Proxy d) -> do+ let twCRTs' :: Kron (GF fp d)+ = fromMaybe (error "internal error: crtSetDec': twCRTs") $ proxyT twCRTs m'p+ zmsToIdx = proxy T.zmsToIndexFact m'p+ elt j i = indexK twCRTs' j (zmsToIdx i)+ trace' = trace :: GF fp d -> fp -- to avoid recomputing powTraces+ cosets <- partitionCosets p+ return $ LP.map (\is -> generate phi+ (\j -> hinv * trace'+ (LP.sum $ LP.map (elt j) is))) cosets
+ Crypto/Lol/Cyclotomic/Tensor/CPP/Instances.hs view
@@ -0,0 +1,30 @@+{-|+Module : Crypto.Lol.Cyclotomic.Tensor.CPP.Instances+Description : CPP Tensor-specific instances.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX++CPP Tensor-specific instances.+-}++{-# LANGUAGE GeneralizedNewtypeDeriving #-}+{-# LANGUAGE PolyKinds #-}+{-# LANGUAGE StandaloneDeriving #-}++module Crypto.Lol.Cyclotomic.Tensor.CPP.Instances where++-- EAC: Do not import Crypto.Lol.Types, because it exports an IrreduciblePoly+-- instance which screw with GHC. Probably #10338.+import Crypto.Lol.Types.Unsafe.Complex+import Crypto.Lol.Types.Unsafe.RRq+import Crypto.Lol.Types.Unsafe.ZqBasic++import Foreign.Storable++deriving instance (Storable a) => Storable (Complex a)+deriving instance (Storable r) => Storable (RRq q r)+deriving instance (Storable i) => Storable (ZqBasic q i)
+ Crypto/Lol/Cyclotomic/Tensor/CPP/common.cpp view
@@ -0,0 +1,27 @@+/* +Module : common.cpp +Description : Shared functions and data. +Copyright : (c) Eric Crockett, 2011-2017 + Chris Peikert, 2011-2017 +License : GPL-2 +Maintainer : ecrockett0@email.com +Stability : experimental +Portability : POSIX +*/ + +#include "types.h" + +hInt_t Zq::q; // should be in zq.cpp; here due to GHC #12152 + +hDim_t ipow(hDim_t base, hShort_t exp) +{ + hDim_t result = 1; + while (exp) { + if (exp & 1) { + result *= base; + } + exp >>= 1; + base *= base; + } + return result; +}
+ Crypto/Lol/Cyclotomic/Tensor/CPP/common.h view
@@ -0,0 +1,20 @@+/*+Module : common.h+Description : Shared functions and data.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX+*/++#ifndef COMMON_H_+#define COMMON_H_++#include "types.h"++// calculates base ** exp+hDim_t ipow(hDim_t base, hShort_t exp);++#endif /* COMMON_H_ */
+ Crypto/Lol/Cyclotomic/Tensor/CPP/crt.cpp view
@@ -0,0 +1,598 @@+/*+Module : crt.cpp+Description : Chinese remainder transform.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX+*/++#include "types.h"+#include "tensor.h"+#include "common.h"++// If this macro is modified, make sure to update all functions below with+// cascading if/else statments so that temp space is allocated when necessary+// (i.e., for all primes >= DFTP_GENERIC_SIZE)+#define DFTP_GENERIC_SIZE 11++hDim_t bitrev (PrimeExponent pe, hDim_t j) {+ hShort_t e;+ hDim_t p = pe.prime;+ hDim_t tempj = j;+ hDim_t acc = 0;++ for(e = pe.exponent-1; e >= 0; e--) {+ div_t qr = div(tempj,p);+ acc += qr.rem * ipow(p,e);+ tempj = qr.quot;+ }+ return acc;+}++template <typename ring> void crtTwiddle (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts,+ PrimeExponent pe, ring* ru)+{+ hDim_t p = pe.prime;+ hShort_t e = pe.exponent;++ pe.exponent -= 1; // used for an argument to bitrev++ if(p == 2) {+ hDim_t mprime = 1<<(e-1);+ hDim_t blockDim = rts*mprime; // size of block in block diagonal tensor matrix++ for(hDim_t i0 = 1; i0 < mprime; i0++) { // loops over i/(p-1) for i = 0..(m'-1), we can skip i0 = 0+ hDim_t temp2 = i0*rts;+ ring twid = ru[bitrev(pe, i0)*tupSize];++ for(hDim_t blockIdx = 0; blockIdx < lts; blockIdx++) {+ hDim_t temp3 = blockIdx*blockDim + temp2;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ hDim_t idx = (temp3 + modOffset)*tupSize;+ y[idx] *= twid;+ }+ }+ }+ }+ else { // This loop is faster, probably due to the division in the loop above.+ // cilk also slows it down+ hDim_t mprime = ipow(p,e-1);+ hDim_t blockDim = rts*(p-1)*mprime; // size of block in block diagonal tensor matrix++ for(hDim_t i0 = 1; i0 < mprime; i0++) { // loops over i/(p-1) for i = 0..(m'-1), we can skip i0 = 0+ hDim_t temp1 = i0*(p-1);+ for(hDim_t i1 = 0; i1 < (p-1); i1++) { // loops over i%(p-1) for i = 0..(m'-1)+ hDim_t temp2 = (temp1+i1)*rts;+ ring twid = ru[bitrev(pe, i0)*(i1+1)*tupSize];++ for(hDim_t blockIdx = 0; blockIdx < lts; blockIdx++) {+ hDim_t temp3 = blockIdx*blockDim + temp2;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ hDim_t idx = (temp3 + modOffset)*tupSize;+ y[idx] *= twid;+ }+ }+ }+ }+ }+}++// dim is power of p+template <typename ring> void dftTwiddle (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts,+ PrimeExponent pe, hDim_t dim, hDim_t rustride, ring* ru)+{+ hDim_t idx;+ hDim_t p = pe.prime;++ pe.exponent -= 1; // used for an argument to bitrev++ if(p == 2) {+ hDim_t mprime = dim>>1; // divides evenly+ hDim_t temp1 = rts*dim; // for use in computing [modified] tensorOffset+ for(hDim_t i0 = 1; i0 < mprime; i0++) { // loops over i/p for i = 0..(dim-1), but we skip i0=0+ hDim_t temp3 = rts*(i0*p+1);+ ring twid = ru[bitrev(pe,i0)*rustride*tupSize];++ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1 + temp3;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ idx = (temp2 + modOffset)*tupSize;+ y[idx] *= twid;+ }+ }+ }+ }+ else {+ hDim_t mprime = dim/p; // divides evenly+ hDim_t temp1 = rts*dim; // for use in computing [modified] tensorOffset+ for(hDim_t i0 = 1; i0 < mprime; i0++) { // loops over i/p for i = 0..(dim-1), but we skip i0=0+ for(hDim_t i1 = 1; i1 < p; i1++) { // loops over i%p for i = 0..(dim-1), but we skip i1=0+ hDim_t temp3 = rts*(i0*p+i1);+ ring twid = ru[bitrev(pe,i0)*i1*rustride*tupSize];++ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1 + temp3;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ idx = (temp2 + modOffset)*tupSize;+ y[idx] *= twid;+ }+ }+ }+ }+ }+}++//implied length of ru is rustride*p+//implied length of tempSpace is p, if p is not a special case+// temp is allowed to be NULL if p < DFTP_GENERIC_SIZE+template <typename ring> void dftp (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts,+ hDim_t p, hDim_t rustride, ring* ru, ring* tempSpace)+{+ hDim_t tensorOffset;++ if(p == 2) {+ hDim_t temp1 = rts<<1;++ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring u = y[tensorOffset*tupSize];+ ring t = y[(tensorOffset+rts)*tupSize];+ y[tensorOffset*tupSize] = u + t;+ y[(tensorOffset+rts)*tupSize] = u - t;+ }+ }+ }+ else if(p == 3) {+ ring ru1 = ru[rustride*tupSize];+ ring ru2 = ru[(rustride<<1)*tupSize];+ hDim_t temp1 = rts*3;++ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring y1, y2, y3;+ y1 = y[tensorOffset*tupSize];+ y2 = y[(tensorOffset+rts)*tupSize];+ y3 = y[(tensorOffset+(rts<<1))*tupSize];+ //q is <32 bits, so we can do 3 additions without overflow+ y[tensorOffset*tupSize] += (y2 + y3);+ y[(tensorOffset+rts)*tupSize] = y1 + (ru1*y2) + (ru2*y3);+ y[(tensorOffset+(rts<<1))*tupSize] = y1 + (ru2*y2) + (ru1*y3);+ }+ }+ }+ else if(p == 5) {+ hDim_t temp1 = rts*5;+ ring ru1 = ru[rustride*tupSize];+ ring ru2 = ru[(rustride<<1)*tupSize];+ ring ru3 = ru[(rustride*3)*tupSize];+ ring ru4 = ru[(rustride<<2)*tupSize];++ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring y1, y2, y3, y4, y5;+ y1 = y[tensorOffset*tupSize];+ y2 = y[(tensorOffset+rts)*tupSize];+ y3 = y[(tensorOffset+(rts<<1))*tupSize];+ y4 = y[(tensorOffset+3*rts)*tupSize];+ y5 = y[(tensorOffset+(rts<<2))*tupSize];+ y[tensorOffset*tupSize] += y2 + y3 + y4 + y5;+ y[(tensorOffset+rts)*tupSize] = y1 + (ru1*y2) + (ru2*y3) + (ru3*y4) + (ru4*y5);+ y[(tensorOffset+(rts<<1))*tupSize] = y1 + (ru2*y2) + (ru4*y3) + (ru1*y4) + (ru3*y5);+ y[(tensorOffset+rts*3)*tupSize] = y1 + (ru3*y2) + (ru1*y3) + (ru4*y4) + (ru2*y5);+ y[(tensorOffset+(rts<<2))*tupSize] = y1 + (ru4*y2) + (ru3*y3) + (ru2*y4) + (ru1*y5);+ }+ }+ }+ else if(p == 7) {+ hDim_t temp1 = rts*7;+ ring ru1 = ru[rustride*tupSize];+ ring ru2 = ru[(rustride<<1)*tupSize];+ ring ru3 = ru[(rustride*3)*tupSize];+ ring ru4 = ru[(rustride<<2)*tupSize];+ ring ru5 = ru[(rustride*5)*tupSize];+ ring ru6 = ru[(rustride*6)*tupSize];++ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring y1, y2, y3, y4, y5, y6, y7;+ y1 = y[tensorOffset*tupSize];+ y2 = y[(tensorOffset+rts)*tupSize];+ y3 = y[(tensorOffset+(rts<<1))*tupSize];+ y4 = y[(tensorOffset+3*rts)*tupSize];+ y5 = y[(tensorOffset+(rts<<2))*tupSize];+ y6 = y[(tensorOffset+rts*5)*tupSize];+ y7 = y[(tensorOffset+rts*6)*tupSize];+ y[tensorOffset*tupSize] += y2 + y3 + y4 + y5 + y6 + y7;+ y[(tensorOffset+rts)*tupSize] = y1 + (ru1*y2) + (ru2*y3) + (ru3*y4) + (ru4*y5) + (ru5*y6) + (ru6*y7);+ y[(tensorOffset+(rts<<1))*tupSize] = y1 + (ru2*y2) + (ru4*y3) + (ru6*y4) + (ru1*y5) + (ru3*y6) + (ru5*y7);+ y[(tensorOffset+rts*3)*tupSize] = y1 + (ru3*y2) + (ru6*y3) + (ru2*y4) + (ru5*y5) + (ru1*y6) + (ru4*y7);+ y[(tensorOffset+(rts<<2))*tupSize] = y1 + (ru4*y2) + (ru1*y3) + (ru5*y4) + (ru2*y5) + (ru6*y6) + (ru3*y7);+ y[(tensorOffset+rts*5)*tupSize] = y1 + (ru5*y2) + (ru3*y3) + (ru1*y4) + (ru6*y5) + (ru4*y6) + (ru2*y7);+ y[(tensorOffset+rts*6)*tupSize] = y1 + (ru6*y2) + (ru5*y3) + (ru4*y4) + (ru3*y5) + (ru2*y6) + (ru1*y7);+ }+ }+ }+ else {+ hDim_t temp1 = rts*p;+ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ for(hDim_t row = 0; row < p; row++) {+ tempSpace[row] = 0;+ //p is small (<< 30 bits), so we can do p additions of mod-q values without overflow+ for(hDim_t col = 0; col < p; col++) {+ tempSpace[row] += (y[(tensorOffset+col*rts)*tupSize]*ru[((col*row) % p)*rustride*tupSize]);+ }+ }++ for(hDim_t row = 0; row < p; row++) {+ y[(tensorOffset+rts*row)*tupSize] = tempSpace[row];+ }+ }+ }+ }+}++template <typename ring> void crtp (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts,+ hDim_t p, hDim_t rustride, ring* ru)+{+ hDim_t tensorOffset;+ if(p == 2) {+ return;+ }+ else if(p == 3) {+ hDim_t temp1 = rts*2;+ ring ru1 = ru[rustride*tupSize];+ ring ru2 = ru[(rustride<<1)*tupSize];++ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring y1, y2;+ y1 = y[tensorOffset*tupSize];+ y2 = y[(tensorOffset+rts)*tupSize];+ y[tensorOffset*tupSize] += (ru1*y2);+ y[(tensorOffset+rts)*tupSize] = y1 + (ru2*y2);+ }+ }+ }+ else if(p == 5) {+ hDim_t temp1 = rts*4;+ ring ru1 = ru[rustride*tupSize];+ ring ru2 = ru[(rustride<<1)*tupSize];+ ring ru3 = ru[(rustride*3)*tupSize];+ ring ru4 = ru[(rustride<<2)*tupSize];++ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring y1, y2, y3, y4;+ y1 = y[tensorOffset*tupSize];+ y2 = y[(tensorOffset+rts)*tupSize];+ y3 = y[(tensorOffset+(rts<<1))*tupSize];+ y4 = y[(tensorOffset+3*rts)*tupSize];++ y[tensorOffset*tupSize] += ((ru1*y2) + (ru2*y3) + (ru3*y4));+ y[(tensorOffset+rts)*tupSize] = y1 + (ru2*y2) + (ru4*y3) + (ru1*y4);+ y[(tensorOffset+(rts<<1))*tupSize] = y1 + (ru3*y2) + (ru1*y3) + (ru4*y4);+ y[(tensorOffset+rts*3)*tupSize] = y1 + (ru4*y2) + (ru3*y3) + (ru2*y4);+ }+ }+ }+ else if(p == 7) {+ hDim_t temp1 = rts*6;+ ring ru1 = ru[rustride*tupSize];+ ring ru2 = ru[(rustride<<1)*tupSize];+ ring ru3 = ru[(rustride*3)*tupSize];+ ring ru4 = ru[(rustride<<2)*tupSize];+ ring ru5 = ru[(rustride*5)*tupSize];+ ring ru6 = ru[(rustride*6)*tupSize];+ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring y1, y2, y3, y4, y5, y6;+ y1 = y[tensorOffset*tupSize];+ y2 = y[(tensorOffset+rts)*tupSize];+ y3 = y[(tensorOffset+(rts<<1))*tupSize];+ y4 = y[(tensorOffset+3*rts)*tupSize];+ y5 = y[(tensorOffset+(rts<<2))*tupSize];+ y6 = y[(tensorOffset+rts*5)*tupSize];+ y[tensorOffset*tupSize] += ((ru1*y2) + (ru2*y3) + (ru3*y4) + (ru4*y5) + (ru5*y6));+ y[(tensorOffset+rts)*tupSize] = y1 + (ru2*y2) + (ru4*y3) + (ru6*y4) + (ru1*y5) + (ru3*y6);+ y[(tensorOffset+(rts<<1))*tupSize] = y1 + (ru3*y2) + (ru6*y3) + (ru2*y4) + (ru5*y5) + (ru1*y6);+ y[(tensorOffset+rts*3)*tupSize] = y1 + (ru4*y2) + (ru1*y3) + (ru5*y4) + (ru2*y5) + (ru6*y6);+ y[(tensorOffset+(rts<<2))*tupSize] = y1 + (ru5*y2) + (ru3*y3) + (ru1*y4) + (ru6*y5) + (ru4*y6);+ y[(tensorOffset+rts*5)*tupSize] = y1 + (ru6*y2) + (ru5*y3) + (ru4*y4) + (ru3*y5) + (ru2*y6);+ }+ }+ }+ else {+ ring* tempSpace = (ring*)malloc((p-1)*sizeof(ring));+ hDim_t temp1 = rts*(p-1);+ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;++ for(hDim_t row = 1; row < p; row++) {+ tempSpace[row-1] = 0;+ for(hDim_t col = 0; col < p-1; col++) {+ tempSpace[row-1] += (y[(tensorOffset+col*rts)*tupSize]*ru[((col*row) % p)*rustride*tupSize]);+ }+ }++ for(hDim_t row = 0; row < p-1; row++) {+ y[(tensorOffset+rts*row)*tupSize] = tempSpace[row];+ }+ }+ }+ free(tempSpace);+ }+}++//takes inverse rus+template <typename ring> void crtpinv (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts,+ hDim_t p, hDim_t rustride, ring* ruinv)+{+ hDim_t tensorOffset;+ if(p == 2) {+ return;+ }+ else if(p == 3) {+ hDim_t temp1 = rts*2;+ ring ru1 = ruinv[rustride*tupSize];+ ring ru2 = ruinv[(rustride<<1)*tupSize];++ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring y1, y2, shift;+ y1 = y[tensorOffset*tupSize];+ y2 = y[(tensorOffset+rts)*tupSize];++ shift = (ru2*y1) + (ru1*y2);++ y[tensorOffset*tupSize] += y2 - shift;+ y[(tensorOffset+rts)*tupSize] = (ru1*y1) + (ru2*y2) - shift;+ }+ }+ }+ else if(p == 5) {+ hDim_t temp1 = rts*4;+ ring ru1 = ruinv[rustride*tupSize];+ ring ru2 = ruinv[(rustride<<1)*tupSize];+ ring ru3 = ruinv[(rustride*3)*tupSize];+ ring ru4 = ruinv[(rustride<<2)*tupSize];++ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring y1, y2, y3, y4, shift;+ y1 = y[tensorOffset*tupSize];+ y2 = y[(tensorOffset+rts)*tupSize];+ y3 = y[(tensorOffset+(rts<<1))*tupSize];+ y4 = y[(tensorOffset+3*rts)*tupSize];++ shift = (ru4*y1) + (ru3*y2) + (ru2*y3) + (ru1*y4);++ y[tensorOffset*tupSize] += y2 + y3 + y4 - shift;+ y[(tensorOffset+rts)*tupSize] = (ru1*y1) + (ru2*y2) + (ru3*y3) + (ru4*y4) - shift;+ y[(tensorOffset+(rts<<1))*tupSize] = (ru2*y1) + (ru4*y2) + (ru1*y3) + (ru3*y4) - shift;+ y[(tensorOffset+rts*3)*tupSize] = (ru3*y1) + (ru1*y2) + (ru4*y3) + (ru2*y4) - shift;+ }+ }+ }+ else if(p == 7) {+ hDim_t temp1 = rts*6;+ ring ru1 = ruinv[rustride*tupSize];+ ring ru2 = ruinv[(rustride<<1)*tupSize];+ ring ru3 = ruinv[(rustride*3)*tupSize];+ ring ru4 = ruinv[(rustride<<2)*tupSize];+ ring ru5 = ruinv[(rustride*5)*tupSize];+ ring ru6 = ruinv[(rustride*6)*tupSize];+ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring y1, y2, y3, y4, y5, y6, shift;+ y1 = y[tensorOffset*tupSize];+ y2 = y[(tensorOffset+rts)*tupSize];+ y3 = y[(tensorOffset+(rts<<1))*tupSize];+ y4 = y[(tensorOffset+3*rts)*tupSize];+ y5 = y[(tensorOffset+(rts<<2))*tupSize];+ y6 = y[(tensorOffset+rts*5)*tupSize];++ shift = (ru6*y1) + (ru5*y2) + (ru4*y3) + (ru3*y4) + (ru2*y5) + (ru1*y6);++ y[tensorOffset*tupSize] += y2 + y3 + y4 + y5 + y6 - shift;+ y[(tensorOffset+rts)*tupSize] = (ru1*y1) + (ru2*y2) + (ru3*y3) + (ru4*y4) + (ru5*y5) + (ru6*y6) - shift;+ y[(tensorOffset+(rts<<1))*tupSize] = (ru2*y1) + (ru4*y2) + (ru6*y3) + (ru1*y4) + (ru3*y5) + (ru5*y6) - shift;+ y[(tensorOffset+rts*3)*tupSize] = (ru3*y1) + (ru6*y2) + (ru2*y3) + (ru5*y4) + (ru1*y5) + (ru4*y6) - shift;+ y[(tensorOffset+(rts<<2))*tupSize] = (ru4*y1) + (ru1*y2) + (ru5*y3) + (ru2*y4) + (ru6*y5) + (ru3*y6) - shift;+ y[(tensorOffset+rts*5)*tupSize] = (ru5*y1) + (ru3*y2) + (ru1*y3) + (ru6*y4) + (ru4*y5) + (ru2*y6) - shift;+ }+ }+ }+ else {+ ring* tempSpace = (ring*)malloc((p-1)*sizeof(ring));+ hDim_t temp1 = rts*(p-1);+ for(hDim_t blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(hDim_t modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ ring shift;+ shift = 0;+ for(hDim_t row = 0; row < p-1; row++) {+ shift += (y[(tensorOffset+row*rts)*tupSize]*ruinv[(p-row-1)*rustride*tupSize]);+ tempSpace[row] = 0;+ for(hDim_t col = 0; col < p-1; col++) {+ tempSpace[row] += (y[(tensorOffset+col*rts)*tupSize]*ruinv[((row*(col+1)) % p)*rustride*tupSize]);+ }+ }++ for(hDim_t row = 0; row < p-1; row++) {+ y[(tensorOffset+rts*row)*tupSize] = tempSpace[row] - shift;+ }+ }+ }+ free(tempSpace);+ }+}++template <typename ring> void ppDFT (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts,+ PrimeExponent pe, hDim_t rustride, ring* ru, ring* temp)+{+ hDim_t p = pe.prime;+ hShort_t e = pe.exponent;++ if(e == 0) {+ return;+ }++ hDim_t primeRuStride = rustride*ipow(p,e-1);++ hShort_t i;++ hDim_t ltsScale = ipow(p,e-1);+ hDim_t rtsScale = 1;+ hDim_t twidRuStride = rustride;+ for(i = 0; i < e; i++) {+ hDim_t rtsDim = rts*rtsScale;+ dftp (y, tupSize, lts*ltsScale, rtsDim, p, primeRuStride, ru, temp);+ dftTwiddle (y, tupSize, lts, rtsDim, pe, ltsScale*p, twidRuStride, ru);++ ltsScale /= p;+ rtsScale *= p;+ twidRuStride *= p;+ pe.exponent -= 1;+ }+}++template <typename ring> void ppDFTInv (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts,+ PrimeExponent pe, hDim_t rustride, ring* ru, ring* temp)+{+ hDim_t p = pe.prime;+ hShort_t e = pe.exponent;++ if(e == 0) {+ return;+ }+ hDim_t primeRuStride = rustride*ipow(p,e-1);++ hShort_t i;++ hDim_t ltsScale = 1;+ hDim_t rtsScale = ipow(p,e-1);+ hDim_t twidRuStride = primeRuStride;+ pe.exponent = 1;+ for(i = 0; i < e; i++) {+ hDim_t rtsDim = rts*rtsScale;+ hDim_t ltsScaleP = ltsScale*p;+ dftTwiddle (y, tupSize, lts, rtsDim, pe, ltsScaleP, twidRuStride, ru);+ dftp (y, tupSize, lts*ltsScale, rtsDim, p, primeRuStride, ru, temp);++ ltsScale = ltsScaleP;+ rtsScale /= p;+ twidRuStride /= p;+ pe.exponent += 1;+ }+}++template <typename ring> void ppcrt (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts,+ PrimeExponent pe, ring* ru)+{+ hDim_t p = pe.prime;+ hDim_t e = pe.exponent;+ hDim_t mprime = ipow(p,e-1);+ ring* temp = 0;+ if(p >= DFTP_GENERIC_SIZE) {+ temp = (ring*)malloc(p*sizeof(ring));+ }++ crtp (y, tupSize, lts*mprime, rts, p, mprime, ru);+ crtTwiddle (y, tupSize, lts, rts, pe, ru);+ pe.exponent -= 1;+ ppDFT (y, tupSize, lts, rts*(p-1), pe, p, ru, temp);+ pe.exponent += 1;++ if(p >= DFTP_GENERIC_SIZE) {+ free(temp);+ }+}++template <typename ring> void ppcrtinv (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts,+ PrimeExponent pe, ring* ru)+{+ hDim_t p = pe.prime;+ hDim_t e = pe.exponent;+ hDim_t mprime = ipow(p,e-1);+ ring* temp = 0;+ if(p >= DFTP_GENERIC_SIZE) {+ temp = (ring*)malloc(p*sizeof(ring));+ }++ pe.exponent -= 1;+ ppDFTInv (y, tupSize, lts, rts*(p-1), pe, p, ru, temp);+ pe.exponent += 1;+ crtTwiddle (y, tupSize, lts, rts, pe, ru);+ crtpinv (y, tupSize, lts*mprime, rts, p, mprime, ru);++ if(p >= DFTP_GENERIC_SIZE) {+ free(temp);+ }+}++extern "C" void tensorCRTRq (hShort_t tupSize, Zq* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE, Zq** ru, hInt_t* qs)+{+ tensorFuserCRT (y, tupSize, ppcrt, totm, peArr, sizeOfPE, ru, qs);+ canonicalizeZq(y,tupSize,totm,qs);+}++//takes inverse rus+extern "C" void tensorCRTInvRq (hShort_t tupSize, Zq* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE,+ Zq** ruinv, Zq* mhatInv, hInt_t* qs)+{+ tensorFuserCRT (y, tupSize, ppcrtinv, totm, peArr, sizeOfPE, ruinv, qs);+ for (hShort_t i = 0; i < tupSize; i++) {+ Zq::q = qs[i];+ for (hDim_t j = 0; j < totm; j++) {+ //careful here! I'm not setting the global q, so I can't rely on Zq multiplication+ y[j*tupSize+i] = y[j*tupSize+i]*mhatInv[i];+ }+ }+ canonicalizeZq(y,tupSize,totm,qs);+}++extern "C" void tensorCRTC (hShort_t tupSize, Complex* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE, Complex** ru)+{+ tensorFuserCRT (y, tupSize, ppcrt, totm, peArr, sizeOfPE, ru, (hInt_t*)0);+}++//takes inverse rus+extern "C" void tensorCRTInvC (hShort_t tupSize, Complex* y, hDim_t totm, PrimeExponent* peArr,+ hShort_t sizeOfPE, Complex** ruinv, Complex* mhatInv)+{+ tensorFuserCRT (y, tupSize, ppcrtinv, totm, peArr, sizeOfPE, ruinv, (hInt_t*)0);+ for (hShort_t i = 0; i < tupSize; i++) {+ for (hDim_t j = 0; j < totm; j++) {+ y[j*tupSize+i] *= mhatInv[i];+ }+ }+}
+ Crypto/Lol/Cyclotomic/Tensor/CPP/g.cpp view
@@ -0,0 +1,273 @@+/* +Module : g.cpp +Description : Multiplication and division by 'g' in different bases. +Copyright : (c) Eric Crockett, 2011-2017 + Chris Peikert, 2011-2017 +License : GPL-2 +Maintainer : ecrockett0@email.com +Stability : experimental +Portability : POSIX +*/ + +#include "types.h" +#include "tensor.h" +#include "common.h" + +template <typename ring> void gPow (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts, hDim_t p) +{ + if (p == 2) {return;} + hDim_t tmp1 = rts*(p-1); + hDim_t tmp2 = tmp1 - rts; + hDim_t blockOffset, modOffset; + hDim_t i; + for (blockOffset = 0; blockOffset < lts; ++blockOffset) { + hDim_t tmp3 = blockOffset * tmp1; + for (modOffset = 0; modOffset < rts; ++modOffset) { + hDim_t tensorOffset = tmp3 + modOffset; + ring last = y[(tensorOffset + tmp2)*tupSize]; + for (i = p-2; i != 0; --i) { + hDim_t idx = tensorOffset + i * rts; + y[idx*tupSize] += (last - y[(idx-rts)*tupSize]); + } + y[tensorOffset*tupSize] += last; + } + } +} + +template <typename ring> void gDec (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts, hDim_t p) +{ + if (p == 2) {return;} + hDim_t tmp1 = rts*(p-1); + hDim_t blockOffset; + hDim_t modOffset; + hDim_t i; + + for (blockOffset = 0; blockOffset < lts; ++blockOffset) { + hDim_t tmp2 = blockOffset * tmp1; + for (modOffset = 0; modOffset < rts; ++modOffset) { + hDim_t tensorOffset = tmp2 + modOffset; + ring acc = y[tensorOffset*tupSize]; + for (i = p-2; i != 0; --i) { + hDim_t idx = tensorOffset + i * rts; + acc += y[idx*tupSize]; + y[idx*tupSize] -= y[(idx-rts)*tupSize]; + } + y[tensorOffset*tupSize] += acc; + } + } +} + +template <typename ring> void gInvPow (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts, hDim_t p) +{ + if (p == 2) {return;} + hDim_t tmp1 = rts * (p-1); + hDim_t blockOffset, modOffset; + hDim_t i; + + for (blockOffset = 0; blockOffset < lts; ++blockOffset) { + hDim_t tmp2 = blockOffset * tmp1; + for (modOffset = 0; modOffset < rts; ++modOffset) { + hDim_t tensorOffset = tmp2 + modOffset; + ring lelts; + lelts = 0; + for (i = 0; i < p-1; ++i) { + lelts += y[(tensorOffset + i*rts)*tupSize]; + } + ring relts; + relts = 0; + for (i = p-2; i >= 0; --i) { + hDim_t idx = tensorOffset + i*rts; + ring z = y[idx*tupSize]; + ring lmul, rmul; + lmul = p-1-i; + rmul = i+1; + y[idx*tupSize] = lmul * lelts - rmul * relts; + lelts -= z; + relts += z; + } + } + } +} + +template <typename ring> void gInvDec (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts, hDim_t p) +{ + if (p == 2) {return;} + hDim_t blockOffset; + hDim_t modOffset; + hDim_t i; + hDim_t tmp1 = rts*(p-1); + + for (blockOffset = 0; blockOffset < lts; ++blockOffset) { + hDim_t tmp2 = blockOffset*tmp1; + for (modOffset = 0; modOffset < rts; ++modOffset) { + hDim_t tensorOffset = tmp2 + modOffset; + ring lastOut; + lastOut = 0; + for (i=1; i < p; ++i) { + ring ri; + ri = i; + lastOut += (ri * y[(tensorOffset + (i-1)*rts)*tupSize]); + } + ring rp; + rp = p; + ring acc = lastOut; + for (i = p-2; i > 0; --i) { + hDim_t idx = tensorOffset + i*rts; + ring tmp = acc; + acc -= y[idx*tupSize]*rp; + y[idx*tupSize] = tmp; + } + y[tensorOffset*tupSize] = acc; + } + } +} + +extern "C" void tensorGPowR (hShort_t tupSize, hInt_t* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, gPow, totm, peArr, sizeOfPE, (hInt_t*)0); +} + +extern "C" void tensorGPowRq (hShort_t tupSize, Zq* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE, hInt_t* qs) +{ + tensorFuserPrime (y, tupSize, gPow, totm, peArr, sizeOfPE, qs); + canonicalizeZq(y,tupSize,totm,qs); +} + +extern "C" void tensorGPowC (hShort_t tupSize, Complex* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, gPow, totm, peArr, sizeOfPE, (hInt_t*)0); +} + +extern "C" void tensorGDecR (hShort_t tupSize, hInt_t* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, gDec, totm, peArr, sizeOfPE, (hInt_t*)0); +} + +extern "C" void tensorGDecRq (hShort_t tupSize, Zq* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE, hInt_t* qs) +{ + tensorFuserPrime (y, tupSize, gDec, totm, peArr, sizeOfPE, qs); + canonicalizeZq(y,tupSize,totm,qs); +} + +extern "C" void tensorGDecC (hShort_t tupSize, Complex* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, gDec, totm, peArr, sizeOfPE, (hInt_t*)0); +} + +hInt_t oddRad(PrimeExponent* peArr, hShort_t sizeOfPE) { + hInt_t oddrad; + oddrad = 1; + for(int i = 0; i < sizeOfPE; i++) { + hShort_t p = peArr[i].prime; + if (p != 2) { + oddrad *= peArr[i].prime; + } + } + return oddrad; +} + +extern "C" hShort_t tensorGInvPowR (hShort_t tupSize, hInt_t* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, gInvPow, totm, peArr, sizeOfPE, (hInt_t*)0); + + hInt_t oddrad = oddRad(peArr, sizeOfPE); + + for(int i = 0; i < tupSize*totm; i++) { + if (y[i] % oddrad) { + y[i] /= oddrad; + } + else { + return 0; + } + } + return 1; +} + +extern "C" hShort_t tensorGInvPowRq (hShort_t tupSize, Zq* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE, hInt_t* qs) +{ + tensorFuserPrime (y, tupSize, gInvPow, totm, peArr, sizeOfPE, qs); + + hInt_t oddrad = oddRad(peArr, sizeOfPE); + + for(int i = 0; i < tupSize; i++) { + Zq::q = qs[i]; // global update + hInt_t ori = reciprocal(Zq::q, oddrad); + Zq oddradInv; + oddradInv = ori; + if (ori == 0) { + return 0; // error condition + } + for(hDim_t j = 0; j < totm; j++) { + y[j*tupSize+i] *= oddradInv; + } + } + + canonicalizeZq(y,tupSize,totm,qs); + return 1; +} + +extern "C" hShort_t tensorGInvPowC (hShort_t tupSize, Complex* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, gInvPow, totm, peArr, sizeOfPE, (hInt_t*)0); + + hInt_t oddrad = oddRad(peArr, sizeOfPE); + Complex oddradInv; + oddradInv = 1 / oddrad; + for(int i = 0; i < tupSize*totm; i++) { + y[i] *= oddradInv; + } + return 1; +} + +extern "C" hShort_t tensorGInvDecR (hShort_t tupSize, hInt_t* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, gInvDec, totm, peArr, sizeOfPE, (hInt_t*)0); + + hInt_t oddrad = oddRad(peArr, sizeOfPE); + + for(int i = 0; i < tupSize*totm; i++) { + if (y[i] % oddrad) { + y[i] /= oddrad; + } + else { + return 0; + } + } + return 1; +} + +extern "C" hShort_t tensorGInvDecRq (hShort_t tupSize, Zq* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE, hInt_t* qs) +{ + tensorFuserPrime (y, tupSize, gInvDec, totm, peArr, sizeOfPE, qs); + + hInt_t oddrad = oddRad(peArr, sizeOfPE); + + for(int i = 0; i < tupSize; i++) { + Zq::q = qs[i]; // global update + hInt_t ori = reciprocal(Zq::q, oddrad); + Zq oddradInv; + oddradInv = ori; + if (ori == 0) { + return 0; // error condition + } + for(hDim_t j = 0; j < totm; j++) { + y[j*tupSize+i] *= oddradInv; + } + } + + canonicalizeZq(y,tupSize,totm,qs); + return 1; +} + +extern "C" hShort_t tensorGInvDecC (hShort_t tupSize, Complex* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, gInvDec, totm, peArr, sizeOfPE, (hInt_t*)0); + + hInt_t oddrad = oddRad(peArr, sizeOfPE); + Complex oddradInv; + oddradInv = 1 / oddrad; + for(int i = 0; i < tupSize*totm; i++) { + y[i] *= oddradInv; + } + return 1; +}
+ Crypto/Lol/Cyclotomic/Tensor/CPP/l.cpp view
@@ -0,0 +1,98 @@+/* +Module : l.cpp +Description : Powerful <-> Decoding basis conversion. +Copyright : (c) Eric Crockett, 2011-2017 + Chris Peikert, 2011-2017 +License : GPL-2 +Maintainer : ecrockett0@email.com +Stability : experimental +Portability : POSIX +*/ + +#include "types.h" +#include "tensor.h" + +template <typename ring> void lp (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts, hDim_t p) +{ + hDim_t blockOffset; + hDim_t modOffset; + int i; + + if(p == 2) {return;} + + hDim_t tmp1 = rts*(p-1); + for (blockOffset = 0; blockOffset < lts; ++blockOffset) { + hDim_t tmp2 = blockOffset*tmp1; + for (modOffset = 0; modOffset < rts; ++modOffset) { + hDim_t idx = tmp2 + modOffset + rts; + for (i = 1; i < p-1; ++i) { + y[idx*tupSize] += y[(idx-rts)*tupSize]; + idx += rts; + } + } + } +} + +template <typename ring> void lpInv (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts, hDim_t p) +{ + hDim_t blockOffset; + hDim_t modOffset; + int i; + + if(p == 2) {return;} + + hDim_t tmp1 = rts*(p-1); + for (blockOffset = 0; blockOffset < lts; ++blockOffset) { + hDim_t tmp2 = blockOffset*tmp1; + for (modOffset = 0; modOffset < rts; ++ modOffset) { + hDim_t tensorOffset = tmp2 + modOffset; + hDim_t idx = tensorOffset + (p-2) * rts; + for (i = p-2; i != 0; --i) { + y[idx*tupSize] -= y[(idx-rts)*tupSize] ; + idx -= rts; + } + } + } +} + +extern "C" void tensorLRq (hShort_t tupSize, Zq* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE, hInt_t* qs) +{ + tensorFuserPrime (y, tupSize, lp, totm, peArr, sizeOfPE, qs); + canonicalizeZq(y,tupSize,totm,qs); +} + +extern "C" void tensorLR (hShort_t tupSize, hInt_t* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, lp, totm, peArr, sizeOfPE, (hInt_t*)0); +} + +extern "C" void tensorLDouble (hShort_t tupSize, double* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, lp, totm, peArr, sizeOfPE, (hInt_t*)0); +} + +extern "C" void tensorLC (hShort_t tupSize, Complex* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, lp, totm, peArr, sizeOfPE, (hInt_t*)0); +} + +extern "C" void tensorLInvRq (hShort_t tupSize, Zq* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE, hInt_t* qs) +{ + tensorFuserPrime (y, tupSize, lpInv, totm, peArr, sizeOfPE, qs); + canonicalizeZq(y,tupSize,totm,qs); +} + +extern "C" void tensorLInvR (hShort_t tupSize, hInt_t* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, lpInv, totm, peArr, sizeOfPE, (hInt_t*)0); +} + +extern "C" void tensorLInvDouble (hShort_t tupSize, double* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, lpInv, totm, peArr, sizeOfPE, (hInt_t*)0); +} + +extern "C" void tensorLInvC (hShort_t tupSize, Complex* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE) +{ + tensorFuserPrime (y, tupSize, lpInv, totm, peArr, sizeOfPE, (hInt_t*)0); +}
+ Crypto/Lol/Cyclotomic/Tensor/CPP/mul.cpp view
@@ -0,0 +1,35 @@+/*+Module : mul.cpp+Description : zipWith (*).+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX+*/++#include "types.h"++template <typename ring> void zipWithStar (ring* a, ring* b, hShort_t tupSize, hDim_t totm, hInt_t* qs)+{+ for(int tupIdx = 0; tupIdx < tupSize; tupIdx++) {+ if(qs) {+ Zq::q = qs[tupIdx];+ }+ for(int i = 0; i < totm; i++) {+ a[i*tupSize+tupIdx] *= b[i*tupSize+tupIdx];+ }+ }+}++//a = zipWith (*) a b+extern "C" void mulRq (hShort_t tupSize, Zq* a, Zq* b, hDim_t totm, hInt_t* qs)+{+ zipWithStar(a, b, tupSize, totm, qs);+}++extern "C" void mulC (hShort_t tupSize, Complex* a, Complex* b, hDim_t totm)+{+ zipWithStar(a, b, tupSize, totm, (hInt_t*)0);+}
+ Crypto/Lol/Cyclotomic/Tensor/CPP/norm.cpp view
@@ -0,0 +1,80 @@+/*+Module : norm.cpp+Description : Compute g*norm(x)^2.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX+*/++#include "types.h"+#include "tensor.h"++template <typename ring> void pNormSq (ring* y, hShort_t tupSize, hDim_t lts, hDim_t rts, hDim_t p)+{+ hDim_t blockOffset;+ hDim_t modOffset;+ hDim_t i;++ if(p==2) {return;}++ hDim_t tmp1 = rts*(p-1);+ for (blockOffset = 0; blockOffset < lts; ++blockOffset) {+ hDim_t tmp2 = blockOffset*tmp1;+ for (modOffset = 0; modOffset < rts; ++modOffset) {+ hDim_t tensorOffset = tmp2 + modOffset;+ ring sum = 0;+ for (i = 0; i < p-1; ++i) {+ sum += y[(tensorOffset + i*rts)*tupSize];+ }+ for (i = 0; i < p-1; ++i) {+ y[(tensorOffset + i*rts)*tupSize] += sum;+ }+ }+ }+}++extern "C" void tensorNormSqR (hShort_t tupSize, hInt_t* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE)+{+ hInt_t* tempSpace = (hInt_t*)malloc(totm*tupSize*sizeof(hInt_t));+ for(hDim_t i = 0; i < totm*tupSize; i++) {+ tempSpace[i]=y[i];+ }++ tensorFuserPrime(y, tupSize, pNormSq, totm, peArr, sizeOfPE, (hInt_t*)0);++ //do dot product and return in index 0+ for(int tupIdx = 0; tupIdx < tupSize; tupIdx++) {+ hInt_t dotprod = 0;+ for(hDim_t i = 0; i < totm; i++) {+ dotprod += (tempSpace[i*tupSize+tupIdx]*y[i*tupSize+tupIdx]);+ }++ y[tupIdx] = dotprod;+ }++ free(tempSpace);+}++extern "C" void tensorNormSqD (hShort_t tupSize, double* y, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE)+{+ double* tempSpace = (double*)malloc(totm*tupSize*sizeof(double));+ for(hDim_t i = 0; i < totm*tupSize; i++) {+ tempSpace[i]=y[i];+ }+ tensorFuserPrime(y, tupSize, pNormSq, totm, peArr, sizeOfPE, (hInt_t*)0);++ //do dot product and return in index 0+ for(int tupIdx = 0; tupIdx < tupSize; tupIdx++) {+ double dotprod = 0;+ for(hDim_t i = 0; i < totm; i++) {+ dotprod += (tempSpace[i*tupSize+tupIdx]*y[i*tupSize+tupIdx]);+ }++ y[tupIdx] = dotprod;+ }++ free(tempSpace);+}
+ Crypto/Lol/Cyclotomic/Tensor/CPP/random.cpp view
@@ -0,0 +1,64 @@+/*+Module : random.cpp+Description : Gaussian sampling on a lattice.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX+*/++#include "types.h"+#include "tensor.h"+#include "common.h"+#include <math.h>++// I had been negating the ru-idx, but this was causing a *negative* mod, resulting in a hard-to-find bug+// current behavior (taking rus, rather than ruInv) matches RT+void primeD (double *y, hShort_t tupSize, hDim_t lts, hDim_t rts, hDim_t p, hDim_t rustride, Complex* ru)+{+ if(p == 2) {+ return;+ }+ hDim_t blockOffset, modOffset, tensorOffset;+ double *tempSpace = (double*)malloc((p-1)*sizeof(double));+ hDim_t temp1 = rts*(p-1);+ for(blockOffset = 0; blockOffset < lts; blockOffset++) {+ hDim_t temp2 = blockOffset*temp1;+ for(modOffset = 0; modOffset < rts; modOffset++) {+ tensorOffset = temp2 + modOffset;+ hDim_t row, col;++ for(row = 0; row < p-1; row++) {+ double acc = 0;+ for(col = 1; col <= (p>>1); col++) {+ acc += 2 * ru[((row*col) % p)*rustride*tupSize].real * y[(tensorOffset+rts*(col-1))*tupSize];+ }+ for(col = (p>>1)+1; col <= p-1; col++) {+ acc += 2 * ru[((row*col) % p)*rustride*tupSize].imag * y[(tensorOffset+rts*(col-1))*tupSize];+ }+ tempSpace[row] = acc/sqrt(2);+ }++ for(row = 0; row < p-1; row++) {+ y[(tensorOffset+rts*row)*tupSize] = tempSpace[row];+ }+ }+ }+ free(tempSpace);+}++void ppD (double *y, hShort_t tupSize, hDim_t lts, hDim_t rts, PrimeExponent pe, Complex *ru)+{+ hDim_t p = pe.prime;+ hDim_t e = pe.exponent;+ hDim_t mprime = ipow(p,e-1);+ primeD (y, tupSize, lts*mprime, rts, p, mprime, ru);+}++//the contents of y will be destroyed, but should be initialized in Haskell-land to independent Guassians over the reals+extern "C" void tensorGaussianDec (hShort_t tupSize, double *y, hDim_t totm, PrimeExponent *peArr, hShort_t sizeOfPE, Complex** ru)+{+ tensorFuserCRT (y, tupSize, ppD, totm, peArr, sizeOfPE, ru, (hInt_t*)0);+}
+ Crypto/Lol/Cyclotomic/Tensor/CPP/tensor.h view
@@ -0,0 +1,67 @@+/*+Module : tensor.h+Description : Templates for the tensor DSL.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX+*/++#ifndef TENSOR_CPP_+#define TENSOR_CPP_++#include "types.h"+#include "common.h"+#ifdef __cplusplus+template <typename ring>+using primeFunc = void (*) (ring*, hShort_t, hDim_t, hDim_t, hDim_t);++template <typename ringy, typename ringru>+using primeCRTFunc = void (*) (ringy*, hShort_t, hDim_t, hDim_t, PrimeExponent, ringru*);++//for square transforms+template <typename ring> void tensorFuserPrime (ring* y, hShort_t tupSize, primeFunc<ring> f, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE, hInt_t* qs)+{+ hDim_t lts = totm;+ hDim_t rts = 1;+ hShort_t i;++ for (i = 0; i < sizeOfPE; ++i) {+ PrimeExponent pe = peArr[i];+ hDim_t ipow_pe = ipow(pe.prime, (pe.exponent-1));+ hDim_t dim = (pe.prime-1) * ipow_pe; // the totient of pe+ lts /= dim;+ for(int tupIdx = 0; tupIdx < tupSize; tupIdx++) {+ if(qs) {+ Zq::q = qs[tupIdx]; // global update+ }+ (*f) (y+tupIdx, tupSize, lts*ipow_pe, rts, pe.prime);+ }+ rts *= dim;+ }+}++template <typename ringy, typename ringru> void tensorFuserCRT (ringy* y, hShort_t tupSize, primeCRTFunc<ringy,ringru> f, hDim_t totm, PrimeExponent* peArr, hShort_t sizeOfPE, ringru** ru, hInt_t* qs)+{+ hDim_t lts = totm;+ hDim_t rts = 1;+ hShort_t i;++ for (i = 0; i < sizeOfPE; ++i) {+ PrimeExponent pe = peArr[i];+ hDim_t ipow_pe = ipow(pe.prime, (pe.exponent-1));+ hDim_t dim = (pe.prime-1) * ipow_pe; // the totient of pe+ lts /= dim;+ for(int tupIdx = 0; tupIdx < tupSize; tupIdx++) {+ if(qs) {+ Zq::q = qs[tupIdx]; // global update+ }+ (*f) (y+tupIdx, tupSize, lts, rts, pe, ru[i]+tupIdx);+ }+ rts *= dim;+ }+}+#endif /* __cplusplus */+#endif /* TENSOR_CPP_ */
+ Crypto/Lol/Cyclotomic/Tensor/CPP/types.h view
@@ -0,0 +1,170 @@+/* +Module : types.h +Description : C++ types and function declarations. +Copyright : (c) Eric Crockett, 2011-2017 + Chris Peikert, 2011-2017 +License : GPL-2 +Maintainer : ecrockett0@email.com +Stability : experimental +Portability : POSIX +*/ + +#ifndef TENSORTYPES_H_ +#define TENSORTYPES_H_ + +#include <inttypes.h> +#include <stdio.h> +#include <stdlib.h> + +typedef int64_t hInt_t ; +typedef int32_t hDim_t ; +typedef int16_t hShort_t ; +typedef int8_t hByte_t ; + +typedef struct +{ + hShort_t prime; + hShort_t exponent; +} PrimeExponent; + + +hInt_t reciprocal (hInt_t a, hInt_t b); + +#define ASSERT(EXP) { \ + if (!(EXP)) { \ + fprintf (stderr, "Assertion in file '%s' line %d : " #EXP " is false\n", __FILE__, __LINE__); \ + exit(-1); \ + } \ +} + +//http://stackoverflow.com/questions/37572628 +#ifdef __cplusplus +//http://stackoverflow.com/a/4421719 +class Zq +{ +public: + hInt_t x; + + static hInt_t q; // declared here, defined in generalfuncs.cpp + + Zq& operator=(const hInt_t& c) + { + this->x = c % q; + return *this; + } + Zq& operator+=(const Zq& b) + { + this->x += b.x; + this->x %= q; + return *this; + } + Zq& operator-=(const Zq& b) + { + this->x -= b.x; + this->x %= q; + return *this; + } + Zq& operator*=(const Zq& b) + { + this->x *= b.x; + this->x %= q; + return *this; + } + Zq& operator/=(const Zq& b) + { + Zq binv; + binv = reciprocal(q,b.x); + ASSERT (binv.x); // binv == 0 indicates that x is not invertible mod q + *this *= binv; + return *this; + } +}; +inline Zq operator+(Zq a, const Zq& b) +{ + a += b; + return a; +} +inline Zq operator-(Zq a, const Zq& b) +{ + a -= b; + return a; +} +inline Zq operator*(Zq a, const Zq& b) +{ + a *= b; + return a; +} +inline Zq operator/(Zq a, const Zq& b) +{ + a /= b; + return a; +} + +void canonicalizeZq (Zq* y, hShort_t tupSize, hDim_t totm, hInt_t* qs); + +class Complex +{ +public: + double real; + double imag; + + Complex& operator=(const hInt_t& c) + { + this->real = c; + this->imag = 0; + return *this; + } + Complex& operator+=(const Complex& b) + { + this->real = this->real+b.real; + this->imag = this->imag+b.imag; + return *this; + } + Complex& operator-=(const Complex& b) + { + this->real = this->real-b.real; + this->imag = this->imag-b.imag; + return *this; + } + Complex& operator*=(const Complex& b) + { + double a = this->real; + this->real = (a*b.real)-(this->imag*b.imag); + this->imag = (a*b.imag)+(this->imag*b.real); + return *this; + } + Complex& operator/=(const Complex& b) + { + Complex bconj; + bconj.real = b.real; + bconj.imag = -b.imag; + *this *= bconj; + double den = (b.real*b.real+b.imag*b.imag); + this->real /= den; + this->imag /= den; + return *this; + } +}; +inline Complex operator+(Complex a, const Complex& b) +{ + a += b; + return a; +} +inline Complex operator-(Complex a, const Complex& b) +{ + a -= b; + return a; +} +inline Complex operator*(Complex a, const Complex& b) +{ + a *= b; + return a; +} +inline Complex operator/(Complex a, const Complex& b) +{ + a /= b; + return a; +} + +#endif /* __cplusplus */ +#endif /* TENSORTYPES_H_ */
+ Crypto/Lol/Cyclotomic/Tensor/CPP/zq.cpp view
@@ -0,0 +1,50 @@+/*+Module : zq.cpp+Description : Implementation of Z_q-specific functions.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX+*/++#include "types.h"+#include "common.h"++// a is the field size. we are looking for reciprocal of b+hInt_t reciprocal (hInt_t a, hInt_t b)+{+ hInt_t fieldSize = a;++ hInt_t y = 1;+ hInt_t lasty = 0;+ while (b != 0) {+ hInt_t quotient = a / b;+ hInt_t tmp = a % b;+ a = b;+ b = tmp;+ tmp = y;+ y = lasty - quotient*y;+ lasty = tmp;+ }+ // if a!=1, then b is not invertible mod a+ if(a!=1) {+ return 0;+ }++ // this actually returns EITHER the reciprocal OR reciprocal + fieldSize+ hInt_t res = lasty + fieldSize;+ return res;+}++void canonicalizeZq (Zq* y, hShort_t tupSize, hDim_t totm, hInt_t* qs) {+ for(int tupIdx = 0; tupIdx<tupSize; tupIdx++) {+ hInt_t q = qs[tupIdx];+ for(hDim_t j = 0; j < totm; j++) {+ if(y[j*tupSize+tupIdx].x<0) {+ y[j*tupSize+tupIdx].x+=q;+ }+ }+ }+}
+ LICENSE view
@@ -0,0 +1,339 @@+ GNU GENERAL PUBLIC LICENSE+ Version 2, June 1991++ Copyright (C) 1989, 1991 Free Software Foundation, Inc.,+ 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA+ Everyone is permitted to copy and distribute verbatim copies+ of this license document, but changing it is not allowed.++ Preamble++ The licenses for most software are designed to take away your+freedom to share and change it. By contrast, the GNU General Public+License is intended to guarantee your freedom to share and change free+software--to make sure the software is free for all its users. This+General Public License applies to most of the Free Software+Foundation's software and to any other program whose authors commit to+using it. (Some other Free Software Foundation software is covered by+the GNU Lesser General Public License instead.) You can apply it to+your programs, too.++ When we speak of free software, we are referring to freedom, not+price. Our General Public Licenses are designed to make sure that you+have the freedom to distribute copies of free software (and charge for+this service if you wish), that you receive source code or can get it+if you want it, that you can change the software or use pieces of it+in new free programs; and that you know you can do these things.++ To protect your rights, we need to make restrictions that forbid+anyone to deny you these rights or to ask you to surrender the rights.+These restrictions translate to certain responsibilities for you if you+distribute copies of the software, or if you modify it.++ For example, if you distribute copies of such a program, whether+gratis or for a fee, you must give the recipients all the rights that+you have. You must make sure that they, too, receive or can get the+source code. And you must show them these terms so they know their+rights.++ We protect your rights with two steps: (1) copyright the software, and+(2) offer you this license which gives you legal permission to copy,+distribute and/or modify the software.++ Also, for each author's protection and ours, we want to make certain+that everyone understands that there is no warranty for this free+software. If the software is modified by someone else and passed on, we+want its recipients to know that what they have is not the original, so+that any problems introduced by others will not reflect on the original+authors' reputations.++ Finally, any free program is threatened constantly by software+patents. We wish to avoid the danger that redistributors of a free+program will individually obtain patent licenses, in effect making the+program proprietary. To prevent this, we have made it clear that any+patent must be licensed for everyone's free use or not licensed at all.++ The precise terms and conditions for copying, distribution and+modification follow.++ GNU GENERAL PUBLIC LICENSE+ TERMS AND CONDITIONS FOR COPYING, DISTRIBUTION AND MODIFICATION++ 0. This License applies to any program or other work which contains+a notice placed by the copyright holder saying it may be distributed+under the terms of this General Public License. The "Program", below,+refers to any such program or work, and a "work based on the Program"+means either the Program or any derivative work under copyright law:+that is to say, a work containing the Program or a portion of it,+either verbatim or with modifications and/or translated into another+language. (Hereinafter, translation is included without limitation in+the term "modification".) Each licensee is addressed as "you".++Activities other than copying, distribution and modification are not+covered by this License; they are outside its scope. The act of+running the Program is not restricted, and the output from the Program+is covered only if its contents constitute a work based on the+Program (independent of having been made by running the Program).+Whether that is true depends on what the Program does.++ 1. You may copy and distribute verbatim copies of the Program's+source code as you receive it, in any medium, provided that you+conspicuously and appropriately publish on each copy an appropriate+copyright notice and disclaimer of warranty; keep intact all the+notices that refer to this License and to the absence of any warranty;+and give any other recipients of the Program a copy of this License+along with the Program.++You may charge a fee for the physical act of transferring a copy, and+you may at your option offer warranty protection in exchange for a fee.++ 2. You may modify your copy or copies of the Program or any portion+of it, thus forming a work based on the Program, and copy and+distribute such modifications or work under the terms of Section 1+above, provided that you also meet all of these conditions:++ a) You must cause the modified files to carry prominent notices+ stating that you changed the files and the date of any change.++ b) You must cause any work that you distribute or publish, that in+ whole or in part contains or is derived from the Program or any+ part thereof, to be licensed as a whole at no charge to all third+ parties under the terms of this License.++ c) If the modified program normally reads commands interactively+ when run, you must cause it, when started running for such+ interactive use in the most ordinary way, to print or display an+ announcement including an appropriate copyright notice and a+ notice that there is no warranty (or else, saying that you provide+ a warranty) and that users may redistribute the program under+ these conditions, and telling the user how to view a copy of this+ License. (Exception: if the Program itself is interactive but+ does not normally print such an announcement, your work based on+ the Program is not required to print an announcement.)++These requirements apply to the modified work as a whole. If+identifiable sections of that work are not derived from the Program,+and can be reasonably considered independent and separate works in+themselves, then this License, and its terms, do not apply to those+sections when you distribute them as separate works. But when you+distribute the same sections as part of a whole which is a work based+on the Program, the distribution of the whole must be on the terms of+this License, whose permissions for other licensees extend to the+entire whole, and thus to each and every part regardless of who wrote it.++Thus, it is not the intent of this section to claim rights or contest+your rights to work written entirely by you; rather, the intent is to+exercise the right to control the distribution of derivative or+collective works based on the Program.++In addition, mere aggregation of another work not based on the Program+with the Program (or with a work based on the Program) on a volume of+a storage or distribution medium does not bring the other work under+the scope of this License.++ 3. You may copy and distribute the Program (or a work based on it,+under Section 2) in object code or executable form under the terms of+Sections 1 and 2 above provided that you also do one of the following:++ a) Accompany it with the complete corresponding machine-readable+ source code, which must be distributed under the terms of Sections+ 1 and 2 above on a medium customarily used for software interchange; or,++ b) Accompany it with a written offer, valid for at least three+ years, to give any third party, for a charge no more than your+ cost of physically performing source distribution, a complete+ machine-readable copy of the corresponding source code, to be+ distributed under the terms of Sections 1 and 2 above on a medium+ customarily used for software interchange; or,++ c) Accompany it with the information you received as to the offer+ to distribute corresponding source code. (This alternative is+ allowed only for noncommercial distribution and only if you+ received the program in object code or executable form with such+ an offer, in accord with Subsection b above.)++The source code for a work means the preferred form of the work for+making modifications to it. For an executable work, complete source+code means all the source code for all modules it contains, plus any+associated interface definition files, plus the scripts used to+control compilation and installation of the executable. However, as a+special exception, the source code distributed need not include+anything that is normally distributed (in either source or binary+form) with the major components (compiler, kernel, and so on) of the+operating system on which the executable runs, unless that component+itself accompanies the executable.++If distribution of executable or object code is made by offering+access to copy from a designated place, then offering equivalent+access to copy the source code from the same place counts as+distribution of the source code, even though third parties are not+compelled to copy the source along with the object code.++ 4. You may not copy, modify, sublicense, or distribute the Program+except as expressly provided under this License. Any attempt+otherwise to copy, modify, sublicense or distribute the Program is+void, and will automatically terminate your rights under this License.+However, parties who have received copies, or rights, from you under+this License will not have their licenses terminated so long as such+parties remain in full compliance.++ 5. You are not required to accept this License, since you have not+signed it. However, nothing else grants you permission to modify or+distribute the Program or its derivative works. These actions are+prohibited by law if you do not accept this License. Therefore, by+modifying or distributing the Program (or any work based on the+Program), you indicate your acceptance of this License to do so, and+all its terms and conditions for copying, distributing or modifying+the Program or works based on it.++ 6. Each time you redistribute the Program (or any work based on the+Program), the recipient automatically receives a license from the+original licensor to copy, distribute or modify the Program subject to+these terms and conditions. You may not impose any further+restrictions on the recipients' exercise of the rights granted herein.+You are not responsible for enforcing compliance by third parties to+this License.++ 7. If, as a consequence of a court judgment or allegation of patent+infringement or for any other reason (not limited to patent issues),+conditions are imposed on you (whether by court order, agreement or+otherwise) that contradict the conditions of this License, they do not+excuse you from the conditions of this License. If you cannot+distribute so as to satisfy simultaneously your obligations under this+License and any other pertinent obligations, then as a consequence you+may not distribute the Program at all. For example, if a patent+license would not permit royalty-free redistribution of the Program by+all those who receive copies directly or indirectly through you, then+the only way you could satisfy both it and this License would be to+refrain entirely from distribution of the Program.++If any portion of this section is held invalid or unenforceable under+any particular circumstance, the balance of the section is intended to+apply and the section as a whole is intended to apply in other+circumstances.++It is not the purpose of this section to induce you to infringe any+patents or other property right claims or to contest validity of any+such claims; this section has the sole purpose of protecting the+integrity of the free software distribution system, which is+implemented by public license practices. Many people have made+generous contributions to the wide range of software distributed+through that system in reliance on consistent application of that+system; it is up to the author/donor to decide if he or she is willing+to distribute software through any other system and a licensee cannot+impose that choice.++This section is intended to make thoroughly clear what is believed to+be a consequence of the rest of this License.++ 8. If the distribution and/or use of the Program is restricted in+certain countries either by patents or by copyrighted interfaces, the+original copyright holder who places the Program under this License+may add an explicit geographical distribution limitation excluding+those countries, so that distribution is permitted only in or among+countries not thus excluded. In such case, this License incorporates+the limitation as if written in the body of this License.++ 9. The Free Software Foundation may publish revised and/or new versions+of the General Public License from time to time. Such new versions will+be similar in spirit to the present version, but may differ in detail to+address new problems or concerns.++Each version is given a distinguishing version number. If the Program+specifies a version number of this License which applies to it and "any+later version", you have the option of following the terms and conditions+either of that version or of any later version published by the Free+Software Foundation. If the Program does not specify a version number of+this License, you may choose any version ever published by the Free Software+Foundation.++ 10. If you wish to incorporate parts of the Program into other free+programs whose distribution conditions are different, write to the author+to ask for permission. For software which is copyrighted by the Free+Software Foundation, write to the Free Software Foundation; we sometimes+make exceptions for this. Our decision will be guided by the two goals+of preserving the free status of all derivatives of our free software and+of promoting the sharing and reuse of software generally.++ NO WARRANTY++ 11. BECAUSE THE PROGRAM IS LICENSED FREE OF CHARGE, THERE IS NO WARRANTY+FOR THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE LAW. EXCEPT WHEN+OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR OTHER PARTIES+PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY OF ANY KIND, EITHER EXPRESSED+OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF+MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE ENTIRE RISK AS+TO THE QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU. SHOULD THE+PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF ALL NECESSARY SERVICING,+REPAIR OR CORRECTION.++ 12. IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING+WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MAY MODIFY AND/OR+REDISTRIBUTE THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES,+INCLUDING ANY GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING+OUT OF THE USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED+TO LOSS OF DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY+YOU OR THIRD PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER+PROGRAMS), EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE+POSSIBILITY OF SUCH DAMAGES.++ END OF TERMS AND CONDITIONS++ How to Apply These Terms to Your New Programs++ If you develop a new program, and you want it to be of the greatest+possible use to the public, the best way to achieve this is to make it+free software which everyone can redistribute and change under these terms.++ To do so, attach the following notices to the program. It is safest+to attach them to the start of each source file to most effectively+convey the exclusion of warranty; and each file should have at least+the "copyright" line and a pointer to where the full notice is found.++ <one line to give the program's name and a brief idea of what it does.>+ Copyright (C) <year> <name of author>++ This program is free software; you can redistribute it and/or modify+ it under the terms of the GNU General Public License as published by+ the Free Software Foundation; either version 2 of the License, or+ (at your option) any later version.++ This program is distributed in the hope that it will be useful,+ but WITHOUT ANY WARRANTY; without even the implied warranty of+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the+ GNU General Public License for more details.++ You should have received a copy of the GNU General Public License along+ with this program; if not, write to the Free Software Foundation, Inc.,+ 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.++Also add information on how to contact you by electronic and paper mail.++If the program is interactive, make it output a short notice like this+when it starts in an interactive mode:++ Gnomovision version 69, Copyright (C) year name of author+ Gnomovision comes with ABSOLUTELY NO WARRANTY; for details type `show w'.+ This is free software, and you are welcome to redistribute it+ under certain conditions; type `show c' for details.++The hypothetical commands `show w' and `show c' should show the appropriate+parts of the General Public License. Of course, the commands you use may+be called something other than `show w' and `show c'; they could even be+mouse-clicks or menu items--whatever suits your program.++You should also get your employer (if you work as a programmer) or your+school, if any, to sign a "copyright disclaimer" for the program, if+necessary. Here is a sample; alter the names:++ Yoyodyne, Inc., hereby disclaims all copyright interest in the program+ `Gnomovision' (which makes passes at compilers) written by James Hacker.++ <signature of Ty Coon>, 1 April 1989+ Ty Coon, President of Vice++This General Public License does not permit incorporating your program into+proprietary programs. If your program is a subroutine library, you may+consider it more useful to permit linking proprietary applications with the+library. If this is what you want to do, use the GNU Lesser General+Public License instead of this License.
+ README view
@@ -0,0 +1,5 @@+This package contains a fast C++ implementation of the 'Tensor' interface for+Lol, via Haskell's FFI.++You can test this package by running `stack test lol-cpp`, and benchmark by+running `stack bench lol-cpp`.
+ Setup.hs view
@@ -0,0 +1,2 @@+import Distribution.Simple+main = defaultMain
+ benchmarks/BenchCPPMain.hs view
@@ -0,0 +1,83 @@+{-|+Module : BenchCPPMain+Description : Main driver for CPP benchmarks.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX++Main driver for CPP benchmarks.+-}++{-# LANGUAGE DataKinds #-}+{-# LANGUAGE NamedFieldPuns #-}+{-# LANGUAGE TypeOperators #-}++module BenchCPPMain where++import Crypto.Lol.Benchmarks+import Crypto.Lol.Benchmarks.Standard+import Crypto.Lol.Cyclotomic.Tensor.CPP+import Crypto.Lol.Factored+import qualified Crypto.Lol.Utils.PrettyPrint.Diagnostic as D+import qualified Crypto.Lol.Utils.PrettyPrint.Table as T+import Crypto.Random.DRBG++import Data.Proxy++-- choose which layers of Lol to benchmark+ls :: [String]+ls = [+ "STensor",+ "Tensor",+ "SUCyc",+ "UCyc",+ "Cyc"+ ]++-- choose which operations to benchmark+bs :: [String]+bs = [+ "unzipPow",+ "unzipDec",+ "unzipCRT",+ "zipWith (*)",+ "crt",+ "crtInv",+ "l",+ "lInv",+ "*g Pow",+ "*g Dec",+ "*g CRT",+ "divg Pow",+ "divg Dec",+ "divg CRT",+ "lift",+ "error",+ "twacePow",+ "twaceDec",+ "twaceCRT",+ "embedPow",+ "embedDec",+ "embedCRT"+ ]++main :: IO ()+main = diagnosticMain++tableMain :: IO ()+tableMain = do+ let opts = (T.defaultOpts $ Just "UCyc"){T.benches=bs}+ g1 <- defaultBenches (Proxy::Proxy CT)+ mapM_ (T.prettyBenches opts) g1++diagnosticMain :: IO ()+diagnosticMain = do+ let opts = D.defaultOpts{D.levels=ls, D.benches=bs}+ b1 <- benchGroup "Single Index"+ [oneIdxBenches (Proxy::Proxy '(F64*F9*F25, Zq 14401)) (Proxy::Proxy CT) (Proxy::Proxy HashDRBG)]+ b2 <- benchGroup "Twace-Embed"+ [twoIdxBenches (Proxy::Proxy '(F64*F9*F25, F64*F9*F25, Zq 14401)) (Proxy::Proxy CT)]+ mapM_ (D.prettyBenches opts) [b1,b2]
+ lol-cpp.cabal view
@@ -0,0 +1,135 @@+name: lol-cpp+-- The package version. See the Haskell package versioning policy (PVP)+-- for standards guiding when and how versions should be incremented.+-- http://www.haskell.org/haskellwiki/Package_versioning_policy+-- PVP summary: +-+------- breaking API changes+-- | | +----- non-breaking API additions+-- | | | +--- code changes with no API change+version: 0.0.0.1+synopsis: A fast C++ backend for <https://hackage.haskell.org/package/lol Λ ∘ λ>.+homepage: https://github.com/cpeikert/Lol+Bug-Reports: https://github.com/cpeikert/Lol/issues+license: GPL-2+license-file: LICENSE+author: Eric Crockett <ecrockett0@gmail.com>, Chris Peikert <cpeikert@alum.mit.edu>+maintainer: Eric Crockett <ecrockett0@gmail.com>+copyright: Eric Crockett, Chris Peikert+category: Crypto+stability: experimental+build-type: Simple+extra-source-files: README, CHANGES.md,+ Crypto/Lol/Cyclotomic/Tensor/CPP/*.h,+ Crypto/Lol/Cyclotomic/Tensor/CPP/*.cpp+cabal-version: >= 1.10+description:+ Λ ∘ λ (Lol) is a general-purpose library for ring-based lattice cryptography.+ This package provides a C++ implementation of Lol's Tensor interface.+source-repository head+ type: git+ location: https://github.com/cpeikert/Lol++-- For information on compiling C with cabal: http://blog.ezyang.com/2010/06/setting-up-cabal-the-ffi-and-c2hs/++Flag llvm+ Description: Compile via LLVM. This produces much better object code,+ but you need to have the LLVM compiler installed.+ -- If you enable this and get errors like "Error: can't resolve `.rodata' {.rodata section}"+ -- then GHC doesn't like your version of LLVM!+ Default: False++Flag opt+ Description: Turn on library optimizations+ Default: True++library+ default-language: Haskell2010+ ghc-options: -fwarn-dodgy-imports+ cc-options: -std=c++11+ Include-dirs: Crypto/Lol/Cyclotomic/Tensor/CPP+ -- Due to #12152, the file containing the definition of `Zq::q` must be linked first,+ -- otherwise dynamic linking (`cabal repl` or `stack ghci`) results in the error:+ -- "Loading temp shared object failed: /tmp/ghc54651_0/libghc_1.so: undefined symbol _ZN2Zq1qE"+ -- For `cabal repl`, we can simply reorder the list so that the file that should be linked+ -- first comes first in the list. However `stack ghci` always links alphabetically,+ -- so we really just have to define `Zq::q` in the first file alphabetically.+ C-sources: Crypto/Lol/Cyclotomic/Tensor/CPP/common.cpp,+ Crypto/Lol/Cyclotomic/Tensor/CPP/crt.cpp,+ Crypto/Lol/Cyclotomic/Tensor/CPP/g.cpp,+ Crypto/Lol/Cyclotomic/Tensor/CPP/l.cpp,+ Crypto/Lol/Cyclotomic/Tensor/CPP/mul.cpp,+ Crypto/Lol/Cyclotomic/Tensor/CPP/norm.cpp,+ Crypto/Lol/Cyclotomic/Tensor/CPP/random.cpp,+ Crypto/Lol/Cyclotomic/Tensor/CPP/zq.cpp++ if flag(llvm)+ ghc-options: -fllvm -optlo-O3++ -- ghc optimizations+ if flag(opt)+ -- makes lift much faster!+ ghc-options: -funfolding-use-threshold1000+ exposed-modules:+ Crypto.Lol.Cyclotomic.Tensor.CPP++ other-modules:+ Crypto.Lol.Cyclotomic.Tensor.CPP.Backend+ Crypto.Lol.Cyclotomic.Tensor.CPP.Extension+ Crypto.Lol.Cyclotomic.Tensor.CPP.Instances++ build-depends:+ arithmoi >= 0.4.1.3,+ base >= 4.9 && < 5,+ bytestring,+ constraints,+ containers >= 0.5.6.2,+ crypto-api,+ data-default >= 0.3.0,+ deepseq >= 1.4.1.1,+ lol >= 0.6.0.0,+ monadcryptorandom,+ MonadRandom >= 0.2,+ mtl >= 2.2.1,+ numeric-prelude >= 0.4.2,+ protocol-buffers,+ protocol-buffers-descriptor,+ random >= 1.1,+ reflection >= 1.5.1,+ repa>=3.4,+ singletons >= 1.1.2.1,+ th-desugar >= 1.5.4,+ tagged-transformer >= 0.7,+ template-haskell >= 2.2.0.0,+ transformers >= 0.4.2.0,+ vector>=0.11,+ vector-th-unbox >= 0.2.1.0++ other-extensions: TemplateHaskell++Benchmark bench-lol-cpp+ type: exitcode-stdio-1.0+ default-language: Haskell2010+ main-is: BenchCPPMain.hs+ ghc-options: -main-is BenchCPPMain+ hs-source-dirs: benchmarks++ ghc-options: -O2 -funfolding-creation-threshold=15000 -funfolding-use-threshold=1000++ build-depends:+ base >= 4.9 && < 5,+ DRBG,+ lol >= 0.6.0.0,+ lol-benches,+ lol-cpp++test-suite test-lol-cpp+ type: exitcode-stdio-1.0+ default-language: Haskell2010+ main-is: TestCPPMain.hs+ ghc-options: -main-is TestCPPMain+ hs-source-dirs: tests+ ghc-options: -threaded -O2++ build-depends:+ base >= 4.9 && < 5,+ lol-cpp,+ lol-tests
+ tests/TestCPPMain.hs view
@@ -0,0 +1,21 @@+{-|+Module : TestCPPMain+Description : Main driver for CPP tests.+Copyright : (c) Eric Crockett, 2011-2017+ Chris Peikert, 2011-2017+License : GPL-2+Maintainer : ecrockett0@email.com+Stability : experimental+Portability : POSIX++Main driver for CPP tests.+-}++module TestCPPMain where++import Crypto.Lol.Cyclotomic.Tensor.CPP+import Crypto.Lol.Tests.Standard+import Data.Proxy++main :: IO ()+main = defaultTestMain (Proxy::Proxy CT)