diff --git a/random-fu/changelog.md b/random-fu/changelog.md
index a22da78..6776e24 100644
--- a/random-fu/changelog.md
+++ b/random-fu/changelog.md
@@ -1,3 +1,7 @@
+* Chnages in 0.3.0.0:
+
+  * Drop usage of `random-source` in favor of `random`
+
 * Changes in 0.2.7.7: Update to random-1.2. Revert 0.2.7.6 changes (which added an extra constraint to `Data.Random.Sample.sampleState` and `Data.Random.Sample.sampleStateT`).
 
 * Changes in 0.2.7.4: Compatibility with ghc 8.8.
diff --git a/random-fu/random-fu.cabal b/random-fu/random-fu.cabal
index e7e3db9..6d7593a 100644
--- a/random-fu/random-fu.cabal
+++ b/random-fu/random-fu.cabal
@@ -1,5 +1,5 @@
 name:                   random-fu
-version:                0.2.7.7
+version:                0.3.0.0
 stability:              provisional
 
 cabal-version:          >= 1.10
@@ -12,21 +12,21 @@ homepage:               https://github.com/mokus0/random-fu
 
 category:               Math
 synopsis:               Random number generation
-description:            Random number generation based on modeling random 
+description:            Random number generation based on modeling random
                         variables in two complementary ways: first, by the
                         parameters of standard mathematical distributions and,
                         second, by an abstract type ('RVar') which can be
                         composed and manipulated monadically and sampled in
                         either monadic or \"pure\" styles.
                         .
-                        The primary purpose of this library is to support 
+                        The primary purpose of this library is to support
                         defining and sampling a wide variety of high quality
                         random variables.  Quality is prioritized over speed,
                         but performance is an important goal too.
                         .
-                        In my testing, I have found it capable of speed 
+                        In my testing, I have found it capable of speed
                         comparable to other Haskell libraries, but still
-                        a fair bit slower than straight C implementations of 
+                        a fair bit slower than straight C implementations of
                         the same algorithms.
 
 tested-with:            GHC == 7.10.3
@@ -83,25 +83,24 @@ Library
   else
     cpp-options:        -Dold_Fixed
     build-depends:      base >= 4 && <4.2
-  
+
   if flag(mtl2)
     build-depends:      mtl == 2.*
     cpp-options:        -DMTL2
   else
     build-depends:      mtl == 1.*
-  
+
   build-depends:        math-functions,
                         monad-loops >= 0.3.0.1,
                         random >= 1.2 && < 1.3,
                         random-shuffle,
-                        random-source == 0.3.*,
-                        rvar == 0.2.*,
+                        rvar >= 0.3,
                         syb,
                         template-haskell,
                         transformers,
                         vector >= 0.7,
                         erf
-  
+
   if impl(ghc == 7.2.1)
     -- Doesn't work under GHC 7.2.1 due to
     -- http://hackage.haskell.org/trac/ghc/ticket/5410
diff --git a/random-fu/src/Data/Random.hs b/random-fu/src/Data/Random.hs
index 58bcc0d..1218132 100644
--- a/random-fu/src/Data/Random.hs
+++ b/random-fu/src/Data/Random.hs
@@ -1,39 +1,39 @@
 -- |Flexible modeling and sampling of random variables.
 --
--- The central abstraction in this library is the concept of a random 
--- variable.  It is not fully formalized in the standard measure-theoretic 
--- language, but rather is informally defined as a \"thing you can get random 
--- values out of\".  Different random variables may have different types of 
+-- The central abstraction in this library is the concept of a random
+-- variable.  It is not fully formalized in the standard measure-theoretic
+-- language, but rather is informally defined as a \"thing you can get random
+-- values out of\".  Different random variables may have different types of
 -- values they can return or the same types but different probabilities for
 -- each value they can return.  The random values you get out of them are
 -- traditionally called \"random variates\".
--- 
--- Most imperative-language random number libraries are all about obtaining 
--- and manipulating random variates.  This one is about defining, manipulating 
--- and sampling random variables.  Computationally, the distinction is small 
--- and mostly just a matter of perspective, but from a program design 
+--
+-- Most imperative-language random number libraries are all about obtaining
+-- and manipulating random variates.  This one is about defining, manipulating
+-- and sampling random variables.  Computationally, the distinction is small
+-- and mostly just a matter of perspective, but from a program design
 -- perspective it provides both a powerfully composable abstraction and a
 -- very useful separation of concerns.
--- 
+--
 -- Abstract random variables as implemented by 'RVar' are composable.  They can
 -- be defined in a monadic / \"imperative\" style that amounts to manipulating
 -- variates, but with strict type-level isolation.  Concrete random variables
 -- are also provided, but they do not compose as generically.  The 'Distribution'
--- type class allows concrete random variables to \"forget\" their concreteness 
--- so that they can be composed.  For examples of both, see the documentation 
--- for 'RVar' and 'Distribution', as well as the code for any of the concrete 
+-- type class allows concrete random variables to \"forget\" their concreteness
+-- so that they can be composed.  For examples of both, see the documentation
+-- for 'RVar' and 'Distribution', as well as the code for any of the concrete
 -- distributions such as 'Uniform', 'Gamma', etc.
--- 
+--
 -- Both abstract and concrete random variables can be sampled (despite the
 -- types GHCi may list for the functions) by the functions in "Data.Random.Sample".
--- 
+--
 -- Random variable sampling is done with regard to a generic basis of primitive
--- random variables defined in "Data.Random.Internal.Primitives".  This basis 
+-- random variables defined in "Data.Random.Internal.Primitives".  This basis
 -- is very low-level and the actual set of primitives is still fairly experimental,
 -- which is why it is in the \"Internal\" sub-heirarchy.  User-defined variables
 -- should use the existing high-level variables such as 'Uniform' and 'Normal'
 -- rather than these basis variables.  "Data.Random.Source" defines classes for
--- entropy sources that provide implementations of these primitive variables. 
+-- entropy sources that provide implementations of these primitive variables.
 -- Several implementations are available in the Data.Random.Source.* modules.
 module Data.Random
     ( -- * Random variables
@@ -43,32 +43,26 @@ module Data.Random
 
       -- ** Concrete ('Distribution')
       Distribution(..), CDF(..), PDF(..),
-      
+
       -- * Sampling random variables
-      Sampleable(..), sample, sampleState, sampleStateT,
-      
+      Sampleable(..), sample, sampleState, samplePure,
+
       -- * A few very common distributions
       Uniform(..), uniform, uniformT,
       StdUniform(..), stdUniform, stdUniformT,
       Normal(..), normal, stdNormal, normalT, stdNormalT,
       Gamma(..), gamma, gammaT,
-      
+
       -- * Entropy Sources
-      MonadRandom, RandomSource, StdRandom(..),
-      
+      StatefulGen, RandomGen,
+
       -- * Useful list-based operations
       randomElement,
       shuffle, shuffleN, shuffleNofM
-      
+
     ) where
 
 import Data.Random.Sample
-import Data.Random.Source (MonadRandom, RandomSource)
-import Data.Random.Source.IO ()
-import Data.Random.Source.MWC ()
-import Data.Random.Source.StdGen ()
-import Data.Random.Source.PureMT ()
-import Data.Random.Source.Std
 import Data.Random.Distribution
 import Data.Random.Distribution.Gamma
 import Data.Random.Distribution.Normal
@@ -78,3 +72,4 @@ import Data.Random.Lift ()
 import Data.Random.List
 import Data.Random.RVar
 
