Repository: PMEAL/porespy

Issue:
add seed argument to all generators
This would be handy for several reason:
- so that we don't have to import numpy just to initialize the rng
- numba doesn't respect the numpy rng so needs special treatment
- the docs could explain how to get a repeatable image 

Discussion:

Patch:
diff --git a/.github/workflows/gh-pages.yml b/.github/workflows/gh-pages.yml
index c5a81e49c..14230c7ac 100644
--- a/.github/workflows/gh-pages.yml
+++ b/.github/workflows/gh-pages.yml
@@ -22,30 +22,21 @@ jobs:
         with:
           python-version: 3.8
 
-      - name: Cache conda
+      - name: Cache pip
         uses: actions/cache@v2
-        env:
-          CACHE_NUMBER: 0  # increase value -> force reset cache
         with:
-          path: ~/conda_pkgs_dir
-          key: ${{ runner.os }}-conda-${{ env.CACHE_NUMBER }}-${{ hashFiles('requirements.txt') }}
-
-      - name: Set up conda/mamba
-        uses: conda-incubator/setup-miniconda@v2
-        with:
-          mamba-version: "*"
-          allow-softlinks: true
-          use-only-tar-bz2: true
-          show-channel-urls: false
-          auto-activate-base: true
-          channel-priority: strict
-          channels: conda-forge,defaults
-          python-version: 3.8
+          # This path is specific to Ubuntu
+          path: ~/.cache/pip
+          # Look to see if there is a cache hit for the corresponding requirements file
+          key: ${{ runner.os }}-pip-${{ hashFiles('requirements.txt') }}
+          restore-keys: |
+            ${{ runner.os }}-pip-
+            ${{ runner.os }}-
 
       - name: Install dependencies (conda)
         run: |
-          mamba install --file requirements/docs.txt
           pip install -r requirements.txt
+          pip install -r requirements/docs.txt
 
       # Build the documentation
       - name: Build the documentation
diff --git a/porespy/__version__.py b/porespy/__version__.py
index f16768b4c..cb1c98223 100644
--- a/porespy/__version__.py
+++ b/porespy/__version__.py
@@ -1,1 +1,1 @@
-__version__ = '2.2.2.dev13'
+__version__ = '2.2.2.dev12'
diff --git a/porespy/generators/__init__.py b/porespy/generators/__init__.py
index 14e86326e..2bef4b9f3 100644
--- a/porespy/generators/__init__.py
+++ b/porespy/generators/__init__.py
@@ -14,6 +14,7 @@ debugging, and illustration purposes.
    :toctree: generated/
 
     generators.rsa
+    generators.random_spheres
     generators.blobs
     generators.borders
     generators.bundle_of_tubes
@@ -44,6 +45,7 @@ from ._imgen import overlapping_spheres
 from ._imgen import polydisperse_spheres
 from ._imgen import RSA
 from ._imgen import rsa
+from ._imgen import random_spheres
 from ._imgen import voronoi_edges
 from ._pseudo_packings import pseudo_gravity_packing
 from ._pseudo_packings import pseudo_electrostatic_packing
diff --git a/porespy/generators/_fractals.py b/porespy/generators/_fractals.py
index 4243fd04a..ec88ec4b8 100644
--- a/porespy/generators/_fractals.py
+++ b/porespy/generators/_fractals.py
@@ -16,7 +16,7 @@ __all__ = [
 ]
 
 
-def random_cantor_dust(shape, n=5, p=2, f=0.8):
+def random_cantor_dust(shape, n=5, p=2, f=0.8, seed=None):
     r"""
     Generates an image of random cantor dust
 
@@ -31,6 +31,10 @@ def random_cantor_dust(shape, n=5, p=2, f=0.8):
         The number of divisions to make on each iteration.
     f : float (default = 0.8)
         The fraction of the set to keep on each iteration.
+    seed : int, optional, default = `None`
+        Initializes numpy's random number generator to the specified state. If not
+        provided, the current global value is used. This means calls to
+        ``np.random.state(seed)`` prior to calling this function will be respected.
 
     Returns
     -------
@@ -44,6 +48,8 @@ def random_cantor_dust(shape, n=5, p=2, f=0.8):
     to view online example.
 
     """
+    if seed is not None:
+        np.random.seed(seed)
     # Parse the given shape and adjust if necessary
     shape = np.array(shape)
     trim = np.mod(shape, (p**n))
diff --git a/porespy/generators/_imgen.py b/porespy/generators/_imgen.py
index 2923959e0..cf968d263 100644
--- a/porespy/generators/_imgen.py
+++ b/porespy/generators/_imgen.py
@@ -110,16 +110,134 @@ def RSA(*args, **kwargs):
     return rsa(*args, **kwargs)
 
 
