packages feed

data-sketches-core (empty) → 0.1.0.0

raw patch · 15 files changed

+1223/−0 lines, 15 filesdep +basedep +data-sketches-coredep +deepseqsetup-changed

Dependencies added: base, data-sketches-core, deepseq, ghc-prim, mwc-random, primitive, vector, vector-algorithms

Files

+ ChangeLog.md view
@@ -0,0 +1,3 @@+# Changelog for data-sketches-core++## Unreleased changes
+ LICENSE view
@@ -0,0 +1,30 @@+Copyright Author name here (c) 2021++All rights reserved.++Redistribution and use in source and binary forms, with or without+modification, are permitted provided that the following conditions are met:++    * Redistributions of source code must retain the above copyright+      notice, this list of conditions and the following disclaimer.++    * Redistributions in binary form must reproduce the above+      copyright notice, this list of conditions and the following+      disclaimer in the documentation and/or other materials provided+      with the distribution.++    * Neither the name of Author name here nor the names of other+      contributors may be used to endorse or promote products derived+      from this software without specific prior written permission.++THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS+"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT+LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR+A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT+OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,+DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY+THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT+(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE+OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ README.md view
@@ -0,0 +1,1 @@+# data-sketches-core
+ Setup.hs view
@@ -0,0 +1,2 @@+import Distribution.Simple+main = defaultMain
+ data-sketches-core.cabal view
@@ -0,0 +1,84 @@+cabal-version: 1.12++-- This file has been generated from package.yaml by hpack version 0.34.4.+--+-- see: https://github.com/sol/hpack++name:           data-sketches-core+version:        0.1.0.0+description:    Please see the README on GitHub at <https://github.com/iand675/datasketches-haskell#readme>+homepage:       https://github.com/iand675/datasketches-haskell#readme+bug-reports:    https://github.com/iand675/datasketches-haskell/issues+author:         Ian Duncan+maintainer:     ian@iankduncan.com+copyright:      2021 Ian Duncan, Rob Bassi, Mercury Technologies+license:        BSD3+license-file:   LICENSE+build-type:     Simple+extra-source-files:+    README.md+    ChangeLog.md++source-repository head+  type: git+  location: https://github.com/iand675/datasketches-haskell++library+  exposed-modules:+      DataSketches.Core.Internal.URef+      DataSketches.Core.Snapshot+      DataSketches.Quantiles.RelativeErrorQuantile.Internal+      DataSketches.Quantiles.RelativeErrorQuantile.Internal.Auxiliary+      DataSketches.Quantiles.RelativeErrorQuantile.Internal.Compactor+      DataSketches.Quantiles.RelativeErrorQuantile.Internal.Constants+      DataSketches.Quantiles.RelativeErrorQuantile.Internal.DoubleBuffer+      DataSketches.Quantiles.RelativeErrorQuantile.Internal.InequalitySearch+      DataSketches.Quantiles.RelativeErrorQuantile.Types+  other-modules:+      Paths_data_sketches_core+  hs-source-dirs:+      src+  default-extensions:+      BangPatterns+      FlexibleInstances+      RecordWildCards+      ScopedTypeVariables+      StandaloneDeriving+      TypeFamilies+      TypeOperators+  build-depends:+      base >=4.7 && <5+    , deepseq+    , ghc-prim+    , mwc-random+    , primitive+    , vector+    , vector-algorithms+  default-language: Haskell2010++test-suite data-sketches-core-test+  type: exitcode-stdio-1.0+  main-is: Spec.hs+  other-modules:+      Paths_data_sketches_core+  hs-source-dirs:+      test+  default-extensions:+      BangPatterns+      FlexibleInstances+      RecordWildCards+      ScopedTypeVariables+      StandaloneDeriving+      TypeFamilies+      TypeOperators+  ghc-options: -threaded -rtsopts -with-rtsopts=-N+  build-depends:+      base >=4.7 && <5+    , data-sketches-core+    , deepseq+    , ghc-prim+    , mwc-random+    , primitive+    , vector+    , vector-algorithms+  default-language: Haskell2010
+ src/DataSketches/Core/Internal/URef.hs view
@@ -0,0 +1,49 @@+module DataSketches.Core.Internal.URef where++import Control.Monad.Primitive+import qualified Data.Vector.Unboxed.Mutable as MUVector+import Data.Vector.Unboxed (Unbox)++-- | An unboxed reference. This works like an 'IORef', but the data is+-- stored in a bytearray instead of a heap object, avoiding+-- significant allocation overhead in some cases. For a concrete+-- example, see this Stack Overflow question:+-- <https://stackoverflow.com/questions/27261813/why-is-my-little-stref-int-require-allocating-gigabytes>.+--+-- The first parameter is the state token type, the same as would be+-- used for the 'ST' monad. If you're using an 'IO'-based monad, you+-- can use the convenience 'IOURef' type synonym instead.+--+-- @since 0.0.2.0+newtype URef s a = URef (MUVector.MVector s a)++-- | Helpful type synonym for using a 'URef' from an 'IO'-based stack.+--+-- @since 0.0.2.0+type IOURef = URef (PrimState IO)++-- | Create a new 'URef'+--+-- @since 0.0.2.0+newURef :: (PrimMonad m, Unbox a) => a -> m (URef (PrimState m) a)+newURef a = fmap URef (MUVector.replicate 1 a)++-- | Read the value in a 'URef'+--+-- @since 0.0.2.0+readURef :: (PrimMonad m, Unbox a) => URef (PrimState m) a -> m a+readURef (URef v) = MUVector.read v 0++-- | Write a value into a 'URef'. Note that this action is strict, and+-- will force evalution of the value.+--+-- @since 0.0.2.0+writeURef :: (PrimMonad m, Unbox a) => URef (PrimState m) a -> a -> m ()+writeURef (URef v) = MUVector.unsafeWrite v 0++-- | Modify a value in a 'URef'. Note that this action is strict, and+-- will force evaluation of the result value.+--+-- @since 0.0.2.0+modifyURef :: (PrimMonad m, Unbox a) => URef (PrimState m) a -> (a -> a) -> m ()+modifyURef u f = readURef u >>= writeURef u . f
+ src/DataSketches/Core/Snapshot.hs view
@@ -0,0 +1,7 @@+module DataSketches.Core.Snapshot where++import Control.Monad.Primitive++class TakeSnapshot a where+  type Snapshot a+  takeSnapshot :: PrimMonad m => a (PrimState m) -> m (Snapshot a)
+ src/DataSketches/Quantiles/RelativeErrorQuantile/Internal.hs view
@@ -0,0 +1,130 @@+{-# LANGUAGE DeriveGeneric #-}+module DataSketches.Quantiles.RelativeErrorQuantile.Internal where++import Data.Primitive (MutVar, readMutVar)+import Data.Word+import qualified Data.Vector as Vector+import GHC.Generics+import System.Random.MWC (Gen)++import DataSketches.Core.Snapshot+import DataSketches.Quantiles.RelativeErrorQuantile.Types+import DataSketches.Core.Internal.URef (URef, readURef)+import DataSketches.Quantiles.RelativeErrorQuantile.Internal.Auxiliary (ReqAuxiliary)+import DataSketches.Quantiles.RelativeErrorQuantile.Internal.Compactor (ReqCompactor)+import Control.DeepSeq (NFData, rnf)+import Control.Monad.Primitive (PrimMonad (PrimState))+import qualified DataSketches.Quantiles.RelativeErrorQuantile.Internal.Compactor as Compactor+import qualified DataSketches.Quantiles.RelativeErrorQuantile.Internal.DoubleBuffer as DoubleBuffer+import qualified DataSketches.Quantiles.RelativeErrorQuantile.Internal.Auxiliary as Auxiliary+import Control.Exception (Exception)+++{- |+This Relative Error Quantiles Sketch is the Haskell implementation based on the paper+"Relative Error Streaming Quantiles", https://arxiv.org/abs/2004.01668, and loosely derived from+a Python prototype written by Pavel Vesely, ported from the Java equivalent.++This implementation differs from the algorithm described in the paper in the following:++The algorithm requires no upper bound on the stream length.+Instead, each relative-compactor counts the number of compaction operations performed+so far (via variable state). Initially, the relative-compactor starts with INIT_NUMBER_OF_SECTIONS.+Each time the number of compactions (variable state) exceeds 2^{numSections - 1}, we double+numSections. Note that after merging the sketch with another one variable state may not correspond+to the number of compactions performed at a particular level, however, since the state variable+never exceeds the number of compactions, the guarantees of the sketch remain valid.++The size of each section (variable k and sectionSize in the code and parameter k in+the paper) is initialized with a value set by the user via variable k.+When the number of sections doubles, we decrease sectionSize by a factor of sqrt(2).+This is applied at each level separately. Thus, when we double the number of sections, the+nominal compactor size increases by a factor of approx. sqrt(2) (+/- rounding).++The merge operation here does not perform "special compactions", which are used in the paper+to allow for a tight mathematical analysis of the sketch.++This implementation provides a number of capabilities not discussed in the paper or provided+in the Python prototype.++The Python prototype only implemented high accuracy for low ranks. This implementation+provides the user with the ability to choose either high rank accuracy or low rank accuracy at+the time of sketch construction.++- The Python prototype only implemented a comparison criterion of "<". This implementation+allows the user to switch back and forth between the "<=" criterion and the "<=" criterion.+-}+data ReqSketch s = ReqSketch+  { k :: !Word32+  , rankAccuracySetting :: !RankAccuracy+  , criterion :: !Criterion+  , sketchRng :: {-# UNPACK #-} !(Gen s)+  , totalN :: {-# UNPACK #-} !(URef s Word64)+  , minValue :: {-# UNPACK #-} !(URef s Double)+  , maxValue :: {-# UNPACK #-} !(URef s Double)+  , sumValue :: {-# UNPACK #-} !(URef s Double)+  , retainedItems :: {-# UNPACK #-} !(URef s Int)+  , maxNominalCapacitiesSize :: {-# UNPACK #-} !(URef s Int)+  , aux :: {-# UNPACK #-} !(MutVar s (Maybe ReqAuxiliary))+  , compactors :: {-# UNPACK #-} !(MutVar s (Vector.Vector (ReqCompactor s)))+  } deriving (Generic)++instance NFData (ReqSketch s) where+  rnf !rs = ()++data ReqSketchSnapshot = ReqSketchSnapshot+    { snapshotRankAccuracySetting :: !RankAccuracy+    , snapshotCriterion :: !Criterion+    , snapshotTotalN :: !Word64+    , snapshotMinValue :: !Double+    , snapshotMaxValue :: !Double+    , snapshotRetainedItems :: !Int+    , snapshotMaxNominalCapacitiesSize :: !Int+    -- , aux :: !(MutVar s (Maybe ()))+    , snapshotCompactors :: !(Vector.Vector (Snapshot ReqCompactor))+    } deriving Show++instance TakeSnapshot ReqSketch where+  type Snapshot ReqSketch = ReqSketchSnapshot +  takeSnapshot ReqSketch{..} = ReqSketchSnapshot rankAccuracySetting criterion+    <$> readURef totalN+    <*> readURef minValue+    <*> readURef maxValue+    <*> readURef retainedItems+    <*> readURef maxNominalCapacitiesSize+    <*> (readMutVar compactors >>= mapM takeSnapshot)++getCompactors :: PrimMonad m => ReqSketch (PrimState m) -> m (Vector.Vector (ReqCompactor (PrimState m)))+getCompactors = readMutVar . compactors++computeTotalRetainedItems :: PrimMonad m => ReqSketch (PrimState m) -> m Int+computeTotalRetainedItems this = do+  compactors <- getCompactors this+  Vector.foldM countBuffer 0 compactors+  where+    countBuffer acc compactor = do+      buff <- Compactor.getBuffer compactor+      buffSize <- DoubleBuffer.getCount buff+      pure $ buffSize + acc++retainedItemCount :: PrimMonad m => ReqSketch (PrimState m) -> m Int+retainedItemCount = readURef . retainedItems++-- | Get the total number of items inserted into the sketch+count :: PrimMonad m => ReqSketch (PrimState m) -> m Word64+count = readURef . totalN++mkAuxiliaryFromReqSketch :: PrimMonad m => ReqSketch (PrimState m) -> m ReqAuxiliary+mkAuxiliaryFromReqSketch this = do+  total <- count this+  retainedItems <- retainedItemCount this+  compactors <- getCompactors this+  Auxiliary.mkAuxiliary (rankAccuracySetting this) total retainedItems compactors++data CumulativeDistributionInvariants+  = CumulativeDistributionInvariantsSplitsAreEmpty+  | CumulativeDistributionInvariantsSplitsAreNotFinite+  | CumulativeDistributionInvariantsSplitsAreNotUniqueAndMontonicallyIncreasing+  deriving (Show, Eq)++instance Exception CumulativeDistributionInvariants
+ src/DataSketches/Quantiles/RelativeErrorQuantile/Internal/Auxiliary.hs view
@@ -0,0 +1,189 @@+module DataSketches.Quantiles.RelativeErrorQuantile.Internal.Auxiliary+  ( ReqAuxiliary(..)+  , MReqAuxiliary (..)+  , mkAuxiliary+  , getQuantile+  -- | Really extra private, just needed for tests+  , mergeSortIn+  ) where++import GHC.TypeLits+import Control.Monad (when)+import Control.Monad.Primitive+import Data.Bits (shiftL)+import Data.Word+import Data.Primitive.MutVar+import Data.Vector.Algorithms.Search+import qualified Data.Vector as Vector+import qualified Data.Vector.Unboxed.Mutable as MUVector+import DataSketches.Quantiles.RelativeErrorQuantile.Types+import DataSketches.Quantiles.RelativeErrorQuantile.Internal.Compactor (ReqCompactor)+import qualified DataSketches.Quantiles.RelativeErrorQuantile.Internal.Compactor as Compactor+import DataSketches.Quantiles.RelativeErrorQuantile.Internal.DoubleBuffer (DoubleBuffer)+import qualified DataSketches.Quantiles.RelativeErrorQuantile.Internal.DoubleBuffer as DoubleBuffer+import qualified Data.Vector.Unboxed as U+import Control.Monad.ST+import DataSketches.Quantiles.RelativeErrorQuantile.Internal.InequalitySearch (find)+import qualified DataSketches.Quantiles.RelativeErrorQuantile.Internal.InequalitySearch as IS+import Debug.Trace+import qualified Data.Vector.Generic.Mutable as MG++data ReqAuxiliary = ReqAuxiliary+  { raWeightedItems :: {-# UNPACK #-} !(U.Vector (Double, Word64))+  , raHighRankAccuracy :: !RankAccuracy+  , raSize :: {-# UNPACK #-} !Word64+  }+  deriving (Show, Eq)++data MReqAuxiliary s = MReqAuxiliary+  { mraWeightedItems :: {-# UNPACK #-} !(MutVar s (MUVector.MVector s (Double, Word64)))+  , mraHighRankAccuracy :: !RankAccuracy+  , mraSize :: {-# UNPACK #-} !Word64+  }++mkAuxiliary :: (PrimMonad m, s ~ PrimState m) => RankAccuracy -> Word64 -> Int -> Vector.Vector (ReqCompactor s) -> m ReqAuxiliary+mkAuxiliary rankAccuracy totalN retainedItems compactors = do+  items <- newMutVar =<< MUVector.replicate retainedItems (0, 0)+  let this = MReqAuxiliary+        { mraWeightedItems = items+        , mraHighRankAccuracy = rankAccuracy+        , mraSize = totalN+        }+  Vector.foldM_ (mergeBuffers this) 0 compactors+  createCumulativeWeights this+  dedup this+  items' <- U.unsafeFreeze =<< readMutVar items+  pure ReqAuxiliary+    { raWeightedItems = items'+    , raHighRankAccuracy = rankAccuracy+    , raSize = totalN+    }+  where+    mergeBuffers this auxCount compactor = do+      buff <- Compactor.getBuffer compactor+      buffSize <-  DoubleBuffer.getCount buff+      let lgWeight = Compactor.getLgWeight compactor+          weight = 1 `shiftL` fromIntegral lgWeight+      mergeSortIn this buff weight auxCount+      pure $ auxCount + buffSize++getWeightedItems :: PrimMonad m => MReqAuxiliary (PrimState m) -> m (MUVector.MVector (PrimState m) (Double, Word64))+getWeightedItems = readMutVar . mraWeightedItems++getItems :: PrimMonad m => MReqAuxiliary (PrimState m) -> m (MUVector.MVector (PrimState m) Double)+getItems = fmap (fst . MUVector.unzip) . getWeightedItems++getWeights :: PrimMonad m => MReqAuxiliary (PrimState m) -> m (MUVector.MVector (PrimState m) Word64)+getWeights = fmap (snd . MUVector.unzip) . getWeightedItems++getQuantile :: ReqAuxiliary -> Double -> Criterion  -> Double+getQuantile this normalRank ltEq = fst (weightedItems U.! ix)+  where+    ix = if searchResult == U.length weightedItems+      then searchResult - 1+      else searchResult+    searchResult = runST $ do+      v <- U.unsafeThaw $ snd $ U.unzip weightedItems+      let search = case ltEq of+            (:<) -> find (IS.:>)+            (:<=) -> find (IS.:>=)+      search v 0 (weightsSize - 1) rank+    weightedItems = raWeightedItems this+    weightsSize = U.length weightedItems+    rank = floor (normalRank * fromIntegral (raSize this))++createCumulativeWeights :: PrimMonad m => MReqAuxiliary (PrimState m) -> m ()+createCumulativeWeights this = do+  weights <- getWeights this+  let size = MUVector.length weights+  let accumulateM i weight = do+        when (i > 0) $ do+          prevWeight <- MUVector.read weights (i - 1)+          MUVector.unsafeWrite weights i (weight + prevWeight)+  forI_ weights (\i -> MUVector.read weights i >>= \x -> accumulateM i x)+  lastWeight <- MUVector.read weights (size - 1)+  when (lastWeight /= mraSize this) $ do+    error "invariant violated: lastWeight does not equal raSize"+  where+    forI_ :: (Monad m, MG.MVector v a) => v (PrimState m) a -> (Int -> m b) -> m ()+    {-# INLINE forI_ #-}+    forI_ v f = loop 0+      where+        loop i +          | i >= n    = return ()+          | otherwise = f i >> loop (i + 1)+        n = MG.length v+++dedup :: PrimMonad m => MReqAuxiliary (PrimState m) -> m ()+dedup this = do+  weightedItems <- getWeightedItems this+  let size = MUVector.length weightedItems+  weightedItemsB <- MUVector.replicate size (0, 0)+  bi <- doDedup weightedItems size weightedItemsB 0 0+  writeMutVar (mraWeightedItems this) $ MUVector.slice 0 bi weightedItemsB+  where+    doDedup weightedItems itemsSize weightedItemsB = go+      where +        go !i !bi+          | i >= itemsSize = pure bi+          | otherwise = do+            let j = i + 1+                hidup = j+                countDups !j !hidup = if j < itemsSize +                  then do+                    (itemI, _) <- MUVector.read weightedItems i+                    (itemJ, _) <- MUVector.read weightedItems j+                    if itemI == itemJ+                      then countDups (j + 1) j+                      else pure (j, hidup)+                  else pure (j, hidup)+            (j', hidup') <- countDups j hidup+            if j' - i == 1 -- no dups+              then do+                (item, weight) <- MUVector.read weightedItems i+                MUVector.unsafeWrite weightedItemsB bi (item, weight)+                go (i + 1) (bi + 1)+              else do+                (item, weight) <- MUVector.read weightedItems hidup'+                MUVector.unsafeWrite weightedItemsB bi (item, weight)+                go j' (bi + 1)++mergeSortIn :: PrimMonad m => MReqAuxiliary (PrimState m) -> DoubleBuffer (PrimState m) -> Word64 -> Int -> m ()+mergeSortIn this bufIn defaultWeight auxCount = do+  DoubleBuffer.sort bufIn+  weightedItems <- getWeightedItems this+  otherItems <- DoubleBuffer.getVector bufIn+  otherBuffSize <- DoubleBuffer.getCount bufIn+  otherBuffCapacity <- DoubleBuffer.getCapacity bufIn+  let totalSize = otherBuffSize + auxCount - 1+      height = case mraHighRankAccuracy this of+        HighRanksAreAccurate -> otherBuffCapacity - 1+        LowRanksAreAccurate -> otherBuffSize - 1+  merge totalSize weightedItems otherItems (auxCount - 1) (otherBuffSize - 1) height +  where+    merge totalSize weightedItems otherItems = go totalSize+      where+        go !k !i !j !h +          | k < 0 = pure ()+          | i >= 0 && j >= 0 = do+            (item, weight) <- MUVector.read weightedItems i+            otherItem <- MUVector.read otherItems h+            if item >= otherItem+               then do+                 MUVector.unsafeWrite weightedItems k (item, weight)+                 continue (i - 1) j h+               else do+                 MUVector.unsafeWrite weightedItems k (otherItem, defaultWeight)+                 continue i (j - 1) (h - 1)+          | i >= 0 = do+            MUVector.read weightedItems i >>= MUVector.write weightedItems k+            continue (i - 1) j h+          | j >= 0 = do+            otherItem <- MUVector.read otherItems h+            MUVector.unsafeWrite weightedItems k (otherItem, defaultWeight)+            continue i (j - 1) (h - 1)+          | otherwise = pure ()+          where+            continue = go (k - 1)+
+ src/DataSketches/Quantiles/RelativeErrorQuantile/Internal/Compactor.hs view
@@ -0,0 +1,238 @@+module DataSketches.Quantiles.RelativeErrorQuantile.Internal.Compactor+  ( ReqCompactor+  , mkReqCompactor+  , CompactorReturn (..)+  , compact+  , getBuffer+  , getCoin+  , getLgWeight+  , getNominalCapacity+  , getNumSections+  , merge+  , nearestEven+  ) where++import GHC.TypeLits+import Data.Bits ((.&.), (.|.), complement, countTrailingZeros, shiftL, shiftR)+import Data.Primitive.MutVar+import Data.Proxy+import Data.Semigroup (Semigroup)+import Data.Word+import DataSketches.Quantiles.RelativeErrorQuantile.Types+import System.Random.MWC (create, Variate(uniform), Gen)+import Control.Exception (assert)+import Control.Monad (when)+import Control.Monad.Primitive+import DataSketches.Quantiles.RelativeErrorQuantile.Internal.Constants+import DataSketches.Quantiles.RelativeErrorQuantile.Internal.DoubleBuffer+import DataSketches.Core.Internal.URef+import DataSketches.Core.Snapshot++data CompactorReturn s = CompactorReturn+  { crDeltaRetItems :: {-# UNPACK #-} !Int+  , crDeltaNominalSize :: {-# UNPACK #-} !Int+  , crDoubleBuffer :: {-# UNPACK #-} !(DoubleBuffer s)+  }++data ReqCompactor s = ReqCompactor+  -- Configuration constants+  { rcRankAccuracy :: !RankAccuracy+  , rcLgWeight :: {-# UNPACK #-} !Word8+  , rcRng :: {-# UNPACK #-} !(Gen s)+  -- State+  , rcState :: {-# UNPACK #-} !(URef s Word64)+  , rcLastFlip :: {-# UNPACK #-} !(URef s Bool)+  , rcSectionSizeFlt :: {-# UNPACK #-} !(URef s Double)+  , rcSectionSize :: {-# UNPACK #-} !(URef s Word32)+  , rcNumSections :: {-# UNPACK #-} !(URef s Word8)+  , rcBuffer :: {-# UNPACK #-} !(MutVar s (DoubleBuffer s))+  }++data ReqCompactorSnapshot = ReqCompactorSnapshot+    { snapshotCompactorRankAccuracy :: !RankAccuracy+    , snapshotCompactorRankAccuracyState :: !Word64+    , snapshotCompactorLastFlip :: !Bool+    , snapshotCompactorSectionSizeFlt :: !Double+    , snapshotCompactorSectionSize :: !Word32+    , snapshotCompactorNumSections :: !Word8+    , snapshotCompactorBuffer :: !(Snapshot DoubleBuffer)+    } deriving (Show)++instance TakeSnapshot ReqCompactor where+  type Snapshot ReqCompactor = ReqCompactorSnapshot+  takeSnapshot ReqCompactor{..} = ReqCompactorSnapshot rcRankAccuracy+    <$> readURef rcState+    <*> readURef rcLastFlip+    <*> readURef rcSectionSizeFlt+    <*> readURef rcSectionSize+    <*> readURef rcNumSections+    <*> (readMutVar rcBuffer >>= takeSnapshot)+++mkReqCompactor+  :: PrimMonad m+  => Gen (PrimState m)+  -> Word8+  -> RankAccuracy+  -> Word32+  -> m (ReqCompactor (PrimState m))+mkReqCompactor g lgWeight rankAccuracy sectionSize = do+  let nominalCapacity = fromIntegral $ nomCapMulti * initNumberOfSections * sectionSize+  buff <- mkBuffer (nominalCapacity * 2) nominalCapacity (rankAccuracy == HighRanksAreAccurate)+  ReqCompactor rankAccuracy lgWeight g+    <$> newURef 0+    <*> newURef False+    <*> newURef (fromIntegral sectionSize)+    <*> newURef sectionSize+    <*> newURef initNumberOfSections+    <*> newMutVar buff++nomCapMult :: Num a => a+nomCapMult = 2++toInt :: Integral a => a -> Int+toInt = fromIntegral++compact :: (PrimMonad m) => ReqCompactor (PrimState m) -> m (CompactorReturn (PrimState m))+compact this = do+  startBuffSize <- getCount =<< getBuffer this+  startNominalCapacity <- getNominalCapacity this+  numSections <- readURef $ rcNumSections this+  sectionSize <- readURef $ rcSectionSize this+  state <- readURef $ rcState this+  let trailingOnes = succ $ countTrailingZeros $ complement state+      sectionsToCompact = min trailingOnes $ fromIntegral numSections+  (compactionStart, compactionEnd) <- computeCompactionRange this sectionsToCompact+  -- TODO, this fails in GHCi but not in tests?+  assert (compactionEnd - compactionStart >= 2) $ do+    coin <- if state .&. 1 == 1+      then fmap not $ readURef $ rcLastFlip this+      else flipCoin this+    writeURef (rcLastFlip this) coin+    buff <- getBuffer this+    promote <- getEvensOrOdds buff compactionStart compactionEnd coin+    trimCount buff $ startBuffSize - (compactionEnd - compactionStart)+    modifyURef (rcState this) (+ 1)+    ensureEnoughSections this+    endBuffSize <- getCount buff+    promoteBuffSize <- getCount promote+    endNominalCapacity <- getNominalCapacity this+    pure $ CompactorReturn+      { crDeltaRetItems = endBuffSize - startBuffSize + promoteBuffSize+      , crDeltaNominalSize = endNominalCapacity - startNominalCapacity+      , crDoubleBuffer = promote+      }++getLgWeight :: ReqCompactor s -> Word8+getLgWeight = rcLgWeight++getBuffer :: PrimMonad m => ReqCompactor (PrimState m) -> m (DoubleBuffer (PrimState m))+getBuffer = readMutVar . rcBuffer++flipCoin :: (PrimMonad m) => ReqCompactor (PrimState m) -> m Bool+flipCoin = uniform . rcRng++getCoin :: PrimMonad m => ReqCompactor (PrimState m) -> m Bool+getCoin = readURef . rcLastFlip++getNominalCapacity :: PrimMonad m => ReqCompactor (PrimState m) -> m Int+getNominalCapacity compactor = do+  numSections <- readURef $ rcNumSections compactor+  sectionSize <- readURef $ rcSectionSize compactor+  pure $ nomCapMult * toInt numSections * toInt sectionSize++getNumSections :: PrimMonad m => ReqCompactor (PrimState m) -> m Word8+getNumSections = readURef . rcNumSections++getSectionSizeFlt :: PrimMonad m => ReqCompactor (PrimState m) -> m Double+getSectionSizeFlt = readURef . rcSectionSizeFlt++getState :: PrimMonad m => ReqCompactor (PrimState m) -> m Word64+getState = readURef . rcState++-- | Merge the other given compactor into this one. They both must have the+-- same @lgWeight@+merge+  :: (PrimMonad m, s ~ PrimState m)+  => ReqCompactor (PrimState m)+  -- ^ The compactor to merge into+  -> ReqCompactor (PrimState m)+  -- ^ The compactor to merge from +  -> m (ReqCompactor s)+merge this otherCompactor = assert (rcLgWeight this == rcLgWeight otherCompactor) $ do+  otherState <- readURef $ rcState otherCompactor+  modifyURef (rcState this) (.|. otherState)+  ensureMaxSections++  buff <- getBuffer this+  sort buff++  otherBuff <- getBuffer otherCompactor+  sort otherBuff++  otherBuffIsBigger <- (>) <$> getCount otherBuff <*> getCount buff+  finalBuff <- if otherBuffIsBigger+     then do+        otherBuff' <- copyBuffer otherBuff+        mergeSortIn otherBuff' buff+        writeMutVar (rcBuffer this) otherBuff'+     else mergeSortIn buff otherBuff+  pure this+  where+    ensureMaxSections = do+      adjusted <- ensureEnoughSections this+      -- loop until no adjustments can be made+      when adjusted ensureMaxSections++-- | Adjust the sectionSize and numSections if possible.+ensureEnoughSections+  :: PrimMonad m+  => ReqCompactor (PrimState m)+  -> m Bool+  -- ^ 'True' if the SectionSize and NumSections were adjusted.+ensureEnoughSections compactor = do+  sectionSizeFlt <- readURef $ rcSectionSizeFlt compactor+  let szf = sectionSizeFlt / sqrt2+      ne = nearestEven szf+  state <- readURef $ rcState compactor+  numSections <- readURef $ rcNumSections compactor+  sectionSize <- readURef $ rcSectionSize compactor+  if state >= (1 `shiftL` toInt (numSections - 1))+     && sectionSize > minK+     && ne >= minK+     then do+       writeURef (rcSectionSizeFlt compactor) szf+       writeURef (rcSectionSize compactor) $ fromIntegral ne+       modifyURef (rcNumSections compactor) (`shiftL` 1)+       buf <- getBuffer compactor+       nomCapacity <- getNominalCapacity compactor+       ensureCapacity buf (2 * nomCapacity)+       pure True+     else pure False++-- | Computes the start and end indices of the compacted region+computeCompactionRange+  :: PrimMonad m+  => ReqCompactor (PrimState m)+  -> Int+  -- ^ secsToCompact the number of contiguous sections to compact+  -> m (Int, Int)+-- ^ the start and end indices of the compacted region in compact form+computeCompactionRange this secsToCompact = do+  buffSize <- getCount =<< getBuffer this+  nominalCapacity <- getNominalCapacity this+  numSections <- readURef $ rcNumSections this+  sectionSize <- readURef $ rcSectionSize this+  let nonCompact = (nominalCapacity `div` 2) + (fromIntegral numSections - secsToCompact) * fromIntegral sectionSize+      nonCompact' = if (buffSize - nonCompact) .&. 1 == 1 then nonCompact - 1 else nonCompact+  pure $ case rcRankAccuracy this of+    HighRanksAreAccurate -> (0, fromIntegral $ buffSize - nonCompact')+    LowRanksAreAccurate -> (nonCompact', buffSize)++-- | Returns the nearest even integer to the given value. Also used by test.+nearestEven+  :: Double+  -- ^ the given value+  -> Int+  -- ^ the nearest even integer to the given value.+nearestEven x = round (x / 2) `shiftL` 1
+ src/DataSketches/Quantiles/RelativeErrorQuantile/Internal/Constants.hs view
@@ -0,0 +1,21 @@+module DataSketches.Quantiles.RelativeErrorQuantile.Internal.Constants where+-- Constants+import Data.Word++sqrt2 :: Double +sqrt2 = sqrt 2++initNumberOfSections :: Num a => a+initNumberOfSections = 3++minK :: Num a => a+minK = 4++nomCapMulti :: Num a => a+nomCapMulti = 2++relRseFactor :: Double+relRseFactor = sqrt (0.0512 / fromIntegral initNumberOfSections)++fixRseFactor :: Double+fixRseFactor = 0.084
+ src/DataSketches/Quantiles/RelativeErrorQuantile/Internal/DoubleBuffer.hs view
@@ -0,0 +1,321 @@+{-# LANGUAGE MagicHash #-}+{-# LANGUAGE UnboxedTuples #-}+module DataSketches.Quantiles.RelativeErrorQuantile.Internal.DoubleBuffer+  ( DoubleBuffer+  , Capacity+  , GrowthIncrement+  , SpaceAtBottom+  , DoubleIsNonFiniteException(..)+  , mkBuffer+  , copyBuffer+  , append+  , ensureCapacity+  , getCountWithCriterion+  , getEvensOrOdds+  , (!) -- getItem+  , growthIncrement+  , spaceAtBottom+  , getCapacity+  , getCount+  , getSpace+  , getVector+  , isEmpty+  , isSorted+  , sort+  , mergeSortIn+  , trimCount+  ) where++import DataSketches.Quantiles.RelativeErrorQuantile.Types+    ( Criterion )+import Control.Monad ( unless, when )+import Control.Monad.Primitive ( PrimMonad(PrimState) )+import Data.Primitive.MutVar+    ( newMutVar, readMutVar, writeMutVar, MutVar )+import qualified Data.Vector.Unboxed as UVector+import qualified Data.Vector.Unboxed.Mutable as MUVector+import DataSketches.Core.Internal.URef+    ( URef, newURef, readURef, writeURef, modifyURef )+import Data.Vector.Algorithms.Intro (sortByBounds)+import GHC.Stack ( HasCallStack )+import System.IO.Unsafe ()+import qualified DataSketches.Quantiles.RelativeErrorQuantile.Internal.InequalitySearch as IS+import Control.Exception ( Exception, throw )+import DataSketches.Core.Snapshot ( TakeSnapshot(..) )++-- | A special buffer of floats specifically designed to support the ReqCompactor class.+data DoubleBuffer s = DoubleBuffer+  { vec :: {-# UNPACK #-} !(MutVar s (MUVector.MVector s Double))+  , count :: {-# UNPACK #-} !(URef s Int)+  , sorted :: {-# UNPACK #-} !(URef s Bool)+  , growthIncrement :: {-# UNPACK #-} !Int+  , spaceAtBottom :: !Bool+  }++data DoubleBufferSnapshot = DoubleBufferSnapshot+    { dbSnapshotVec :: UVector.Vector Double+    , dbSnapshotCount :: !Int+    , dbSnapshotSorted :: !Bool+    , dbSnapshotGrowthIncrement :: !Int+    , dbSnapshotSpaceAtBottom :: !Bool+    } deriving (Show)++instance TakeSnapshot DoubleBuffer where+  type Snapshot DoubleBuffer = DoubleBufferSnapshot++  takeSnapshot DoubleBuffer{..} = DoubleBufferSnapshot+    <$> (readMutVar vec >>= UVector.freeze)+    <*> readURef count+    <*> readURef sorted+    <*> pure growthIncrement+    <*> pure spaceAtBottom++type Capacity = Int+type GrowthIncrement = Int+type SpaceAtBottom = Bool++-- | Constructs an new empty FloatBuffer with an initial capacity specified by+-- the <code>capacity</code> argument.+mkBuffer :: PrimMonad m => Capacity -> GrowthIncrement -> SpaceAtBottom -> m (DoubleBuffer (PrimState m))+mkBuffer capacity_ growthIncrement spaceAtBottom = do+  vec <- newMutVar =<< MUVector.new capacity_+  count <- newURef 0+  sorted <- newURef True+  pure $ DoubleBuffer{..}++copyBuffer :: PrimMonad m => DoubleBuffer (PrimState m) -> m (DoubleBuffer (PrimState m))+copyBuffer buf@DoubleBuffer{..} = do+  vec <- newMutVar =<< MUVector.clone =<< getVector buf+  count <- newURef =<< getCount buf+  sorted <- newURef =<< readURef sorted+  pure $ DoubleBuffer {..}++-- | Appends the given item to the active array and increments the active count.+-- This will expand the array if necessary.+append :: PrimMonad m => DoubleBuffer (PrimState m) -> Double -> m ()+append buf@DoubleBuffer{..} x = do+  ensureSpace buf 1+  index <- if spaceAtBottom+    then+      (\capacity_ count_ -> capacity_ - count_ - 1)+        <$> getCapacity buf+        <*> getCount buf+    else readURef count+  modifyURef count (+ 1)+  getVector buf >>= \vec -> MUVector.unsafeWrite vec index x+  writeURef sorted False+{-# SCC append #-}++-- | Ensures that the capacity of this FloatBuffer is at least newCapacity.+-- If newCapacity &lt; capacity(), no action is taken.+ensureSpace :: PrimMonad m => DoubleBuffer (PrimState m) -> Int -> m ()+ensureSpace buf@DoubleBuffer{..} space = do+  count_ <- readURef count+  capacity_ <- getCapacity buf+  let notEnoughSpace = count_ + space > capacity_+  when notEnoughSpace $ do+    let newCap = count_ + space + growthIncrement+    ensureCapacity buf newCap++getVector :: (PrimMonad m, PrimState m ~ s) => DoubleBuffer s -> m (MUVector.MVector s Double)+getVector = readMutVar . vec+{-# INLINE getVector #-}++getCapacity :: PrimMonad m => DoubleBuffer (PrimState m) -> m Int+getCapacity = fmap MUVector.length . getVector+{-# INLINE getCapacity #-}++ensureCapacity :: PrimMonad m => DoubleBuffer (PrimState m) -> Int -> m ()+ensureCapacity buf@DoubleBuffer{..} newCapacity = do+  capacity_ <- getCapacity buf+  when (newCapacity > capacity_) $ do+    count_ <- getCount buf+    (srcPos, destPos) <- if spaceAtBottom+      then do+        pure (capacity_ - count_, newCapacity - count_)+      else pure (0, 0)+    oldVec <- getVector buf+    newVec <- MUVector.new newCapacity+    MUVector.unsafeCopy+      (MUVector.slice destPos count_ newVec)+      (MUVector.slice srcPos count_ oldVec)+    writeMutVar vec newVec+{-# SCC ensureCapacity #-}++newtype DoubleIsNonFiniteException = DoubleIsNonFiniteException Double+  deriving (Show, Eq)++instance Exception DoubleIsNonFiniteException++getCountWithCriterion :: PrimMonad m => DoubleBuffer (PrimState m) -> Double -> Criterion -> m Int+getCountWithCriterion buf@DoubleBuffer{..} value criterion = do+  when (isNaN value || isInfinite value) $ throw $ DoubleIsNonFiniteException value+  sort buf+  count_ <- getCount buf+  vec <- getVector buf+  (low, high) <- if spaceAtBottom+    then do+      capacity_ <- getCapacity buf+      pure (capacity_ - count_, capacity_ - 1)+    else pure (0, count_)++  ix <- IS.find criterion vec low high value+  pure $! if ix == MUVector.length vec+    then 0+    else ix - low + 1++-- data EvensOrOdds = Evens | Odds++getEvensOrOdds :: PrimMonad m => DoubleBuffer (PrimState m) -> Int -> Int -> Bool -> m (DoubleBuffer (PrimState m))+getEvensOrOdds buf@DoubleBuffer{..} startOffset endOffset odds = do+  (start, end) <- if spaceAtBottom+    then do+      basis <- (-) <$> getCapacity buf <*> getCount buf+      pure (basis + startOffset, basis + endOffset)+    else pure (startOffset, endOffset)+  sort buf+  let range = endOffset - startOffset+  vec <- getVector buf+  out <- MUVector.new (range `div` 2)+  go vec out start 0+  where+    odd = if odds then 1 else 0+    go vec !out !i !j = if j < MUVector.length out+      then do+        MUVector.unsafeWrite out j =<< MUVector.unsafeRead vec (i + odd)+        go vec out (i + 2) (j + 1)+      else do+        count <- newURef (MUVector.length out)+        sorted <- newURef True+        vec <- newMutVar out+        pure DoubleBuffer+          { vec = vec+          , count = count+          , sorted = sorted+          , growthIncrement = 0+          , spaceAtBottom = spaceAtBottom+          }+{-# SCC getEvensOrOdds #-}+++(!) :: PrimMonad m => DoubleBuffer (PrimState m) -> Int -> m Double+(!) buf offset = do+  index <- if spaceAtBottom buf+    then do+      capacity_ <- getCapacity buf+      count_ <- getCount buf+      pure $! capacity_ - count_ + offset+    else pure offset+  vec <- getVector buf+  MUVector.read vec index++getCount :: PrimMonad m => DoubleBuffer (PrimState m) -> m Int+getCount = readURef . count++getSpace :: PrimMonad m => DoubleBuffer (PrimState m) -> m Int+getSpace buf@DoubleBuffer{..} = (-) <$> getCapacity buf <*> getCount buf++isEmpty :: PrimMonad m => DoubleBuffer (PrimState m) -> m Bool+isEmpty buf = (== 0) <$> getCount buf++isSorted :: PrimMonad m => DoubleBuffer (PrimState m) -> m Bool+isSorted = readURef . sorted++-- | Sorts the active region+sort :: PrimMonad m => DoubleBuffer (PrimState m) -> m ()+sort buf@DoubleBuffer{..} = do+  sorted_ <- isSorted buf+  unless sorted_ $ do+    capacity_ <- getCapacity buf+    count_ <- getCount buf+    let (start, end) = if spaceAtBottom+          then (capacity_ - count_, capacity_)+          else (0, count_)+    vec <- getVector buf+    sortByBounds compare vec start end+    writeURef sorted True+{-# SCC sort #-}++-- | Merges the incoming sorted buffer into this sorted buffer.+mergeSortIn :: (PrimMonad m, HasCallStack) => DoubleBuffer (PrimState m) -> DoubleBuffer (PrimState m) -> m ()+mergeSortIn this bufIn = do+  sort this+  sort bufIn++  thatBuf <- getVector bufIn+  bufInLen <- getCount bufIn++  ensureSpace this bufInLen+  count_ <- getCount this+  let totalLength = count_ + bufInLen++  thisBuf <- getVector this++  if spaceAtBottom this+    then do -- scan up, insert at bottom+      capacity_ <- getCapacity this+      bufInCapacity_ <- getCapacity bufIn+      inSs <- takeSnapshot bufIn+      let i = capacity_ - count_+      let j = bufInCapacity_ - bufInLen+      let targetStart = capacity_ - totalLength+      let k = targetStart+      mergeUpwards thisBuf thatBuf capacity_ bufInCapacity_ i j k+    else do -- scan down, insert at top+      let i = count_ - 1+      let j = bufInLen - 1+      let k = totalLength+      mergeDownwards thisBuf thatBuf i j (k - 1)++  modifyURef (count this) (+ bufInLen)+  writeURef (sorted this) True+  pure ()+  where+    mergeUpwards thisBuf thatBuf capacity_ bufInCapacity_ = go+      where+        go !i !j !k+          -- for loop ended+          | k >= capacity_ = pure ()+          -- both valid+          | i < capacity_ && j < bufInCapacity_ = do+            iVal <- MUVector.read thisBuf i+            jVal <- MUVector.read thatBuf j+            if iVal <= jVal+              then MUVector.unsafeWrite thisBuf k iVal >> go (i + 1) j (k + 1)+              else MUVector.unsafeWrite thisBuf k jVal >> go i (j + 1) (k + 1)+          -- i is valid+          | i < capacity_ = do+            MUVector.unsafeWrite thisBuf k =<< MUVector.read thisBuf i+            go (i + 1) j (k + 1)+          -- j is valid+          | j < bufInCapacity_ = do+            MUVector.unsafeWrite thisBuf k =<< MUVector.read thatBuf j+            go i (j + 1) (k + 1)+          -- neither is valid, break;+          | otherwise = pure ()+    mergeDownwards thisBuf thatBuf !i !j !k+      -- for loop ended+      | k < 0 = pure ()+      -- both valid+      | i >= 0 && j >= 0 = do+        iVal <- MUVector.read thisBuf i+        jVal <- MUVector.read thatBuf j+        if iVal >= jVal+          then do+            MUVector.unsafeWrite thisBuf k iVal >> continue (i - 1) j (k - 1)+          else do+            MUVector.unsafeWrite thisBuf k jVal >> continue i (j - 1) (k - 1)+      | i >= 0 = do+        MUVector.unsafeWrite thisBuf k =<< MUVector.read thisBuf i+        continue (i - 1) j (k - 1)+      | j >= 0 = do+        MUVector.unsafeWrite thisBuf k =<< MUVector.read thatBuf j+        continue i (j - 1) (k - 1)+      -- neither is valid, break;+      | otherwise = pure ()+      where+        continue = mergeDownwards thisBuf thatBuf+{-# SCC mergeSortIn #-}++trimCount :: PrimMonad m => DoubleBuffer (PrimState m) -> Int -> m ()+trimCount DoubleBuffer{..} newCount = modifyURef count (\oldCount -> if newCount < oldCount then newCount else oldCount)
+ src/DataSketches/Quantiles/RelativeErrorQuantile/Internal/InequalitySearch.hs view
@@ -0,0 +1,121 @@+module DataSketches.Quantiles.RelativeErrorQuantile.Internal.InequalitySearch where++import Control.Monad.Primitive+import Data.Vector.Generic.Mutable (MVector)+import qualified Data.Vector.Generic.Mutable as MV++data (:<) = (:<)+data (:<=) = (:<=)+data (:>) = (:>)+data (:>=) = (:>=)+-- JavaDoc copypasta+-- +-- This provides efficient, unique and unambiguous binary searching for inequality comparison criteria+--  for ordered arrays of values that may include duplicate values. The inequality criteria include+--  <, >, ==, >=, <=. All the inequality criteria use the same search algorithm.+--  (Although == is not an inequality, it is included for convenience.)++--  In order to make the searching unique and unambiguous, we modified the traditional binary+--  search algorithm to search for adjacent pairs of values <i>{A, B}</i> in the values array+--  instead of just a single value, where <i>A</i> and <i>B</i> are the array indicies of two+--  adjacent values in the array. For all the search criteria, if the algorithm reaches the ends of+--  the search range, the algorithm calls the <i>resolve()</i> method to determine what to+--   return to the caller. If the key value cannot be resolved, it returns a -1 to the caller.++--  Given an array of values <i>arr[]</i> and the search key value <i>v</i>, the algorithms for+--  the searching criteria are as follows:</p>+-- +--  <li><b>LT:</b> Find the highest ranked adjacent pair <i>{A, B}</i> such that:<br>+--  <i>arr[A] &lt; v &le; arr[B]</i>. The normal return is the index <i>A</i>.+--  </li>+--  <li><b>LE:</b>  Find the highest ranked adjacent pair <i>{A, B}</i> such that:<br>+--  <i>arr[A] &le; v &lt; arr[B]</i>. The normal return is the index <i>A</i>.+--  </li>+--  <li><b>EQ:</b>  Find the adjacent pair <i>{A, B}</i> such that:<br>+--  <i>arr[A] &le; v &le; arr[B]</i>. The normal return is the index <i>A</i> or <i>B</i> whichever+--  equals <i>v</i>, otherwise it returns -1.+--  </li>+--  <li><b>GE:</b>  Find the lowest ranked adjacent pair <i>{A, B}</i> such that:<br>+--  <i>arr[A] &lt; v &le; arr[B]</i>. The normal return is the index <i>B</i>.+--  </li>+--  <li><b>GT:</b>  Find the lowest ranked adjacent pair <i>{A, B}</i> such that:<br>+--  <i>arr[A] &le; v &lt; arr[B]</i>. The normal return is the index <i>B</i>.+--  </li>+--  </ul>+class InequalitySearch s where+  inequalityCompare :: Ord a+    => s+    -> a +    -- ^ V+    -> a+    -- ^ A+    -> a+    -- ^ B+    -> Ordering+    -- ^ 'GT' means we must search higher in the array, 'LT' means we must+    -- search lower in the array, or `EQ`, which means we have found +    -- the correct bounding pair.+  getIndex +    :: (PrimMonad m, MVector v a, Ord a) +    => s +    -> v (PrimState m) a +    -> Int +    -> Int +    -> a +    -> m Int+  resolve +    :: s +    -> Int -- Vector length+    -> (Int, Int) +    -- ^ Final low index, high index (lo, hi)+    -> (Int, Int) +    -- ^ Initial search region (low, high)+    -> Int+    -- ^ A thing++instance InequalitySearch (:<) where+  inequalityCompare _ v a b +    | v <= a = LT+    | b < v = GT+    | otherwise = EQ+  getIndex _ _ a _ _ = pure a+  resolve _ vl (lo, hi) (low, high) = if lo >= high then high else vl++instance InequalitySearch (:<=) where+  inequalityCompare _ v a b+    | v < a = LT+    | b <= v = GT+    | otherwise = EQ+  getIndex _ _ a _ _ = pure a+  resolve _ vl (lo, hi) (low, high) = if lo >= high then high else vl ++instance InequalitySearch (:>) where+  inequalityCompare _ v a b+    | v < a = LT+    | b <= v = GT+    | otherwise = EQ+  getIndex _ _ _ b _ = pure b+  resolve _ vl (lo, hi) (low, high) = if hi <= low then low else vl++instance InequalitySearch (:>=) where+  inequalityCompare _ v a b+    | v <= a = LT+    | b < v = GT+    | otherwise = EQ+  getIndex _ _ _ b _ = pure b+  resolve _ vl (lo, hi) (low, high) = if hi <= low then low else vl++find :: (InequalitySearch s, PrimMonad m, MVector v a, Ord a) => s -> v (PrimState m) a -> Int -> Int -> a -> m Int+find strat v low high x = go low (high - 1)+  where+    go lo hi +      | lo <= hi && lo < high = do+          let mid = lo + ((hi - lo) `div` 2)+          midV <- MV.read v mid+          midV' <- MV.read v (mid + 1)+          case inequalityCompare strat x midV midV' of+            LT -> go lo (mid - 1)+            EQ -> getIndex strat v mid (mid + 1) x+            GT -> go (mid + 1) hi+      | otherwise = pure $! resolve strat (MV.length v) (lo, hi) (low, high)+{-# INLINE find #-}
+ src/DataSketches/Quantiles/RelativeErrorQuantile/Types.hs view
@@ -0,0 +1,25 @@+module DataSketches.Quantiles.RelativeErrorQuantile.Types where+import Control.Monad.Primitive+import qualified DataSketches.Quantiles.RelativeErrorQuantile.Internal.InequalitySearch as IS++data Criterion = (:<) | (:<=)+  deriving (Show, Eq)++instance IS.InequalitySearch Criterion where+  inequalityCompare c = case c of+    (:<) -> IS.inequalityCompare (IS.:<)+    (:<=) -> IS.inequalityCompare(IS.:<=)+  resolve c = case c of+    (:<) -> IS.resolve (IS.:<)+    (:<=) -> IS.resolve (IS.:<=)+  getIndex c = case c of+    (:<) -> IS.getIndex (IS.:<)+    (:<=) -> IS.getIndex (IS.:<=)+++data RankAccuracy +  = HighRanksAreAccurate +  -- ^ High ranks are prioritized for better accuracy.+  | LowRanksAreAccurate+  -- ^ Low ranks are prioritized for better accuracy+  deriving (Show, Eq)
+ test/Spec.hs view
@@ -0,0 +1,2 @@+main :: IO ()+main = putStrLn "Test suite not yet implemented"