+import System.Random.Stateful (StatefulGen, RandomGen)
diff --git a/random-fu/src/Data/Random/Distribution/Categorical.hs b/random-fu/src/Data/Random/Distribution/Categorical.hs
index 575a0a2..9980a4c 100644
--- a/random-fu/src/Data/Random/Distribution/Categorical.hs
+++ b/random-fu/src/Data/Random/Distribution/Categorical.hs
@@ -23,9 +23,7 @@ import Data.Random.Distribution.Uniform
 import Control.Arrow
 import Control.Monad
 import Control.Monad.ST
-import Data.Foldable (Foldable(foldMap))
 import Data.STRef
-import Data.Traversable (Traversable(traverse, sequenceA))
 
 import Data.List
 import Data.Function
diff --git a/random-fu/src/Data/Random/Distribution/Normal.hs b/random-fu/src/Data/Random/Distribution/Normal.hs
index 90b9e11..9aa041b 100644
--- a/random-fu/src/Data/Random/Distribution/Normal.hs
+++ b/random-fu/src/Data/Random/Distribution/Normal.hs
@@ -22,14 +22,13 @@ module Data.Random.Distribution.Normal
     , knuthPolarNormalPair
     ) where
 
-import Data.Random.Internal.Words
 import Data.Bits
 
-import Data.Random.Source
 import Data.Random.Distribution
 import Data.Random.Distribution.Uniform
 import Data.Random.Distribution.Ziggurat
 import Data.Random.RVar
+import Data.Word
 
 import Data.Vector.Generic (Vector)
 import qualified Data.Vector as V
@@ -37,6 +36,8 @@ import qualified Data.Vector.Unboxed as UV
 
 import Data.Number.Erf
 
+import qualified System.Random.Stateful as Random
+
 -- |A random variable that produces a pair of independent
 -- normally-distributed values.
 normalPair :: (Floating a, Distribution StdUniform a) => RVar (a,a)