-def rsa(im_or_shape: np.array,
-        r: int,
-        volume_fraction: int = 1,
-        clearance: int = 0,
-        protrusion: int = 0,
-        n_max: int = 100000,
-        mode: str = "contained",
-        return_spheres: bool = False,
-        smooth: bool = True,
-        seed: int = None):
+def random_spheres(
+    im_or_shape: np.array,
+    r: int,
+    volume_fraction: int = 1,
+    clearance: int = 0,
+    protrusion: int = 0,
+    n_max: int = 100000,
+    edges: str = "contained",
+    return_spheres: bool = False,
+    smooth: bool = True,
+    seed: int = None,
+):
+    r"""
+    Generates a sphere or disk packing using Random Sequential Addition
+
+    Parameters
+    ----------
+    im_or_shape : ndarray or list
+        To provide flexibility, this argument accepts either an image into
+        which the spheres are inserted, or a shape which is used to create an
+        empty image.  In both cases the spheres are added as ``True`` values
+        to the background.  Since ``True`` is considered the pore space, then
+        the added spheres represent holes.
+    r : int
+        The radius of the disk or sphere to insert.
+    volume_fraction : scalar (default is 1.0)
+        The fraction of the image that should be filled with spheres.  The
+        spheres are added as ``True``'s, so each sphere addition increases the
+        ``volume_fraction`` until the specified limit is reached.  Note that if
+        ``n_max`` is reached first, then ``volume_fraction`` will not be
+        achieved.  Also, ``volume_fraction`` is not counted correctly if the
+        ``mode`` is ``'extended'``.
+    clearance : int (optional, default = 0)
+        The amount of space to put between each sphere. Negative values are
+        acceptable to create overlaps, so long as ``abs(clearance) < r``.
+    protrusion : int (optional, default = 0)
+        The amount by which inserted spheres are allowed to protrude outside of
+        the given background.  If set to 0 (the default) then all spheres will
+        be fully inside the region marked ``False`` in the input image. If > 0, then
+        spheres will extend into the region marked ``True`` in the input image.
+    n_max : int (default is 100,000)
+        The maximum number of spheres to add.  Using a low value may halt
+        the addition process prior to reaching the specified
+        ``volume_fraction``.  If ``None`` is given, then no limit is used.
+    edges : string (default is 'contained')
+        Controls how the edges of the image are handled.  Options are:
+
+        ============ ===============================================================
+        edges        description
+        ============ ===============================================================
+        'contained'  Spheres are all completely within the image
+        'extended'   Spheres are allowed to extend beyond the edge of the
+                     image.  In this mode the volume fraction will be less than
+                     requested since some spheres extend beyond the image, but their
+                     entire volume is counted as added for computational efficiency.
+        ============ ===============================================================
+
+    return_spheres : bool
+        If ``True`` then an image containing only the spheres is returned
+        rather than the input image with the spheres added, which is the
+        default behavior.
+    smooth : bool
+        Indicates whether balls should have smooth faces (``True``) or should
+        include the bumps on the extremities (``False``).
+    seed : int
+        The seed to supply to the random number generators. Because this function
+        uses ``numba`` for speed, calling the normal ``numpy.random.seed(<seed>)``
+        has no effect. To get a repeatable image, the seed must be passed to the
+        function so it can be initialized the way ``numba`` requires. The default
+        is ``None``, which means each call will produce a new realization.
+
+    Returns
+    -------
+    image : ndarray
+        An image with spheres of specified radius *added* to the background.
+
+    See Also
+    --------
+    pseudo_gravity_packing
+    pseudo_electrostatic_packing
+
+    Notes
+    -----
+    This algorithm ensures that spheres do not overlap but does not
+    guarantee they are tightly packed.
+
+    This function adds spheres to the background of the received ``im``, which
+    allows iteratively adding spheres of different radii to the unfilled space
+    by repeatedly passing in the result of previous calls to RSA.
+
+    References
+    ----------
+    [1] Random Heterogeneous Materials, S. Torquato (2001)
+
+    Examples
+    --------
+    `Click here
+    <https://porespy.org/examples/generators/reference/random_spheres.html>`_
+    to view online example.
+
+    """
+    im = rsa(
+        im_or_shape=im_or_shape,
+        r=r,
+        volume_fraction=volume_fraction,
+        clearance=clearance,
+        protrusion=protrusion,
+        n_max=n_max,
+        mode=edges,
+        return_spheres=return_spheres,
+        smooth=smooth,
+        seed=seed)
+    return im
+
+
+@deprecated("This function will be renamed random_spheres in a future version")
+def rsa(
+    im_or_shape: np.array,
+    r: int,
+    volume_fraction: int = 1,
+    clearance: int = 0,
+    protrusion: int = 0,
+    n_max: int = 100000,
+    mode: str = "contained",
+    return_spheres: bool = False,
+    smooth: bool = True,
+    seed: int = None,
+):
     r"""
     Generates a sphere or disk packing using Random Sequential Addition
 
@@ -355,7 +473,13 @@ def _make_choice(options_im, free_sites):
     return coords, count
 
 
-def bundle_of_tubes(shape: List[int], spacing: int, distribution=None, smooth=True):
+def bundle_of_tubes(
+    shape: List[int],
+    spacing: int,
+    distribution=None,
+    smooth=True,
+    seed=None,
+):
     r"""
     Create a 3D image of a bundle of tubes, in the form of a rectangular
     plate with randomly sized holes through it.