@@ -141,7 +142,7 @@ realFloatStdNormal = runZiggurat (normalZ p getIU `asTypeOf` (undefined :: Ziggu
 
         getIU :: (Num a, Distribution Uniform a) => RVarT m (Int, a)
         getIU = do
-            i <- getRandomWord8
+            i <- Random.uniformWord8 RGen
             u <- uniformT (-1) 1
             return (fromIntegral i .&. (2^p-1), u)
 
@@ -167,10 +168,24 @@ doubleStdNormalZ = mkZiggurat_ True
     where
         getIU :: RVarT m (Int, Double)
         getIU = do
-            !w <- getRandomWord64
+            !w <- Random.uniformWord64 RGen
             let (u,i) = wordToDoubleWithExcess w
             return $! (fromIntegral i .&. (doubleStdNormalC-1), u+u-1)
 
+-- NOTE: inlined from random-source
+{-# INLINE wordToDouble #-}
+-- |Pack the low 52 bits from a 'Word64' into a 'Double' in the range [0,1).
+-- Used to convert a 'stdUniform' 'Word64' to a 'stdUniform' 'Double'.
+wordToDouble :: Word64 -> Double
+wordToDouble x = (encodeFloat $! toInteger (x .&. 0x000fffffffffffff {- 2^52-1 -})) $ (-52)
+
+{-# INLINE wordToDoubleWithExcess #-}
+-- |Same as wordToDouble, but also return the unused bits (as the 12
+-- least significant bits of a 'Word64')
+wordToDoubleWithExcess :: Word64 -> (Double, Word64)
+wordToDoubleWithExcess x = (wordToDouble x, x `shiftR` 52)
+
+
 -- |A random variable sampling from the standard normal distribution
 -- over the 'Float' type.
 floatStdNormal :: RVarT m Float
@@ -193,10 +208,24 @@ floatStdNormalZ = mkZiggurat_ True
     where
         getIU :: RVarT m (Int, Float)
         getIU = do
-            !w <- getRandomWord32
+            !w <- Random.uniformWord32 RGen
             let (u,i) = word32ToFloatWithExcess w
             return (fromIntegral i .&. (floatStdNormalC-1), u+u-1)
 
+-- NOTE: inlined from random-source
+{-# INLINE word32ToFloat #-}
+-- |Pack the low 23 bits from a 'Word32' into a 'Float' in the range [0,1).
+-- Used to convert a 'stdUniform' 'Word32' to a 'stdUniform' 'Double'.
+word32ToFloat :: Word32 -> Float
+word32ToFloat x = (encodeFloat $! toInteger (x .&. 0x007fffff {- 2^23-1 -} )) $ (-23)
+
+{-# INLINE word32ToFloatWithExcess #-}
+-- |Same as word32ToFloat, but also return the unused bits (as the 9
+-- least significant bits of a 'Word32')
+word32ToFloatWithExcess :: Word32 -> (Float, Word32)
+word32ToFloatWithExcess x = (word32ToFloat x, x `shiftR` 23)
+
+
 normalCdf :: (Real a) => a -> a -> a -> Double
 normalCdf m s x = normcdf ((realToFrac x - realToFrac m) / realToFrac s)
 
diff --git a/random-fu/src/Data/Random/Distribution/Uniform.hs b/random-fu/src/Data/Random/Distribution/Uniform.hs
index a580fe8..8eb5669 100644
--- a/random-fu/src/Data/Random/Distribution/Uniform.hs
+++ b/random-fu/src/Data/Random/Distribution/Uniform.hs
@@ -38,10 +38,8 @@ module Data.Random.Distribution.Uniform
     ) where
 
 import Data.Random.Internal.TH
-import Data.Random.Internal.Words
 import Data.Random.Internal.Fixed
 
-import Data.Random.Source
 import Data.Random.Distribution
 import Data.Random.RVar
 
@@ -51,39 +49,14 @@ import Data.Int
 
 import Control.Monad.Loops
 
+import qualified System.Random.Stateful as Random
+
 -- |Compute a random 'Integral' value between the 2 values provided (inclusive).
 {-# INLINE integralUniform #-}
-integralUniform :: (Integral a) => a -> a -> RVarT m a
-integralUniform !x !y = if x < y then integralUniform' x y else integralUniform' y x
-
-{-# SPECIALIZE integralUniform' :: Int     -> Int     -> RVarT m Int   #-}
-{-# SPECIALIZE integralUniform' :: Int8    -> Int8    -> RVarT m Int8  #-}
-{-# SPECIALIZE integralUniform' :: Int16   -> Int16   -> RVarT m Int16 #-}
-{-# SPECIALIZE integralUniform' :: Int32   -> Int32   -> RVarT m Int32 #-}
-{-# SPECIALIZE integralUniform' :: Int64   -> Int64   -> RVarT m Int64 #-}
-{-# SPECIALIZE integralUniform' :: Word    -> Word    -> RVarT m Word   #-}
-{-# SPECIALIZE integralUniform' :: Word8   -> Word8   -> RVarT m Word8  #-}
-{-# SPECIALIZE integralUniform' :: Word16  -> Word16  -> RVarT m Word16 #-}
-{-# SPECIALIZE integralUniform' :: Word32  -> Word32  -> RVarT m Word32 #-}
-{-# SPECIALIZE integralUniform' :: Word64  -> Word64  -> RVarT m Word64 #-}
-{-# SPECIALIZE integralUniform' :: Integer -> Integer -> RVarT m Integer #-}
-integralUniform' :: (Integral a) => a -> a -> RVarT m a
-integralUniform' !l !u
-    | nReject == 0  = fmap shift prim
-    | otherwise     = fmap shift loop
-    where
-        m = 1 + toInteger u - toInteger l
-        (bytes, nPossible) = bytesNeeded m
-        nReject = nPossible `mod` m
-
-        !prim = getRandomNByteInteger bytes
-        !shift = \(!z) -> l + (fromInteger $! (z `mod` m))
-
-        loop = do
-            z <- prim
-            if z < nReject
-                then loop
-                else return z
+integralUniform :: Random.UniformRange a => a -> a -> RVarT m a
+integralUniform !x !y = Random.uniformRM (x, y) RGen
+  -- Maybe switch to uniformIntegralM (requires exposing from `random` internals):
+  -- Random.uniformIntegralM (x, y) RGen
 
 integralUniformCDF :: (Integral a, Fractional b) => a -> a -> a -> b
 integralUniformCDF a b x
@@ -92,15 +65,6 @@ integralUniformCDF a b x
     | x > b     = 1
     | otherwise = (fromIntegral x - fromIntegral a) / (fromIntegral b - fromIntegral a)
 
--- TODO: come up with a decent, fast heuristic to decide whether to return an extra
--- byte.  May involve moving calculation of nReject into this function, and then
--- accepting first if 4*nReject < nPossible or something similar.
-bytesNeeded :: Integer -> (Int, Integer)
-bytesNeeded x = head (dropWhile ((<= x).snd) powersOf256)
-
-powersOf256 :: [(Int, Integer)]
-powersOf256 = zip [0..] (iterate (256 *) 1)
-
 -- |Compute a random value for a 'Bounded' type, between 'minBound' and 'maxBound'
 -- (inclusive for 'Integral' or 'Enum' types, in ['minBound', 'maxBound') for Fractional types.)
 boundedStdUniform :: (Distribution Uniform a, Bounded a) => RVar a
@@ -120,13 +84,17 @@ boundedEnumStdUniformCDF = enumUniformCDF minBound maxBound
 -- |Compute a uniform random 'Float' value in the range [0,1)
 floatStdUniform :: RVarT m Float
 floatStdUniform = do
-    x <- getRandomWord32
-    return (word32ToFloat x)
+    x <- uniformRangeRVarT (0, 1)
+    -- exclude 1. TODO: come up with something smarter
+    if x == 1 then floatStdUniform else pure x
 
 -- |Compute a uniform random 'Double' value in the range [0,1)
 {-# INLINE doubleStdUniform #-}
 doubleStdUniform :: RVarT m Double
-doubleStdUniform = getRandomDouble
+doubleStdUniform = do
+    x <- uniformRangeRVarT (0, 1)
+    -- exclude 1. TODO: come up with something smarter
+    if x == 1 then doubleStdUniform else pure x
 
 -- |Compute a uniform random value in the range [0,1) for any 'RealFloat' type
 realFloatStdUniform :: RealFloat a => RVarT m a
@@ -289,27 +257,19 @@ $( replicateInstances ''Int integralTypes [d|
         instance CDF Uniform Int            where cdf   (Uniform a b) = integralUniformCDF a b
     |])
 
-instance Distribution StdUniform Word8      where rvarT _ = getRandomWord8
-instance Distribution StdUniform Word16     where rvarT _ = getRandomWord16
-instance Distribution StdUniform Word32     where rvarT _ = getRandomWord32
-instance Distribution StdUniform Word64     where rvarT _ = getRandomWord64
+instance Distribution StdUniform Word8      where rvarT _ = Random.uniformWord8 RGen
+instance Distribution StdUniform Word16     where rvarT _ = Random.uniformWord16 RGen
+instance Distribution StdUniform Word32     where rvarT _ = Random.uniformWord32 RGen
+instance Distribution StdUniform Word64     where rvarT _ = Random.uniformWord64 RGen
+instance Distribution StdUniform Word       where rvarT _ = uniformRVarT
 
-instance Distribution StdUniform Int8       where rvarT _ = fromIntegral `fmap` getRandomWord8
-instance Distribution StdUniform Int16      where rvarT _ = fromIntegral `fmap` getRandomWord16
-instance Distribution StdUniform Int32      where rvarT _ = fromIntegral `fmap` getRandomWord32
-instance Distribution StdUniform Int64      where rvarT _ = fromIntegral `fmap` getRandomWord64
+instance Distribution StdUniform Int8       where rvarT _ = uniformRVarT
+instance Distribution StdUniform Int16      where rvarT _ = uniformRVarT
+instance Distribution StdUniform Int32      where rvarT _ = uniformRVarT
+instance Distribution StdUniform Int64      where rvarT _ = uniformRVarT
 
-instance Distribution StdUniform Int where
-    rvar _ =
-        $(if toInteger (maxBound :: Int) > toInteger (maxBound :: Int32)
-            then [|fromIntegral `fmap` getRandomWord64 :: RVar Int|]
-            else [|fromIntegral `fmap` getRandomWord32 :: RVar Int|])
+instance Distribution StdUniform Int        where rvarT _ = uniformRVarT
 
-instance Distribution StdUniform Word where
-    rvar _ =
-        $(if toInteger (maxBound :: Word) > toInteger (maxBound :: Word32)
-            then [|fromIntegral `fmap` getRandomWord64 :: RVar Word|]
-            else [|fromIntegral `fmap` getRandomWord32 :: RVar Word|])
 
 -- Integer has no StdUniform...
 
@@ -331,7 +291,7 @@ instance CDF Uniform Float                  where cdf   (Uniform a b) = realUnif
 instance CDF Uniform Double                 where cdf   (Uniform a b) = realUniformCDF a b
 
 instance Distribution StdUniform Float      where rvarT _ = floatStdUniform
-instance Distribution StdUniform Double     where rvarT _ = getRandomDouble
+instance Distribution StdUniform Double     where rvarT _ = uniformRangeRVarT (0, 1)
 instance CDF StdUniform Float               where cdf   _ = realStdUniformCDF
 instance CDF StdUniform Double              where cdf   _ = realStdUniformCDF
 instance PDF StdUniform Float               where pdf   _ = realStdUniformPDF
@@ -357,7 +317,7 @@ $( replicateInstances ''Char [''Char, ''Bool, ''Ordering] [d|
 
 instance Distribution StdUniform ()         where rvarT ~StdUniform = return ()
 instance CDF StdUniform ()                  where cdf   ~StdUniform = return 1
-instance Distribution StdUniform Bool       where rvarT ~StdUniform = fmap even (getRandomWord8)
+instance Distribution StdUniform Bool       where rvarT ~StdUniform = uniformRVarT
 instance CDF StdUniform Bool                where cdf   ~StdUniform = boundedEnumStdUniformCDF
 
 instance Distribution StdUniform Char       where rvarT ~StdUniform = boundedEnumStdUniform
diff --git a/random-fu/src/Data/Random/Lift.hs b/random-fu/src/Data/Random/Lift.hs
index 658ac72..e62eb4b 100644
--- a/random-fu/src/Data/Random/Lift.hs
+++ b/random-fu/src/Data/Random/Lift.hs
@@ -5,7 +5,6 @@ module Data.Random.Lift where
 import Data.RVar
 import qualified Data.Functor.Identity as T
 import qualified Control.Monad.Trans.Class as T
-import Data.Random.Source.Std
 
 #ifndef MTL2
 import qualified Control.Monad.Identity as MTL
@@ -41,7 +40,7 @@ instance Monad m => Lift T.Identity m where
     lift = return . T.runIdentity
 
 instance Lift (RVarT T.Identity) (RVarT m) where
-    lift x = runRVar x StdRandom
+    lift x = runRVar x RGen
 
 -- | This instance is again incoherent with the others, but provides a
 -- more-specific instance to resolve the overlap between the
@@ -58,7 +57,7 @@ instance Monad m => Lift MTL.Identity m where
     lift = return . MTL.runIdentity
 
 instance Lift (RVarT MTL.Identity) (RVarT m) where
-    lift x = runRVarTWith (return . MTL.runIdentity) x StdRandom
+    lift x = runRVarTWith (return . MTL.runIdentity) x RGen
 
 -- | This instance is again incoherent with the others, but provides a
 -- more-specific instance to resolve the overlap between the
diff --git a/random-fu/src/Data/Random/RVar.hs b/random-fu/src/Data/Random/RVar.hs
index 1e69672..029ff49 100644
--- a/random-fu/src/Data/Random/RVar.hs
+++ b/random-fu/src/Data/Random/RVar.hs
@@ -2,13 +2,14 @@
 module Data.Random.RVar
     ( RVar, runRVar
     , RVarT, runRVarT, runRVarTWith
+    , RGen(..), uniformRVarT, uniformRangeRVarT
     ) where
 
 import Data.Random.Lift
-import Data.Random.Internal.Source
 import Data.RVar hiding (runRVarT)
+import System.Random.Stateful
 
 -- |Like 'runRVarTWith', but using an implicit lifting (provided by the
 -- 'Lift' class)
-runRVarT :: (Lift n m, RandomSource m s) => RVarT n a -> s -> m a
+runRVarT :: (Lift n m, StatefulGen g m) => RVarT n a -> g -> m a
 runRVarT = runRVarTWith lift
diff --git a/random-fu/src/Data/Random/Sample.hs b/random-fu/src/Data/Random/Sample.hs
index c1c730c..b7c0f93 100644
--- a/random-fu/src/Data/Random/Sample.hs
+++ b/random-fu/src/Data/Random/Sample.hs
@@ -9,11 +9,12 @@
 module Data.Random.Sample where
 
 import Control.Monad.State
+import Control.Monad.Reader
 import Data.Random.Distribution
 import Data.Random.Lift
 import Data.Random.RVar
-import Data.Random.Source
-import Data.Random.Source.Std
+
+import System.Random.Stateful
 
 -- |A typeclass allowing 'Distribution's and 'RVar's to be sampled.  Both may
 -- also be sampled via 'runRVar' or 'runRVarT', but I find it psychologically
@@ -21,26 +22,29 @@ import Data.Random.Source.Std
 -- separate abstractions for one base concept: a random variable.
 class Sampleable d m t where
     -- |Directly sample from a distribution or random variable, using the given source of entropy.
-    sampleFrom :: RandomSource m s => s -> d t -> m t
+    sampleFrom :: StatefulGen g m => g -> d t -> m t
 
 instance Distribution d t => Sampleable d m t where
-    sampleFrom src d = runRVarT (rvar d) src
+    sampleFrom gen d = runRVarT (rvar d) gen
 
 -- This instance overlaps with the other, but because RVarT is not a Distribution there is no conflict.
 instance Lift m n => Sampleable (RVarT m) n t where
-    sampleFrom src x = runRVarT x src
+    sampleFrom gen x = runRVarT x gen
 
 -- |Sample a random variable using the default source of entropy for the
 -- monad in which the sampling occurs.
-sample :: (Sampleable d m t, MonadRandom m) => d t -> m t
-sample = sampleFrom StdRandom
+sample :: (Sampleable d m t, StatefulGen g m, MonadReader g m) => d t -> m t
+sample thing = ask >>= \gen -> sampleFrom gen thing
 
 -- |Sample a random variable in a \"functional\" style.  Typical instantiations
 -- of @s@ are @System.Random.StdGen@ or @System.Random.Mersenne.Pure64.PureMT@.
-sampleState :: (Sampleable d (State s) t, MonadRandom (State s)) => d t -> s -> (t, s)
-sampleState thing = runState (sample thing)
+-- sample :: (Distribution d a, StatefulGen g m, MonadReader g m) => d t -> m t
+-- sample thing gen = runStateGen gen (\stateGen -> sampleFrom stateGen thing)
 
--- |Sample a random variable in a \"semi-functional\" style.  Typical instantiations
--- of @s@ are @System.Random.StdGen@ or @System.Random.Mersenne.Pure64.PureMT@.
-sampleStateT :: (Sampleable d (StateT s m) t, MonadRandom (StateT s m)) => d t -> s -> m (t, s)
-sampleStateT thing = runStateT (sample thing)
+sampleState :: (Distribution d t, RandomGen g, MonadState g m) => d t -> m t
+sampleState thing = sampleFrom StateGenM thing
+
+-- |Sample a random variable in a \"functional\" style.  Typical instantiations
+-- of @g@ are @System.Random.StdGen@ or @System.Random.Mersenne.Pure64.PureMT@.
+samplePure :: (Distribution d t, RandomGen g) => d t -> g -> (t, g)
+samplePure thing gen = runStateGen gen (\stateGen -> sampleFrom stateGen thing)
diff --git a/random-source/src/Data/Random/Source/Internal/Prim.hs b/random-source/src/Data/Random/Source/Internal/Prim.hs
index e4bc67e..8d060d4 100644
--- a/random-source/src/Data/Random/Source/Internal/Prim.hs
+++ b/random-source/src/Data/Random/Source/Internal/Prim.hs
@@ -9,7 +9,7 @@ import Data.Typeable
 
 -- |A 'Prompt' GADT describing a request for a primitive random variate.
 -- Random variable definitions will request their entropy via these prompts,
--- and entropy sources will satisfy those requests.  The functions in
+-- and entropy sources will satisfy those requests.  The 'Integer' 0 <= U < 256^n'Integer' 0 <= U < 256^n'Integer' 0 <= U < 256^nfunctions in
 -- "Data.Random.Source.Internal.TH" extend incomplete entropy-source definitions
 -- to complete ones, essentially defining a very flexible 
 -- implementation-defaulting system.
diff --git a/rvar/changelog.md b/rvar/changelog.md
index 8d02357..e9b65c2 100644
--- a/rvar/changelog.md
+++ b/rvar/changelog.md
@@ -1,3 +1,9 @@
+* Changes in 0.3.0.0:
+
+  * Drop usage of `random-source` in favor of `random`
+  * Add `Prim` type that resembles one from `random-source`
+  * Add `RGen` type that serves the same purpose as `StdRandom` in `random-source`
+
 * Changes in 0.2.0.6: None. (Pacify Hackage.)
 
 * Changes in 0.2.0.4: Update for GHC 8.8.
diff --git a/rvar/rvar.cabal b/rvar/rvar.cabal
index 34d110b..67d79fb 100644
--- a/rvar/rvar.cabal
+++ b/rvar/rvar.cabal
@@ -1,5 +1,5 @@
 name:                   rvar
-version:                0.2.0.6
+version:                0.3.0.0
 stability:              stable
 
 cabal-version:          >= 1.10
@@ -45,6 +45,7 @@ Library
   hs-source-dirs:       src
   default-language:     Haskell2010
   exposed-modules:      Data.RVar
+  other-modules:        Data.RVar.Prim
 
   if flag(mtl2)
     build-depends:      mtl == 2.*
@@ -53,6 +54,7 @@ Library
     build-depends:      mtl == 1.1.*
 
   build-depends:        base            >= 3 && <5,
+                        bytestring,
                         MonadPrompt     == 1.0.*,
-                        random-source   == 0.3.*,
-                        transformers    >= 0.2 && < 0.6
+                        transformers    >= 0.2 && < 0.6,
+                        random          >= 1.2.0
\ No newline at end of file
diff --git a/rvar/src/Data/RVar.hs b/rvar/src/Data/RVar.hs
index cdcba6c..af5ee15 100644
--- a/rvar/src/Data/RVar.hs
+++ b/rvar/src/Data/RVar.hs
@@ -1,49 +1,45 @@
 {-
  -      ``Data/Random/RVar''
  -}
-{-# LANGUAGE
-    RankNTypes,
-    MultiParamTypeClasses,
-    FlexibleInstances,
-    GADTs,
-    ScopedTypeVariables,
-    CPP
-  #-}
+{-# LANGUAGE CPP #-}
+{-# LANGUAGE FlexibleInstances #-}
+{-# LANGUAGE GADTs #-}
+{-# LANGUAGE MultiParamTypeClasses #-}
+{-# LANGUAGE RankNTypes #-}
+{-# LANGUAGE ScopedTypeVariables #-}
 
 -- |Random variables.  An 'RVar' is a sampleable random variable.  Because
 -- probability distributions form a monad, they are quite easy to work with
 -- in the standard Haskell monadic styles.  For examples, see the source for
 -- any of the 'Distribution' instances - they all are defined in terms of
 -- 'RVar's.
+{-# LANGUAGE FlexibleContexts #-}
+
 module Data.RVar
-    ( RandomSource
-    , MonadRandom
-        ( getRandomWord8
-        , getRandomWord16
-        , getRandomWord32
-        , getRandomWord64
-        , getRandomDouble
-        , getRandomNByteInteger
-        )
-
-    , RVar
-    , runRVar, sampleRVar
+    ( RVar
+    , runRVar, sampleReaderRVar, sampleStateRVar
+    , pureRVar
 
     , RVarT
-    , runRVarT, sampleRVarT
-    , runRVarTWith, sampleRVarTWith
-    ) where
+    , runRVarT, sampleReaderRVarT, sampleStateRVarT
+    , runRVarTWith, sampleReaderRVarTWith, sampleStateRVarTWith
 
+    , RGen(..)
+    , uniformRVarT
+    , uniformRangeRVarT
+
+    , Prim(..)
+    ) where
 
-import Data.Random.Internal.Source (Prim(..), MonadRandom(..), RandomSource(..))
-import Data.Random.Source ({-instances-})
 
-import qualified Control.Monad.Trans.Class as T
-import Control.Monad (liftM, ap)
-import Control.Monad.Prompt (MonadPrompt(..), PromptT, runPromptT)
 import qualified Control.Monad.IO.Class as T
-import qualified Control.Monad.Trans as MTL
+import Control.Monad.Prompt (MonadPrompt(..), PromptT, runPromptT)
+import Control.Monad.Reader as MTL
+import Control.Monad.State as MTL
+import qualified Control.Monad.Trans.Class as T
 import qualified Data.Functor.Identity as T
+import Data.RVar.Prim
+import System.Random.Stateful
 
 -- |An opaque type modeling a \"random variable\" - a value
 -- which depends on the outcome of some random event.  'RVar's
@@ -80,14 +76,21 @@ import qualified Data.Functor.Identity as T
 -- (where @sampleState = runState . sampleRVar@)
 type RVar = RVarT T.Identity
 
+-- | Sample random variable using `RandomGen` generator as source of entropy
+pureRVar :: RandomGen g => RVar a -> g -> (a, g)
+pureRVar rvar g = runStateGen g (runRVar rvar)
+
 -- |\"Run\" an 'RVar' - samples the random variable from the provided
 -- source of entropy.
-runRVar :: RandomSource m s => RVar a -> s -> m a
+runRVar :: StatefulGen g m => RVar a -> g -> m a
 runRVar = runRVarTWith (return . T.runIdentity)
 
 -- |@sampleRVar x@ is equivalent to @runRVar x 'StdRandom'@.
-sampleRVar :: MonadRandom m => RVar a -> m a
-sampleRVar = sampleRVarTWith (return . T.runIdentity)
+sampleReaderRVar :: (StatefulGen g m, MonadReader g m) => RVar a -> m a
+sampleReaderRVar = sampleReaderRVarTWith (return . T.runIdentity)
+
+sampleStateRVar :: (RandomGen g, MonadState g m) => RVar a -> m a
+sampleStateRVar = sampleStateRVarTWith (return . T.runIdentity)
 
 -- |A random variable with access to operations in an underlying monad.  Useful
 -- examples include any form of state for implementing random processes with hysteresis,
@@ -139,11 +142,15 @@ sampleRVar = sampleRVarTWith (return . T.runIdentity)
 -- >                 replicateM count rwalkState
 newtype RVarT m a = RVarT { unRVarT :: PromptT Prim m a }
 
-runRVarT :: RandomSource m s => RVarT m a -> s -> m a
+runRVarT :: StatefulGen g m => RVarT m a -> g -> m a
 runRVarT = runRVarTWith id
 
-sampleRVarT :: MonadRandom m => RVarT m a -> m a
-sampleRVarT = sampleRVarTWith id
+
+sampleStateRVarT :: (RandomGen g, MonadState g m) => RVarT m a -> m a
+sampleStateRVarT rvar = runRVarT rvar StateGenM
+
+sampleReaderRVarT :: (StatefulGen g m, MonadReader g m) => RVarT m a -> m a
+sampleReaderRVarT rvar = ask >>= runRVarT rvar
 
 -- | \"Runs\" an 'RVarT', sampling the random variable it defines.
 --
@@ -171,21 +178,53 @@ sampleRVarT = sampleRVarTWith id
 -- or functions manipulating 'RVar's would have to use higher-ranked
 -- types to enforce the same kind of isolation and polymorphism.
 {-# INLINE runRVarTWith #-}
-runRVarTWith :: forall m n s a. RandomSource m s => (forall t. n t -> m t) -> RVarT n a -> s -> m a
-runRVarTWith liftN (RVarT m) src = runPromptT return bindP bindN m
+runRVarTWith :: forall m n g a. StatefulGen g m => (forall t. n t -> m t) -> RVarT n a -> g -> m a
+runRVarTWith liftN (RVarT m) gen = runPromptT return bindP bindN m
+    where
+        bindP :: forall t. (Prim t -> (t -> m a) -> m a)
+        bindP prim cont = uniformPrimM prim gen >>= cont
+
+        bindN :: forall t. n t -> (t -> m a) -> m a
+        bindN nExp cont = liftN nExp >>= cont
+
+{-# INLINE uniformPrimM #-}
+uniformPrimM :: StatefulGen g m => Prim t -> g -> m t
+uniformPrimM prim g =
+    case prim of
+        PrimWord8             -> uniformWord8 g
+        PrimWord16            -> uniformWord16 g
+        PrimWord32            -> uniformWord32 g
+        PrimWord64            -> uniformWord64 g
+        PrimShortByteString n -> uniformShortByteString n g
+
+
+-- |@sampleRVarTWith lift x@ is equivalent to @runRVarTWith lift x 'StdRandom'@.
+{-# INLINE sampleReaderRVarTWith #-}
+sampleReaderRVarTWith ::
+       forall m n a g. (StatefulGen g m, MonadReader g m)
+    => (forall t. n t -> m t)
+    -> RVarT n a
+    -> m a
+sampleReaderRVarTWith liftN (RVarT m) = runPromptT return bindP bindN m
     where
         bindP :: forall t. (Prim t -> (t -> m a) -> m a)
-        bindP prim cont = getRandomPrimFrom src prim >>= cont
+        bindP prim cont = ask >>= uniformPrimM prim >>= cont
 
         bindN :: forall t. n t -> (t -> m a) -> m a
         bindN nExp cont = liftN nExp >>= cont
 
+
 -- |@sampleRVarTWith lift x@ is equivalent to @runRVarTWith lift x 'StdRandom'@.
-sampleRVarTWith :: forall m n a. MonadRandom m => (forall t. n t -> m t) -> RVarT n a -> m a
-sampleRVarTWith liftN (RVarT m) = runPromptT return bindP bindN m
+{-# INLINE sampleStateRVarTWith #-}
+sampleStateRVarTWith ::
+       forall m n a g. (RandomGen g, MonadState g m)
+    => (forall t. n t -> m t)
+    -> RVarT n a
+    -> m a
+sampleStateRVarTWith liftN (RVarT m) = runPromptT return bindP bindN m
     where
         bindP :: forall t. (Prim t -> (t -> m a) -> m a)
-        bindP prim cont = getRandomPrim prim >>= cont
+        bindP prim cont = uniformPrimM prim StateGenM >>= cont
 
         bindN :: forall t. n t -> (t -> m a) -> m a
         bindN nExp cont = liftN nExp >>= cont
@@ -197,9 +236,6 @@ instance Monad (RVarT n) where
     return x = RVarT (return $! x)
     (RVarT m) >>= k = RVarT (m >>= \x -> x `seq` unRVarT (k x))
 
-instance MonadRandom (RVarT n) where
-    getRandomPrim = RVarT . prompt
-
 instance Applicative (RVarT n) where
     pure  = return
     (<*>) = ap
@@ -222,3 +258,26 @@ instance MTL.MonadIO m => MTL.MonadIO (RVarT m) where
     liftIO = MTL.lift . MTL.liftIO
 
 #endif
+
+data RGen = RGen
+
+instance StatefulGen RGen (RVarT m) where
+    uniformWord8 RGen = RVarT $ prompt PrimWord8
+    {-# INLINE uniformWord8 #-}
+    uniformWord16 RGen = RVarT $ prompt PrimWord16
+    {-# INLINE uniformWord16 #-}
+    uniformWord32 RGen = RVarT $ prompt PrimWord32
+    {-# INLINE uniformWord32 #-}
+    uniformWord64 RGen = RVarT $ prompt PrimWord64
+    {-# INLINE uniformWord64 #-}
+    uniformShortByteString n RGen = RVarT $ prompt (PrimShortByteString n)
+    {-# INLINE uniformShortByteString #-}
+
+
+uniformRVarT :: Uniform a => RVarT m a
+uniformRVarT = uniformM RGen
+{-# INLINE uniformRVarT #-}
+
+uniformRangeRVarT :: UniformRange a => (a, a) -> RVarT m a
+uniformRangeRVarT r = uniformRM r RGen
+{-# INLINE uniformRangeRVarT #-}
diff --git a/rvar/src/Data/RVar/Prim.hs b/rvar/src/Data/RVar/Prim.hs
new file mode 100644
index 0000000..7399bb9
--- /dev/null
+++ b/rvar/src/Data/RVar/Prim.hs
@@ -0,0 +1,37 @@
+{-# LANGUAGE DeriveDataTypeable #-}
+{-# LANGUAGE GADTs #-}
+{-# LANGUAGE RankNTypes #-}
+-- |This is an internal interface to support the 'RVar' abstraction.  It
+-- reifies the operations provided by `System.Random.Stateful.StatefulGen` in a
+-- uniform and efficient way, as functions of type @Prim a -> m a@.
+module Data.RVar.Prim (Prim(..)) where
+
+import Data.Typeable
+import Data.Word
+import Data.ByteString.Short
+
+-- |A 'Prompt' GADT describing a request for a primitive random variate.  Random variable
+-- definitions will request their entropy via these prompts, and entropy sources will
+-- satisfy those requests. This data type is needed for creating
+-- `System.Random.Stateful.StatefulGen` instance for `Data.RVar.RVarT`
+--
+data Prim a where
+    -- | An unsigned byte, uniformly distributed from 0 to 0xff
+    PrimWord8           :: Prim Word8
+    -- | An unsigned 16-bit word, uniformly distributed from 0 to 0xffff
+    PrimWord16          :: Prim Word16
+    -- | An unsigned 32-bit word, uniformly distributed from 0 to 0xffffffff
+    PrimWord32          :: Prim Word32
+    -- | An unsigned 64-bit word, uniformly distributed from 0 to 0xffffffffffffffff
+    PrimWord64          :: Prim Word64
+    -- | A uniformly distributed `ShortByteString` of length @n@ bytes
+    PrimShortByteString :: !Int -> Prim ShortByteString
+    deriving (Typeable)
+
+instance Show (Prim a) where
+    showsPrec _p PrimWord8         = showString "PrimWord8"
+    showsPrec _p PrimWord16        = showString "PrimWord16"
+    showsPrec _p PrimWord32        = showString "PrimWord32"
+    showsPrec _p PrimWord64        = showString "PrimWord64"
+    showsPrec  p (PrimShortByteString n) =
+      showParen (p > 10) (showString "PrimShortByteString " . showsPrec 11 n)
diff --git a/stack.yaml b/stack.yaml
index d19a64c..8b446ce 100644
--- a/stack.yaml
+++ b/stack.yaml
@@ -35,7 +35,7 @@ resolver: lts-16.31
 #  - wai
 packages:
 - random-fu
-- random-source
+# - random-source
 - rvar
 - tests/speed
 # Dependency packages to be pulled from upstream that are not in the resolver
@@ -46,7 +46,6 @@ extra-deps:
 - flexible-defaults-0.0.3@sha256:6a7ab000561e1075003cb1053dfbbb4020ae2b02916776d1479c9c3fc82f5d0d,2508
 - splitmix-0.1.0.3@sha256:fc3aae74c467f4b608050bef53aec17904a618731df9407e655d8f3bf8c32d5c,6049
 - mwc-random-0.15.0.1@sha256:841c86f132c45cb81116e1d3a8a150cecc27820c2b4e38f8cf86e3fe7735c2ab,3370
-
 # Override default flag values for local packages and extra-deps
 # flags: {}
 
diff --git a/tests/speed/Bench.hs b/tests/speed/Bench.hs
index 1ec329e..2465278 100644
--- a/tests/speed/Bench.hs
+++ b/tests/speed/Bench.hs
@@ -1,4 +1,5 @@
-{-# LANGUAGE FlexibleContexts, RankNTypes #-}
+{-# LANGUAGE FlexibleContexts #-}
+{-# LANGUAGE RankNTypes #-}
 
 {-# OPTIONS_GHC -fno-warn-missing-signatures #-}
 {-# OPTIONS_GHC -fno-warn-type-defaults      #-}
@@ -6,7 +7,6 @@
 module Main where
 
 import Data.Random
-import Data.Random.Source.DevRandom
 import Data.Random.Distribution.Beta
 import Data.Random.Distribution.Binomial
 import Data.Random.Distribution.Dirichlet
@@ -16,35 +16,32 @@ import Data.Random.Distribution.Poisson
 import Data.Random.Distribution.Rayleigh
 import Data.Random.Distribution.Triangular
 
-import Control.Applicative ((<$>))
+import System.Random.Stateful (runStateGen_)
+
 import Control.DeepSeq (NFData)
 import Control.Monad
+import qualified Control.Monad.Random as CMR
 import Control.Monad.ST
 import Control.Monad.State
-import qualified Control.Monad.Random as CMR
-import Control.Monad.Trans (lift)
 import Foreign
 import System.Random hiding (uniform)
 import qualified System.Random.MWC as MWC
 import TestSupport
-import qualified Data.Vector.Unboxed as Vec
 
 import Criterion.Main
 
 main = do
     let count = 64000
 
-    mwcSrc <- newGenIO
-    mtSrc  <- newMTSrc
-    stdSrc <- newStdSrc
+    mwcGen <- newMWCGenIO
+    mtGen  <- newMTGenM
+    stdGen <- newStdGenM
 
     defaultMain
         [ bgroup "dists"
-            [ bgroup "MWC"          (dists mwcSrc       count)
-            , bgroup "PureMT"       (dists mtSrc        count)
-            , bgroup "StdGen"       (dists stdSrc       count)
-            , bgroup "DevRandom"    (dists DevRandom    count)
-            , bgroup "DevURandom"   (dists DevURandom   count)
+            [ bgroup "MWC"          (dists mwcGen       count)
+            , bgroup "PureMT"       (dists mtGen        count)
+            , bgroup "StdGen"       (dists stdGen       count)
             ]
 
         , bgroup "IO StdGen"
@@ -52,39 +49,39 @@ main = do
                 sum' <$> replicateM count (randomRIO (10,50)) :: IO Double
 
             , bench "uniform A" $ nfIO $ do
-                sum' <$> replicateM count (sampleFrom stdSrc (uniform 10 50)) :: IO Double
+                sum' <$> replicateM count (sampleFrom stdGen (uniform 10 50)) :: IO Double
 
             , bench "uniform B" $ nfIO $ do
-                sum' <$> sampleFrom stdSrc (replicateM count (uniform 10 50)) :: IO Double
+                sum' <$> sampleFrom stdGen (replicateM count (uniform 10 50)) :: IO Double
             ]
 
         , bgroup "pure StdGen"
             [ bgroup "Double"
                 [ bench "getRandomRs" $ nfIO $ do
-                    src <- newStdGen
-                    let xs = CMR.evalRand (CMR.getRandomRs (10,50)) src
+                    gen <- newStdGen
+                    let xs = CMR.evalRand (CMR.getRandomRs (10,50)) gen
                     return $ sum' (take count xs :: [Double])
                 , bench "RVarT, State - sample replicateM" $ nfIO $ do
-                    src <- newStdGen
-                    let xs = evalState (runRVar (replicateM count (uniform 10 50)) StdRandom) src
+                    gen <- newStdGen
+                    let xs = runStateGen_ gen $ runRVar (replicateM count (uniform 10 50))
                     return (sum' xs :: Double)
                 , bench "RVarT, State - replicateM sample" $ nfIO $ do
-                    src <- newStdGen
-                    let xs = evalState (replicateM count (runRVar (uniform 10 50) StdRandom)) src
+                    gen <- newStdGen
+                    let xs = runStateGen_ gen $ replicateM count . runRVar (uniform 10 50)
                     return (sum' xs :: Double)
                 ]
             , bgroup "Int"
                 [ bench "getRandomRs" $ nfIO $ do
-                    src <- newStdGen
-                    let xs = CMR.evalRand (CMR.getRandomRs (10,50)) src
+                    gen <- newStdGen
+                    let xs = CMR.evalRand (CMR.getRandomRs (10,50)) gen
                     return $ sum' (take count xs :: [Int])
                 , bench "RVarT, State - sample replicateM" $ nfIO $ do
-                    src <- newStdGen
-                    let xs = evalState (runRVar (replicateM count (uniform 10 50)) StdRandom) src
+                    gen <- newStdGen
+                    let xs = runStateGen_ gen $ runRVar (replicateM count (uniform 10 50))
                     return $ (sum' xs :: Int)
                 , bench "RVarT, State - replicateM sample" $ nfIO $ do
-                    src <- newStdGen
-                    let xs = evalState (replicateM count (runRVar (uniform 10 50) StdRandom)) src
+                    gen <- newStdGen
+                    let xs = runStateGen_ gen $ replicateM count . runRVar (uniform 10 50)
                     return $ (sum' xs :: Int)
                 ]
             ]
@@ -92,31 +89,32 @@ main = do
         , bgroup "MWC"
             [ bgroup "stdUniform"
                 [ bench "Double" $ nfIO $ do
-                    src <- newGenIO
-                    xs <- stToIO $ replicateM count (MWC.uniform src)
+                    gen <- newMWCGenIO
+                    xs <- stToIO $ replicateM count (MWC.uniform gen)
                     return $ sum' (xs :: [Double])
                 , bench "Int" $ nfIO $ do
-                    src <- newGenIO
-                    xs <- stToIO $ replicateM count (MWC.uniform src)
+                    gen <- newMWCGenIO
+                    xs <- stToIO $ replicateM count (MWC.uniform gen)
                     return $ sum' (xs :: [Int])
                 ]
             , bgroup "uniform"
                 [ bench "Double" $ nfIO $ do
-                    src <- newGenIO
-                    us <- stToIO $ replicateM count (MWC.uniform src)
+                    gen <- newMWCGenIO
+                    us <- stToIO $ replicateM count (MWC.uniform gen)
                     let xs = [(u - 0.5) * 20| u <- us]
                     return $ sum' (xs :: [Double])
                 ]
 --            , bgroup "normal"
 --                [ bench "Double" $ nfIO $ do
---                    src <- newGenIO
---                    xs <- stToIO $ replicateM count (MWC.normal src)
+--                    gen <- newGenIO
+--                    xs <- stToIO $ replicateM count (MWC.normal gen)
 --                    return $ sum' (xs :: [Double])
 --                ]
+            -- FIXME: uniformVector no longer works on Double
             -- , bgroup "uniformVector"
             --     [ bench "Double" $ nfIO $ do
-            --         src <- newGenIO
-            --         xs <- stToIO $ MWC.uniformVector src count
+            --         gen <- newGenIO
+            --         xs <- stToIO $ MWC.uniformVector gen count
             --         -- unboxed, so don't need to force it, but we sum it
             --         -- anyway to make the comparison fair between other tests
             --         return $ (Vec.sum xs :: Double)
@@ -127,71 +125,71 @@ main = do
 
         ]
 
-dists src count =
-    [ multiSuite src count "uniform"     (Uniform 10 50)
-    , multiSuite src count "stdUniform"  StdUniform
-    , multiSuite src count "poisson"     (Poisson 3 :: Num t => Poisson Double t)
-    , multiSuite src count "binomial 10" (Binomial 10 (0.5 :: Float))
-    , doubleSuite src count "stdNormal"   StdNormal
-    , doubleSuite src count "exponential" (Exp  2)
-    , doubleSuite src count "beta"        (Beta  2 5)
-    , doubleSuite src count "gamma"       (Gamma 2 5)
-    , doubleSuite src count "triangular"  (Triangular 2 5 14)
-    , doubleSuite src count "rayleigh"    (Rayleigh 1.6)
+dists gen count =
+    [ multiSuite gen count "uniform"     (Uniform 10 50)
+    , multiSuite gen count "stdUniform"  StdUniform
+    , multiSuite gen count "poisson"     (Poisson 3 :: Num t => Poisson Double t)
+    , multiSuite gen count "binomial 10" (Binomial 10 (0.5 :: Float))
+    , doubleSuite gen count "stdNormal"   StdNormal
+    , doubleSuite gen count "exponential" (Exp  2)
+    , doubleSuite gen count "beta"        (Beta  2 5)
+    , doubleSuite gen count "gamma"       (Gamma 2 5)
+    , doubleSuite gen count "triangular"  (Triangular 2 5 14)
+    , doubleSuite gen count "rayleigh"    (Rayleigh 1.6)
 
     , bench "dirichlet" $ nfIO $ do
-        sum' <$> sampleFrom src (dirichlet [1..fromIntegral count :: Double])
+        sum' <$> sampleFrom gen (dirichlet [1..fromIntegral count :: Double])
 
     , bgroup "multinomial"
         [ bgroup "many p"
             [ bench desc $ nfIO $ do
-                sum' <$> sampleFrom src (multinomial [1..1e4 :: Double] (n :: Int))
+                sum' <$> sampleFrom gen (multinomial [1..1e4 :: Double] (n :: Int))
             | (desc, n) <- [("small n", 10), ("medium n", 10^4), ("large n", 10^8)]
             ]
         , bgroup "few p"
             [bench desc $ nfIO $ do
                 replicateM_ 1000 $ do
-                    sum' <$> sampleFrom src (multinomial [1..10 :: Double] (n :: Int))
+                    sum' <$> sampleFrom gen (multinomial [1..10 :: Double] (n :: Int))
             | (desc, n) <- [("small n", 10), ("medium n", 10^4), ("large n", 10^8)]
             ]
         ]
 
     , bench "shuffle" $ nfIO $ do
-        sum' <$> sampleFrom src (shuffle [1..count])
+        sum' <$> sampleFrom gen (shuffle [1..count])
     ]
 
-multiSuite :: (Distribution d Double, Distribution d Int, RandomSource IO s) => s -> Int -> String -> (forall t. Num t => d t) -> Benchmark
-multiSuite src count name dist = bgroup name
-    [ doubleSuite src count "Double" dist
-    , intSuite    src count "Int"    dist
+multiSuite :: (Distribution d Double, Distribution d Int, StatefulGen g IO) => g -> Int -> String -> (forall t. Num t => d t) -> Benchmark
+multiSuite gen count name dist = bgroup name
+    [ doubleSuite gen count "Double" dist
+    , intSuite    gen count "Int"    dist
     ]
 
-doubleSuite :: (Distribution d Double, RandomSource IO s) => s -> Int -> String -> d Double -> Benchmark
+doubleSuite :: (Distribution d Double, StatefulGen g IO) => g -> Int -> String -> d Double -> Benchmark
 doubleSuite = suite
 
-intSuite :: (Distribution d Int, RandomSource IO s) => s -> Int -> String -> d Int -> Benchmark
+intSuite :: (Distribution d Int, StatefulGen g IO) => g -> Int -> String -> d Int -> Benchmark
 intSuite = suite
 
-suite :: (Storable t, Num t, Distribution d t, NFData t, RandomSource IO s) => s -> Int -> String -> d t -> Benchmark
-suite src count name var = bgroup name
+suite :: (Storable t, Num t, Distribution d t, NFData t, StatefulGen g IO) => g -> Int -> String -> d t -> Benchmark
+suite gen count name var = bgroup name
     [ bench "single sample" $ nfIO $ do
-        sampleFrom src var
+        sampleFrom gen var
 
     -- Ideally, these would all be the same speed:
     , bench "sum of samples (implicit rvar)" $ nfIO $ do
-        sumM count (sampleFrom src var)
+        sumM count (sampleFrom gen var)
 
     , bench "sum of samples (explicit rvar)" $ nfIO $ do
-        sumM count (sampleFrom src (rvar var))
+        sumM count (sampleFrom gen (rvar var))
 
     , bench "sample of sum" $ nfIO $ do
-        sampleFrom src (sumM count (rvar var))
+        sampleFrom gen (sumM count (rvar var))
 
     , bench "array of samples" $ nfIO $ do
         allocaArray count $ \ptr -> do
             sequence_
                 [ do
-                    x <- sampleFrom src var
+                    x <- sampleFrom gen var
 
                     pokeElemOff ptr offset x
                 | offset <- [0 .. count-1]
@@ -199,7 +197,7 @@ suite src count name var = bgroup name
             sumBuf count ptr
 
     , bench "RVarT IO arrays" $ nfIO $ do
-        allocaArray count $ \ptr -> flip runRVarT src $ do
+        allocaArray count $ \ptr -> flip runRVarT gen $ do
             sequence_
                 [ do
                     x <- rvarT var
diff --git a/tests/speed/TestSupport.hs b/tests/speed/TestSupport.hs
index 71aae7f..20d5e68 100644
--- a/tests/speed/TestSupport.hs
+++ b/tests/speed/TestSupport.hs
@@ -2,41 +2,25 @@
 
 {-# OPTIONS_GHC -fno-warn-simplifiable-class-constraints #-}
 {-# OPTIONS_GHC -fno-warn-type-defaults                  #-}
-{-# OPTIONS_GHC -fno-warn-missing-signatures              #-}
+{-# OPTIONS_GHC -fno-warn-missing-signatures             #-}
 
 module TestSupport where
 
 import System.Random.Mersenne.Pure64
 import System.Random.MWC
 import Data.List
-import Data.StateRef
 import Foreign
-import System.Random
+import System.Random.Stateful
+import Control.Monad.ST (RealWorld)
 
--- type Src = IORef PureMT
--- getTestSource = newMTSrc
+newMTGenM :: IO (IOGenM PureMT)
+newMTGenM = newIOGenM =<< newPureMT
 
--- type Src = IORef StdGen
--- getTestSource = newStdSrc
+newStdGenM :: IO (IOGenM StdGen)
+newStdGenM = newIOGenM =<< newStdGen
 
-type Src = Gen RealWorld
-getTestSource :: IO Src
-getTestSource = newGenIO
-
-newMTSrc :: IO (IORef PureMT)
-newMTSrc = do
-    mt <- newPureMT
-    newReference mt
-
-newStdSrc :: IO (IORef StdGen)
-newStdSrc = do
-    mt <- newStdGen
-    newReference mt
-
-newGenIO :: IO (Gen RealWorld)
-newGenIO = do
-    seed <- withSystemRandom (save :: Gen RealWorld -> IO Seed)
-    restore seed
+newMWCGenIO :: IO (Gen RealWorld)
+newMWCGenIO = createSystemRandom
 
 
 sum' xs = foldl' (+) 0 xs
diff --git a/tests/speed/speed-tests.cabal b/tests/speed/speed-tests.cabal
index 1d1ea8b..93e7e24 100644
--- a/tests/speed/speed-tests.cabal
+++ b/tests/speed/speed-tests.cabal
@@ -13,15 +13,12 @@ description:            Various benchmarks for the random-fu library.
 flag split-random-fu
 
 benchmark random-fu-bench
-  type:             exitcode-stdio-1.0
+  type:                 exitcode-stdio-1.0
   main-is:              Bench.hs
   build-depends:        base >= 4, criterion, MonadRandom, mtl,
-                        stateref, mersenne-random-pure64, mwc-random >= 0.5,
-                        random, vector, deepseq
-  if flag(split-random-fu)
-    build-depends:      random-fu, random-source
-  else
-    build-depends:      random-fu
+                        stateref, mersenne-random-pure64, mwc-random >= 0.15.0.1,
+                        random, vector, deepseq,
+                        random-fu
  other-modules:         TestSupport