@@ -371,6 +495,10 @@ def bundle_of_tubes(shape: List[int], spacing: int, distribution=None, smooth=Tr
         be distributed between this values down to 3 voxels.
     distribution : scipy.stats object
         A handle to a scipy stats object with the desired parameters.
+    seed : int, optional, default = `None`
+        Initializes numpy's random number generator to the specified state. If not
+        provided, the current global value is used. This means calls to
+        ``np.random.state(seed)`` prior to calling this function will be respected.
 
     Returns
     -------
@@ -384,6 +512,8 @@ def bundle_of_tubes(shape: List[int], spacing: int, distribution=None, smooth=Tr
     to view online example.
 
     """
+    if seed is not None:
+        np.random.seed(seed)
     shape = np.array(shape)
     if len(shape) == 2:
         shape = np.hstack((shape, [1]))
@@ -415,7 +545,8 @@ def polydisperse_spheres(shape: List[int],
                          porosity: float,
                          dist,
                          nbins: int = 5,
-                         r_min: int = 5):
+                         r_min: int = 5,
+                         seed=None):
     r"""
     Create an image of randomly placed, overlapping spheres with a
     distribution of radii.
@@ -443,6 +574,10 @@ def polydisperse_spheres(shape: List[int],
         monodisperse spheres that span 0.05 and 0.95 of the possible
         values produced by the provided distribution, then overlays them
         to get polydispersivity.
+    seed : int, optional, default = `None`
+        Initializes numpy's random number generator to the specified state. If not
+        provided, the current global value is used. This means calls to
+        ``np.random.state(seed)`` prior to calling this function will be respected.
 
     Returns
     -------
@@ -456,6 +591,8 @@ def polydisperse_spheres(shape: List[int],
     to view online example.
 
     """
+    if seed is not None:
+        np.random.seed(seed)
     shape = np.array(shape)
     if np.size(shape) == 1:
         shape = np.full((3,), int(shape))
@@ -473,8 +610,13 @@ def polydisperse_spheres(shape: List[int],
     return im
 
 
-def voronoi_edges(shape: List[int], ncells: int, r: int = 0,
-                  flat_faces: bool = True):
+def voronoi_edges(
+    shape: List[int],
+    ncells: int,
+    r: int = 0,
+    flat_faces: bool = True,
+    seed=None,
+):
     r"""
     Create an image from the edges of a Voronoi tessellation.
 
@@ -493,6 +635,10 @@ def voronoi_edges(shape: List[int], ncells: int, r: int = 0,
         Whether the Voronoi edges should lie on the boundary of the
         image (``True``), or if edges outside the image should be removed
         (``False``).
+    seed : int, optional, default = `None`
+        Initializes numpy's random number generator to the specified state. If not
+        provided, the current global value is used. This means calls to
+        ``np.random.state(seed)`` prior to calling this function will be respected.
 
     Returns
     -------
@@ -506,6 +652,8 @@ def voronoi_edges(shape: List[int], ncells: int, r: int = 0,
     to view online example.
 
     """
+    if seed is not None:
+        np.random.seed(seed)
     logger.info(f"Generating {ncells} cells")
     shape = np.array(shape)
     if np.size(shape) == 1:
@@ -709,11 +857,14 @@ def lattice_spheres(shape: List[int],
     return im
 
 
-def overlapping_spheres(shape: List[int],
-                        r: int,
-                        porosity: float,
-                        maxiter: int = 10,
-                        tol: float = 0.01):
+def overlapping_spheres(
+    shape: List[int],
+    r: int,
+    porosity: float,
+    maxiter: int = 10,
+    tol: float = 0.01,
+    seed=None,
+):
     r"""
     Generate a packing of overlapping mono-disperse spheres
 
@@ -731,6 +882,10 @@ def overlapping_spheres(shape: List[int],
         the porosity of the final image to match the given value.
     tol : float
         Tolerance for porosity of the final image compared to the given value.
+    seed : int, optional, default = `None`
+        Initializes numpy's random number generator to the specified state. If not
+        provided, the current global value is used. This means calls to
+        ``np.random.state(seed)`` prior to calling this function will be respected.
 
     Returns
     -------
@@ -750,6 +905,8 @@ def overlapping_spheres(shape: List[int],
     to view online example.
 
     """
+    if seed is not None:
+        np.random.seed(seed)
     shape = np.array(shape)
     if np.size(shape) == 1:
         shape = np.full((3, ), int(shape))
@@ -796,8 +953,13 @@ def overlapping_spheres(shape: List[int],
     return ~f(N)
 
 
-def blobs(shape: List[int], porosity: float = 0.5, blobiness: int = 1,
-          divs: int = 1):
+def blobs(
+    shape: List[int],
+    porosity: float = 0.5,
+    blobiness: int = 1,
+    divs: int = 1,
+    seed=None,
+):
     """
     Generates an image containing amorphous blobs
 
@@ -821,6 +983,10 @@ def blobs(shape: List[int], porosity: float = 0.5, blobiness: int = 1,
         equivalent to ``[2, 2, 2]`` for a 3D image.  The number of cores
         used is specified in ``porespy.settings.ncores`` and defaults to
         all cores.
+    seed : int, optional, default = `None`
+        Initializes numpy's random number generator to the specified state. If not
+        provided, the current global value is used. This means calls to
+        ``np.random.state(seed)`` prior to calling this function will be respected.
 
     Returns
     -------
@@ -848,6 +1014,8 @@ def blobs(shape: List[int], porosity: float = 0.5, blobiness: int = 1,
     to view online example.
 
     """
+    if seed is not None:
+        np.random.seed(seed)
     if isinstance(shape, int):
         shape = [shape]*3
     if len(shape) == 1:
@@ -877,13 +1045,16 @@ def blobs(shape: List[int], porosity: float = 0.5, blobiness: int = 1,
     return im
 
 
-def _cylinders(shape: List[int],
-               r: int,
-               ncylinders: int,
-               phi_max: float = 0,
-               theta_max: float = 90,
-               length: float = None,
-               verbose: bool = True):
+def _cylinders(
+    shape: List[int],
+    r: int,
+    ncylinders: int,
+    phi_max: float = 0,
+    theta_max: float = 90,
+    length: float = None,
+    verbose: bool = True,
+    seed=None,
+):
     r"""
     Generates a binary image of overlapping cylinders.
 
@@ -918,6 +1089,10 @@ def _cylinders(shape: List[int],
         cylinder. Note that one or both of the ends *may* still lie
         outside the domain, depending on the randomly chosen center point
         of the cylinder.
+    seed : int, optional, default = `None`
+        Initializes numpy's random number generator to the specified state. If not
+        provided, the current global value is used. This means calls to
+        ``np.random.state(seed)`` prior to calling this function will be respected.
 
     Returns
     -------
@@ -925,6 +1100,8 @@ def _cylinders(shape: List[int],
         A boolean array with ``True`` values denoting the pore space
 
     """
+    if seed is not None:
+        np.random.seed(seed)
     shape = np.array(shape)
     if np.size(shape) == 1:
         shape = np.full((3, ), int(shape))
@@ -972,14 +1149,17 @@ def _cylinders(shape: List[int],
     return ~dt
 
 
-def cylinders(shape: List[int],
-              r: int,
-              ncylinders: int = None,
-              porosity: float = None,
-              phi_max: float = 0,
-              theta_max: float = 90,
-              length: float = None,
-              maxiter: int = 3):
+def cylinders(
+    shape: List[int],
+    r: int,
+    ncylinders: int = None,
+    porosity: float = None,
+    phi_max: float = 0,
+    theta_max: float = 90,
+    length: float = None,
+    maxiter: int = 3,
+    seed=None,
+):
     r"""
     Generates a binary image of overlapping cylinders given porosity OR
     number of cylinders.
@@ -1026,6 +1206,10 @@ def cylinders(shape: List[int],
         typically effective in getting very close to the targeted
         porosity), but a greater number can be input to improve the
         achieved porosity.
+    seed : int, optional, default = `None`
+        Initializes numpy's random number generator to the specified state. If not
+        provided, the current global value is used. This means calls to
+        ``np.random.state(seed)`` prior to calling this function will be respected.
 
     Returns
     -------
@@ -1058,6 +1242,9 @@ def cylinders(shape: List[int],
     to view online example.
 
     """
+    if seed is not None:
+        np.random.seed(seed)
+
     if ncylinders is not None:
         im = _cylinders(
             shape=shape,
diff --git a/porespy/generators/_noise.py b/porespy/generators/_noise.py
index a86ec45ad..de6cfc179 100644
--- a/porespy/generators/_noise.py
+++ b/porespy/generators/_noise.py
@@ -3,8 +3,16 @@ from porespy.tools import norm_to_uniform
 import psutil
 
 
-def fractal_noise(shape, frequency=0.05, octaves=4, gain=0.5, mode='simplex',
-                  seed=None, cores=None, uniform=True):
+def fractal_noise(
+    shape,
+    frequency=0.05,
+    octaves=4,
+    gain=0.5,
+    mode='simplex',
+    seed=None,
+    cores=None,
+    uniform=True,
+):
     r"""
     Generate fractal noise which can be thresholded to create binary
     images with realistic structures across scales.
diff --git a/porespy/generators/_pseudo_packings.py b/porespy/generators/_pseudo_packings.py
index 7cca1551e..619e42842 100644
--- a/porespy/generators/_pseudo_packings.py
+++ b/porespy/generators/_pseudo_packings.py
@@ -1,5 +1,4 @@
 import logging
-import numba
 import numpy as np
 import scipy.ndimage as spim
 from edt import edt
@@ -8,15 +7,28 @@ from porespy import settings
 from porespy.tools import get_tqdm, ps_round, get_border
 from porespy.tools import _insert_disks_at_points
 from porespy.filters import trim_disconnected_blobs, fftmorphology
-import random
+# import random
+from numba import njit
 
 
 tqdm = get_tqdm()
 logger = logging.getLogger(__name__)
 
 
-def pseudo_gravity_packing(im, r, clearance=0, axis=0, maxiter=1000,
-                           edges='contained'):
+@njit
+def _set_seed(a):
+    np.random.seed(a)
+
+
+def pseudo_gravity_packing(
+    im,
+    r,
+    clearance=0,
+    axis=0,
+    maxiter=1000,
+    edges='contained',
+    seed=None,
+):
     r"""
     Iteratively inserts spheres at the lowest accessible point in an image,
     mimicking a gravity packing.
@@ -37,15 +49,24 @@ def pseudo_gravity_packing(im, r, clearance=0, axis=0, maxiter=1000,
     maxiter : int (default is 1000)
         The maximum number of spheres to add
     edges : string (default is 'contained')
-        Controls how the edges of the image are handled.  Options are:
+        Controls how spheres at the edges of the image are handled.  Options are:
+
+        ============ ================================================================
+        edges        description
+        ============ ================================================================
+        'contained'  Spheres are all completely within the image
+        'extended'   Spheres are allowed to extend beyond the edge of the
+                     image.  In this mode the volume fraction will be less that
+                     requested since some spheres extend beyond the image, but their
+                     entire volume is counted as added for computational efficiency.
+        ============ ================================================================
 
-        'contained'
-            Spheres are all completely within the image
-        'extended'
-            Spheres are allowed to extend beyond the edge of the
-            image.  In this mode the volume fraction will be less that
-            requested since some spheres extend beyond the image, but their
-            entire volume is counted as added for computational efficiency.
+    seed : int, optional, default = `None`
+        The seed to supply to the random number generator. Because this function
+        uses ``numba`` for speed, calling the normal ``numpy.random.seed(<seed>)``
+        has no effect. To get a repeatable image, the seed must be passed to the
+        function so it can be initialized the way ``numba`` requires. The default
+        is ``None``, which means each call will produce a new realization.
 
     Returns
     -------
@@ -63,6 +84,10 @@ def pseudo_gravity_packing(im, r, clearance=0, axis=0, maxiter=1000,
     """
     logger.debug(f'Adding spheres of radius {r}')
 
+    if seed is not None:  # Initialize rng so numba sees it
+        _set_seed(seed)
+        np.random.seed(seed)
+
     im = np.swapaxes(im, 0, axis)
     im_temp = np.zeros_like(im, dtype=bool)
     r = r - 1
@@ -102,11 +127,16 @@ def pseudo_gravity_packing(im, r, clearance=0, axis=0, maxiter=1000,
     return im_temp
 
 
-def pseudo_electrostatic_packing(im, r, sites=None,
-                                 clearance=0,
-                                 protrusion=0,
-                                 edges='extended',
-                                 maxiter=1000):
+def pseudo_electrostatic_packing(
+    im,
+    r,
+    sites=None,
+    clearance=0,
+    protrusion=0,
+    edges='extended',
+    maxiter=1000,
+    seed=None,
+):
     r"""
     Iterativley inserts spheres as close to the given sites as possible.
 
@@ -129,15 +159,24 @@ def pseudo_electrostatic_packing(im, r, sites=None,
     maxiter : int (optional, default=1000)
         The maximum number of spheres to insert.
     edges : string (default is 'contained')
-        Controls how the edges of the image are handled.  Options are:
+        Controls how spheres at the edges of the image are handled.  Options are:
+
+        ============ ================================================================
+        edges        description
+        ============ ================================================================
+        'contained'  Spheres are all completely within the image
+        'extended'   Spheres are allowed to extend beyond the edge of the
+                     image.  In this mode the volume fraction will be less that
+                     requested since some spheres extend beyond the image, but their
+                     entire volume is counted as added for computational efficiency.
+        ============ ================================================================
 
-        'contained'
-            Spheres are all completely within the image
-        'extended'
-            Spheres are allowed to extend beyond the edge of the
-            image.  In this mode the volume fraction will be less that
-            requested since some spheres extend beyond the image, but their
-            entire volume is counted as added for computational efficiency.
+    seed : int, optional, default = `None`
+        The seed to supply to the random number generator. Because this function
+        uses ``numba`` for speed, calling the normal ``numpy.random.seed(<seed>)``
+        has no effect. To get a repeatable image, the seed must be passed to the
+        function so it can be initialized the way ``numba`` requires. The default
+        is ``None``, which means each call will produce a new realization.
 
     Returns
     -------
@@ -151,7 +190,10 @@ def pseudo_electrostatic_packing(im, r, sites=None,
     to view online example.
 
     """
-    random.seed(0)
+    if seed is not None:  # Initialize rng so numba sees it
+        _set_seed(seed)
+        np.random.seed(seed)
+
     im_temp = np.zeros_like(im, dtype=bool)
     dt_im = edt(im)
     if sites is None:
@@ -180,7 +222,7 @@ def pseudo_electrostatic_packing(im, r, sites=None,
             hits = dt[options] < dtmax
         if hits.size == 0:
             break
-        choice = random.choice(np.where(hits)[0])
+        choice = np.random.choice(np.where(hits)[0])
         cen = np.vstack([options[i][choice] for i in range(im.ndim)])
         im_temp = _insert_disks_at_points(im_temp, coords=cen,
                                           radii=np.array([r-clearance]),
@@ -201,7 +243,7 @@ if __name__ == "__main__":
     fig, ax = plt.subplots(1, 2)
     im = ps.generators.pseudo_gravity_packing(im=np.ones(shape, dtype=bool),
                                               r=7, clearance=3,
-                                              edges='contained')
+                                              edges='contained', seed=0)
     ax[0].imshow(im, origin='lower')
 
     sites = np.zeros(shape, dtype=bool)
@@ -210,5 +252,5 @@ if __name__ == "__main__":
                                                                dtype=bool),
                                                     r=5, sites=sites,
                                                     clearance=4,
-                                                    maxiter=50)
+                                                    maxiter=50, seed=0)
     ax[1].imshow(im, origin='lower')
diff --git a/requirements/docs.txt b/requirements/docs.txt
index 6edb595a5..3bfbf7d38 100644
--- a/requirements/docs.txt
+++ b/requirements/docs.txt
@@ -5,7 +5,7 @@ jupyterlab_widgets
 mock
 myst-nb
 pandoc
-pydata-sphinx-theme=0.9
+pydata-sphinx-theme==0.9
 sphinx
 sphinx-copybutton
 sphinx-design
diff --git a/setup.cfg b/setup.cfg
index bdeca8ab8..6d7a9b872 100644
--- a/setup.cfg
+++ b/setup.cfg
@@ -1,5 +1,5 @@
 [bumpversion]
-current_version = 2.2.2.dev13
+current_version = 2.2.2.dev12
 parse = (?P<major>\d+)\.(?P<minor>\d+)\.(?P<patch>\d+)\.(?P<release>\D+)(?P<build>\d+)?
 serialize = {major}.{minor}.{patch}.{release}{build}
 


Test patch:
diff --git a/test/unit/test_generators.py b/test/unit/test_generators.py
index a015cc785..6e817ec93 100644
--- a/test/unit/test_generators.py
+++ b/test/unit/test_generators.py
@@ -33,6 +33,13 @@ class GeneratorTest():
         porosity = im.sum() / im.size
         np.testing.assert_allclose(porosity, 0.5, rtol=1e-2)
 
+    def test_cylinders_w_seed(self):
+        im1 = ps.generators.cylinders(shape=[50, 50, 50], r=1, ncylinders=20, seed=0)
+        im2 = ps.generators.cylinders(shape=[50, 50, 50], r=1, ncylinders=20, seed=0)
+        im3 = ps.generators.cylinders(shape=[50, 50, 50], r=1, ncylinders=20, seed=1)
+        assert np.all(im1 == im2)
+        assert not np.all(im1 == im3)
+
     def test_insert_shape_center_defaults(self):
         im = np.zeros([11, 11])
         shape = np.ones([3, 3])
@@ -145,6 +152,13 @@ class GeneratorTest():
         assert N == 100
         assert im.shape == (101, 101, 1)
 
+    def test_bundle_of_tubes_w_seed(self):
+        im1 = ps.generators.bundle_of_tubes(shape=[101, 101], spacing=10, seed=0)
+        im2 = ps.generators.bundle_of_tubes(shape=[101, 101], spacing=10, seed=0)
+        im3 = ps.generators.bundle_of_tubes(shape=[101, 101], spacing=10, seed=1)
+        assert np.all(im1 == im2)
+        assert not np.all(im1 == im3)
+
     def test_overlapping_spheres_2d(self):
         phis = np.arange(0.1, 0.9, 0.2)
         for phi in phis:
@@ -161,6 +175,16 @@ class GeneratorTest():
             phi_actual = im.sum() / np.size(im)
             assert abs(phi_actual - phi) < 0.02
 
+    def test_overlapping_spheres_w_seed(self):
+        im1 = ps.generators.overlapping_spheres(shape=[50, 50, 50], r=5,
+                                                porosity=0.5, seed=0)
+        im2 = ps.generators.overlapping_spheres(shape=[50, 50, 50], r=5,
+                                                porosity=0.5, seed=0)
+        im3 = ps.generators.overlapping_spheres(shape=[50, 50, 50], r=5,
+                                                porosity=0.5, seed=1)
+        assert np.all(im1 == im2)
+        assert not np.all(im1 == im3)
+
     def test_polydisperse_spheres(self):
         phis = np.arange(0.1, 0.5, 0.2)
         dist = spst.norm(loc=7, scale=2)
@@ -178,6 +202,16 @@ class GeneratorTest():
         top_slice = im[:, :, 0]
         assert np.sum(top_slice) == 1398
 
+    def test_voronoi_edges_w_seed(self):
+        im1 = ps.generators.voronoi_edges(
+            shape=[50, 50, 50], r=2, ncells=25, seed=0)
+        im2 = ps.generators.voronoi_edges(
+            shape=[50, 50, 50], r=2, ncells=25, seed=0)
+        im3 = ps.generators.voronoi_edges(
+            shape=[50, 50, 50], r=2, ncells=25, seed=1)
+        assert np.all(im1 == im2)
+        assert not np.all(im1 == im3)
+
     def test_lattice_spheres_square(self):
         im = ps.generators.lattice_spheres(
             shape=[101, 101], r=5, spacing=10, lattice='sc')
@@ -212,6 +246,13 @@ class GeneratorTest():
         im = ps.generators.blobs(shape=[101])
         assert len(list(im.shape)) == 3
 
+    def test_blobs_w_seed(self):
+        im1 = ps.generators.blobs(shape=[101, 101], seed=0)
+        im2 = ps.generators.blobs(shape=[101, 101], seed=0)
+        im3 = ps.generators.blobs(shape=[101, 101], seed=1)
+        assert np.all(im1 == im2)
+        assert not np.all(im1 == im3)
+
     def test_rsa_2d_contained(self):
         im = np.zeros([100, 100], dtype=int)
         im = ps.generators.rsa(im, r=10, volume_fraction=0.5, mode='contained')
@@ -290,6 +331,13 @@ class GeneratorTest():
         rsa2p = ps.generators.rsa(im_or_shape=[200, 200], r=1, clearance=2)
         assert rsa0.sum() > rsa2p.sum()
 
+    def test_rsa_w_seed(self):
+        im1 = ps.generators.rsa([50, 50], r=5, seed=0)
+        im2 = ps.generators.rsa([50, 50], r=5, seed=0)
+        im3 = ps.generators.rsa([50, 50], r=5, seed=1)
+        assert np.all(im1 == im2)
+        assert not np.all(im1 == im3)
+
     def test_line_segment(self):
         X0 = [3, 4]
         X1 = [5, 9]
@@ -324,16 +372,27 @@ class GeneratorTest():
         e4 = im.sum()/im.size
         assert e4 > e3
 
-    def test_pseudo_gravity_packing_values(self):
+    def test_pseudo_gravity_packing_2D(self):
         np.random.seed(0)
-        # 2d
-        im = np.ones([50, 50], dtype=bool)
-        im = ps.generators.pseudo_gravity_packing(im=im, r=5, clearance=0)
-        assert_allclose(np.linalg.norm(im), 37.3497, rtol=1e-5)
-        # 3d
-        im = np.ones([50, 50, 50], dtype=bool)
-        im = ps.generators.pseudo_gravity_packing(im=im, r=5, clearance=0)
-        assert_allclose(np.linalg.norm(im), 218.3804, rtol=1e-5)
+        im = np.ones([100, 100], dtype=bool)
+        im = ps.generators.pseudo_gravity_packing(im=im, r=8, clearance=1)
+        assert im.sum() == 4578
+
+    def test_pseudo_gravity_packing_3D(self):
+        np.random.seed(0)
+        im = np.ones([100, 100, 100], dtype=bool)
+        im = ps.generators.pseudo_gravity_packing(im=im, r=8, clearance=1)
+        assert im.sum() == 279240
+
+    def test_pseudo_gravity_packing_w_seed(self):
+        im1 = ps.generators.pseudo_gravity_packing(
+            im=np.ones([50, 50], dtype=bool), r=5, seed=0)
+        im2 = ps.generators.pseudo_gravity_packing(
+            im=np.ones([50, 50], dtype=bool), r=5, seed=0)
+        im3 = ps.generators.pseudo_gravity_packing(
+            im=np.ones([50, 50], dtype=bool), r=5, seed=1)
+        assert np.all(im1 == im2)
+        assert not np.all(im1 == im3)
 
     def test_pseudo_electrostatic_packing(self):
         im1 = ps.generators.blobs(shape=[100, 100])
@@ -342,19 +401,35 @@ class GeneratorTest():
         assert (im1.sum() > im2.sum())
         assert im2.sum() > 0
 
-    def test_pseudo_electrostatic_packing_values(self):
-        np.random.seed(0)
-        # 2d
-        im1 = ps.generators.blobs(shape=[100, 100])
+    def test_pseudo_electrostatic_packing_w_seed(self):
+        im1 = ps.generators.pseudo_electrostatic_packing(
+            im=ps.generators.blobs(shape=[100, 100], seed=0), r=5, seed=0)
         im2 = ps.generators.pseudo_electrostatic_packing(
-            im=im1, r=3, clearance=1, protrusion=1)
-        assert_allclose(np.linalg.norm(im2), 46.604721, rtol=1e-5)
-        # 3d
+            im=ps.generators.blobs(shape=[100, 100], seed=0), r=5, seed=0)
+        im3 = ps.generators.pseudo_electrostatic_packing(
+            im=ps.generators.blobs(shape=[100, 100], seed=0), r=5, seed=1)
+        assert np.all(im1 == im2)
+        assert not np.all(im1 == im3)
+
+    def test_pseudo_electrostatic_packing_2D(self):
         np.random.seed(0)
-        im1 = ps.generators.blobs(shape=[50, 50, 50])
-        im2 = ps.generators.pseudo_electrostatic_packing(
-            im=im1, r=3, clearance=1, protrusion=1)
-        assert_allclose(np.linalg.norm(im2), 130.16528, rtol=1e-5)
+        im = np.ones([100, 100], dtype=bool)
+        sites = np.zeros_like(im)
+        sites[50, 50] = True
+        im = ps.generators.pseudo_electrostatic_packing(im=im, r=8,
+                                                        sites=sites,
+                                                        maxiter=10)
+        assert im.sum() == 1930
+
+    def test_pseudo_electrostatic_packing_3D(self):
+        np.random.seed(0)
+        im = np.ones([100, 100, 100], dtype=bool)
+        sites = np.zeros_like(im)
+        sites[50, 50, 50] = True
+        im = ps.generators.pseudo_electrostatic_packing(im=im, r=8,
+                                                        sites=sites,
+                                                        maxiter=10)
+        assert im.sum() == 21030
 
     def test_faces(self):
         im = ps.generators.faces(shape=[10, 10], inlet=0)
@@ -403,6 +478,13 @@ class GeneratorTest():
         assert im3D.shape == (128, 128, 128)
         np.testing.assert_allclose(im3D.sum()/im3D.size, 0.31614160537720)
 
+    def test_cantor_dust_w_seed(self):
+        im1 = ps.generators.random_cantor_dust([128, 128], 6, f=0.8, seed=0)
+        im2 = ps.generators.random_cantor_dust([128, 128], 6, f=0.8, seed=0)
+        im3 = ps.generators.random_cantor_dust([128, 128], 6, f=0.8, seed=1)
+        assert np.all(im1 == im2)
+        assert not np.all(im1 == im3)
+
     def test_sierpinski_foam(self):
         im2D = ps.generators.sierpinski_foam(4, 4, 2)
         assert im2D.shape == (324, 324)
@@ -470,38 +552,6 @@ class GeneratorTest():
         assert np.all(p.shape == s)
         assert p.sum() == 13690
 
-    def test_pseudo_gravity_packing_2D(self):
-        np.random.seed(0)
-        im = np.ones([100, 100], dtype=bool)
-        im = ps.generators.pseudo_gravity_packing(im=im, r=8, clearance=1)
-        assert im.sum() == 4578
-
-    def test_pseudo_gravity_packing_3D(self):
-        np.random.seed(0)
-        im = np.ones([100, 100, 100], dtype=bool)
-        im = ps.generators.pseudo_gravity_packing(im=im, r=8, clearance=1)
-        assert im.sum() == 279240
-
-    def test_pseudo_electrostatic_packing_2D(self):
-        np.random.seed(0)
-        im = np.ones([100, 100], dtype=bool)
-        sites = np.zeros_like(im)
-        sites[50, 50] = True
-        im = ps.generators.pseudo_electrostatic_packing(im=im, r=8,
-                                                        sites=sites,
-                                                        maxiter=10)
-        assert im.sum() == 1930
-
-    def test_pseudo_electrostatic_packing_3D(self):
-        np.random.seed(0)
-        im = np.ones([100, 100, 100], dtype=bool)
-        sites = np.zeros_like(im)
-        sites[50, 50, 50] = True
-        im = ps.generators.pseudo_electrostatic_packing(im=im, r=8,
-                                                        sites=sites,
-                                                        maxiter=10)
-        assert im.sum() == 21030
-
     def test_spheres_from_coords(self):
         df = pd.DataFrame({'X': [10, 20, 40, 40],
                            'Y': [10, 30, 50, 10],