From 54482d08414d87b96b3f6d92dc349ffe016dda27 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Fri, 12 Apr 2024 23:59:36 -0600 Subject: [PATCH 01/16] First steps at using random generators. --- discretize/tests.py | 52 +++++++++++++-------- discretize/utils/mesh_utils.py | 9 ++-- tests/base/test_interpolation.py | 7 ++- tests/base/test_operators.py | 16 ++++++- tests/base/test_tensor.py | 3 ++ tests/base/test_tensor_innerproduct.py | 11 ++++- tests/boundaries/test_boundary_integrals.py | 5 ++ tests/boundaries/test_boundary_poisson.py | 12 +++++ tests/cyl/test_cyl3D.py | 2 +- tests/tree/test_tree_operators.py | 3 -- 10 files changed, 88 insertions(+), 32 deletions(-) diff --git a/discretize/tests.py b/discretize/tests.py index 501520e12..f5ebbe2b5 100644 --- a/discretize/tests.py +++ b/discretize/tests.py @@ -78,11 +78,10 @@ "You break it, you fix it.", ] -# Initiate random number generator -rng = np.random.default_rng() +_happiness_rng = np.random.default_rng() -def setup_mesh(mesh_type, nC, nDim): +def setup_mesh(mesh_type, nC, nDim, rng=None): """Generate arbitrary mesh for testing. For the mesh type, number of cells along each axis and dimension specified, @@ -106,13 +105,15 @@ def setup_mesh(mesh_type, nC, nDim): discretize.base.BaseMesh A discretize mesh of class specified by the input argument *mesh_type* """ + if "random" in mesh_type: + rng = np.random.default_rng(rng) if "TensorMesh" in mesh_type: if "uniform" in mesh_type: h = [nC, nC, nC] elif "random" in mesh_type: - h1 = np.random.rand(nC) * nC * 0.5 + nC * 0.5 - h2 = np.random.rand(nC) * nC * 0.5 + nC * 0.5 - h3 = np.random.rand(nC) * nC * 0.5 + nC * 0.5 + h1 = rng.random(nC) * nC * 0.5 + nC * 0.5 + h2 = rng.random(nC) * nC * 0.5 + nC * 0.5 + h3 = rng.random(nC) * nC * 0.5 + nC * 0.5 h = [hi / np.sum(hi) for hi in [h1, h2, h3]] # normalize else: raise Exception("Unexpected mesh_type") @@ -127,14 +128,14 @@ def setup_mesh(mesh_type, nC, nDim): else: h = [nC, nC, nC] elif "random" in mesh_type: - h1 = np.random.rand(nC) * nC * 0.5 + nC * 0.5 + h1 = rng.random(nC) * nC * 0.5 + nC * 0.5 if "symmetric" in mesh_type: h2 = [ 2 * np.pi, ] else: - h2 = np.random.rand(nC) * nC * 0.5 + nC * 0.5 - h3 = np.random.rand(nC) * nC * 0.5 + nC * 0.5 + h2 = rng.random(nC) * nC * 0.5 + nC * 0.5 + h3 = rng.random(nC) * nC * 0.5 + nC * 0.5 h = [hi / np.sum(hi) for hi in [h1, h2, h3]] # normalize h[1] = h[1] * 2 * np.pi else: @@ -179,9 +180,9 @@ def setup_mesh(mesh_type, nC, nDim): if "uniform" in mesh_type or "notatree" in mesh_type: h = [nC, nC, nC] elif "random" in mesh_type: - h1 = np.random.rand(nC) * nC * 0.5 + nC * 0.5 - h2 = np.random.rand(nC) * nC * 0.5 + nC * 0.5 - h3 = np.random.rand(nC) * nC * 0.5 + nC * 0.5 + h1 = rng.random(nC) * nC * 0.5 + nC * 0.5 + h2 = rng.random(nC) * nC * 0.5 + nC * 0.5 + h3 = rng.random(nC) * nC * 0.5 + nC * 0.5 h = [hi / np.sum(hi) for hi in [h1, h2, h3]] # normalize else: raise Exception("Unexpected mesh_type") @@ -302,6 +303,7 @@ class for the given operator. Within the test class, the user sets the parameter meshTypes = ["uniformTensorMesh"] _meshType = meshTypes[0] meshDimension = 3 + rng = np.random.default_rng() def setupMesh(self, nC): """Generate mesh and set as current mesh for testing. @@ -316,7 +318,7 @@ def setupMesh(self, nC): Float Maximum cell width for the mesh """ - mesh, max_h = setup_mesh(self._meshType, nC, self.meshDimension) + mesh, max_h = setup_mesh(self._meshType, nC, self.meshDimension, rng=self.rng) self.M = mesh return max_h @@ -335,7 +337,7 @@ def getError(self): """ return 1.0 - def orderTest(self): + def orderTest(self, rng=None): """Perform an order test. For number of cells specified in meshSizes setup mesh, call getError @@ -496,9 +498,9 @@ class in older versions of `discretize`. np.testing.assert_allclose(orders[-1], expected_order, rtol=rtol) elif test_type == "all": np.testing.assert_allclose(orders, expected_order, rtol=rtol) - print(np.random.choice(happiness)) + print(_happiness_rng.choice(happiness)) except AssertionError as err: - print(np.random.choice(sadness)) + print(_happiness_rng.choice(sadness)) raise err return orders @@ -552,6 +554,7 @@ def check_derivative( tolerance=0.85, eps=1e-10, ax=None, + rng=None, ): """Perform a basic derivative check. @@ -580,6 +583,9 @@ def check_derivative( ax : matplotlib.pyplot.Axes, optional An axis object for the convergence plot if *plotIt = True*. Otherwise, the function will create a new axis. + rng : numpy.random.Generator, int, optional + The random number generator to use for the adjoint test, if an integer or None + it used to seed a new `numpy.random.default_rng`. Returns ------- @@ -611,6 +617,8 @@ def check_derivative( __tracebackhide__ = True # matplotlib is a soft dependencies for discretize, # lazy-loaded to decrease load time of discretize. + rng = np.random.default_rng(rng) + try: import matplotlib import matplotlib.pyplot as plt @@ -627,7 +635,7 @@ def check_derivative( x0 = mkvc(x0) if dx is None: - dx = np.random.randn(len(x0)) + dx = rng.standard_normal(len(x0)) h = np.logspace(-1, -num, num) E0 = np.ones(h.shape) @@ -700,7 +708,7 @@ def _plot_it(axes, passed): f" {tolerance} of the expected order {expectedOrder}." ) print("{0!s} PASS! {1!s}".format("=" * 25, "=" * 25)) - print(np.random.choice(happiness) + "\n") + print(_happiness_rng.choice(happiness) + "\n") if plotIt: _plot_it(ax, True) except AssertionError as err: @@ -709,7 +717,7 @@ def _plot_it(axes, passed): "*" * 57, "<" * 25, ">" * 25, "*" * 57 ) ) - print(np.random.choice(sadness) + "\n") + print(_happiness_rng.choice(sadness) + "\n") if plotIt: _plot_it(ax, False) raise err @@ -773,6 +781,7 @@ def assert_isadjoint( rtol=1e-6, atol=0.0, assert_error=True, + rng=None, ): r"""Do a dot product test for the forward operator and its adjoint operator. @@ -823,6 +832,9 @@ def assert_isadjoint( assertion error if failed). If set to False, the result of the test is returned as boolean and a message is printed. + rng : numpy.random.Generator, int, optional + The random number generator to use for the adjoint test, if an integer or None + it used to seed a new `numpy.random.default_rng`. Returns ------- @@ -837,6 +849,8 @@ def assert_isadjoint( """ __tracebackhide__ = True + rng = np.random.default_rng(rng) + def random(size, iscomplex): """Create random data of size and dtype of .""" out = rng.standard_normal(size) diff --git a/discretize/utils/mesh_utils.py b/discretize/utils/mesh_utils.py index 2f5a73891..5c7530298 100644 --- a/discretize/utils/mesh_utils.py +++ b/discretize/utils/mesh_utils.py @@ -28,7 +28,7 @@ def random_model(shape, seed=None, anisotropy=None, its=100, bounds=None): ---------- shape : (dim) tuple of int shape of the model. - seed : int, optional + seed : numpy.random.Generator, int, optional pick which model to produce, prints the seed if you don't choose anisotropy : numpy.ndarray, optional this is the kernel that is convolved with the model @@ -67,15 +67,14 @@ def random_model(shape, seed=None, anisotropy=None, its=100, bounds=None): if bounds is None: bounds = [0, 1] + rng = np.random.default_rng(seed) if seed is None: - seed = np.random.randint(1e3) - print("Using a seed of: ", seed) + print("Using a seed of: ", rng.bit_generator.seed_seq) if type(shape) in num_types: shape = (shape,) # make it a tuple for consistency - np.random.seed(seed) - mr = np.random.rand(*shape) + mr = rng.random(*shape) if anisotropy is None: if len(shape) == 1: smth = np.array([1, 10.0, 1], dtype=float) diff --git a/tests/base/test_interpolation.py b/tests/base/test_interpolation.py index c9826128a..cea6abd02 100644 --- a/tests/base/test_interpolation.py +++ b/tests/base/test_interpolation.py @@ -3,7 +3,7 @@ import discretize -np.random.seed(182) +gen = np.random.default_rng(182) MESHTYPES = ["uniformTensorMesh", "randomTensorMesh"] TOLERANCES = [0.9, 0.5, 0.5] @@ -49,6 +49,7 @@ class TestInterpolation1D(discretize.tests.OrderTest): tolerance = TOLERANCES meshDimension = 1 meshSizes = [8, 16, 32, 64, 128] + rng = gen def getError(self): funX = lambda x: np.cos(2 * np.pi * x) @@ -100,6 +101,7 @@ class TestInterpolation2d(discretize.tests.OrderTest): tolerance = TOLERANCES meshDimension = 2 meshSizes = [8, 16, 32, 64] + rng = gen def getError(self): funX = lambda x, y: np.cos(2 * np.pi * y) @@ -185,6 +187,7 @@ class TestInterpolationSymCyl(discretize.tests.OrderTest): tolerance = 0.6 meshDimension = 3 meshSizes = [32, 64, 128, 256] + rng = gen def getError(self): funX = lambda x, y: np.cos(2 * np.pi * y) @@ -252,6 +255,7 @@ class TestInterpolationCyl(discretize.tests.OrderTest): meshTypes = ["uniformCylMesh", "randomCylMesh"] # MESHTYPES + meshDimension = 3 meshSizes = [8, 16, 32, 64] + rng = gen def getError(self): func = lambda x, y, z: np.cos(2 * np.pi * x) + np.cos(y) + np.cos(2 * np.pi * z) @@ -320,6 +324,7 @@ class TestInterpolation3D(discretize.tests.OrderTest): tolerance = TOLERANCES meshDimension = 3 meshSizes = [8, 16, 32, 64] + rng = gen def getError(self): funX = lambda x, y, z: np.cos(2 * np.pi * y) diff --git a/tests/base/test_operators.py b/tests/base/test_operators.py index aa84d76f5..1679c2a25 100644 --- a/tests/base/test_operators.py +++ b/tests/base/test_operators.py @@ -5,7 +5,7 @@ # Tolerance TOL = 1e-14 -np.random.seed(1) +gen = np.random.default_rng(1) MESHTYPES = [ "uniformTensorMesh", @@ -45,6 +45,7 @@ class TestCurl(discretize.tests.OrderTest): name = "Curl" meshTypes = MESHTYPES + rng = gen def getError(self): # fun: i (cos(y)) + j (cos(z)) + k (cos(x)) @@ -82,6 +83,7 @@ class TestCurl2D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 2 meshSizes = [8, 16, 32, 64] + rng = gen def getError(self): # Test function @@ -108,6 +110,7 @@ def test_order(self): # 1 # because of the averaging involved in the ghost point. u_b = (u_n + u_g)/2 # ) # meshSizes = [8, 16, 32, 64] +# rng = gen # # def getError(self): # # Test function @@ -135,6 +138,7 @@ class TestCellGrad2D_Dirichlet(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 2 meshSizes = [8, 16, 32, 64] + rng = gen def getError(self): # Test function @@ -163,6 +167,7 @@ class TestCellGrad3D_Dirichlet(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 3 meshSizes = [8, 16, 32] + rng = gen def getError(self): # Test function @@ -214,6 +219,7 @@ class TestCellGrad2D_Neumann(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 2 meshSizes = [8, 16, 32, 64] + rng = gen def getError(self): # Test function @@ -242,6 +248,7 @@ class TestCellGrad3D_Neumann(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 3 meshSizes = [8, 16, 32] + rng = gen def getError(self): # Test function @@ -292,6 +299,7 @@ class TestFaceDiv3D(discretize.tests.OrderTest): name = "Face Divergence 3D" meshTypes = MESHTYPES meshSizes = [8, 16, 32, 64] + rng = gen def getError(self): # Test function @@ -327,6 +335,7 @@ class TestFaceDiv2D(discretize.tests.OrderTest): meshTypes = MESHTYPES meshDimension = 2 meshSizes = [8, 16, 32, 64] + rng = gen def getError(self): # Test function @@ -351,6 +360,7 @@ def test_order(self): class TestNodalGrad(discretize.tests.OrderTest): name = "Nodal Gradient" meshTypes = MESHTYPES + rng = gen def getError(self): # Test function @@ -378,6 +388,7 @@ class TestNodalGrad2D(discretize.tests.OrderTest): name = "Nodal Gradient 2D" meshTypes = MESHTYPES meshDimension = 2 + rng = gen def getError(self): # Test function @@ -405,6 +416,7 @@ class TestAveraging1D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 1 meshSizes = [16, 32, 64] + rng = gen def getError(self): num = self.getAve(self.M) * self.getHere(self.M) @@ -526,6 +538,7 @@ class TestAveraging2D(discretize.tests.OrderTest): meshTypes = MESHTYPES meshDimension = 2 meshSizes = [16, 32, 64] + rng = gen def getError(self): num = self.getAve(self.M) * self.getHere(self.M) @@ -645,6 +658,7 @@ class TestAveraging3D(discretize.tests.OrderTest): meshTypes = MESHTYPES meshDimension = 3 meshSizes = [16, 32, 64] + rng = gen def getError(self): num = self.getAve(self.M) * self.getHere(self.M) diff --git a/tests/base/test_tensor.py b/tests/base/test_tensor.py index 9349cafa7..cf42f226f 100644 --- a/tests/base/test_tensor.py +++ b/tests/base/test_tensor.py @@ -6,6 +6,8 @@ TOL = 1e-10 +gen = np.random.default_rng(123) + class BasicTensorMeshTests(unittest.TestCase): def setUp(self): @@ -276,6 +278,7 @@ def test_cell_nodes(self, mesh): class TestPoissonEqn(discretize.tests.OrderTest): name = "Poisson Equation" meshSizes = [10, 16, 20] + rng = gen def getError(self): # Create some functions to integrate diff --git a/tests/base/test_tensor_innerproduct.py b/tests/base/test_tensor_innerproduct.py index 156403255..8ca27543d 100644 --- a/tests/base/test_tensor_innerproduct.py +++ b/tests/base/test_tensor_innerproduct.py @@ -4,7 +4,7 @@ from discretize import TensorMesh from discretize.utils import sdinv -np.random.seed(50) +gen = np.random.default_rng(50) class TestInnerProducts(discretize.tests.OrderTest): @@ -14,6 +14,7 @@ class TestInnerProducts(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh", "uniformCurv", "rotateCurv"] meshDimension = 3 meshSizes = [16, 32] + rng = gen def getError(self): call = lambda fun, xyz: fun(xyz[:, 0], xyz[:, 1], xyz[:, 2]) @@ -173,6 +174,7 @@ class TestInnerProductsFaceProperties3D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 3 meshSizes = [16, 32] + rng = gen def getError(self): call = lambda fun, xyz: fun(xyz[:, 0], xyz[:, 1], xyz[:, 2]) @@ -319,6 +321,7 @@ class TestInnerProductsEdgeProperties3D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 3 meshSizes = [16, 32] + rng = gen def getError(self): call = lambda fun, xyz: fun(xyz[:, 0], xyz[:, 1], xyz[:, 2]) @@ -407,6 +410,7 @@ class TestInnerProducts2D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh", "uniformCurv", "rotateCurv"] meshDimension = 2 meshSizes = [4, 8, 16, 32, 64, 128] + rng = gen def getError(self): z = 5 # Because 5 is just such a great number. @@ -552,6 +556,7 @@ class TestInnerProductsFaceProperties2D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 2 meshSizes = [8, 16, 32] + rng = gen def getError(self): call = lambda fun, xy: fun(xy[:, 0], xy[:, 1]) @@ -642,6 +647,7 @@ class TestInnerProductsEdgeProperties2D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 2 meshSizes = [8, 16, 32] + rng = gen def getError(self): call = lambda fun, xy: fun(xy[:, 0], xy[:, 1]) @@ -704,6 +710,7 @@ class TestInnerProducts1D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 1 meshSizes = [4, 8, 16, 32, 64, 128] + rng = gen def getError(self): y = 12 # Because 12 is just such a great number. @@ -750,7 +757,7 @@ class TestTensorSizeErrorRaises(unittest.TestCase): def setUp(self): self.mesh3D = TensorMesh([4, 4, 4]) - self.model = np.random.rand(self.mesh3D.nC) + self.model = np.ones(self.mesh3D.nC) def test_edge_inner_product_surface(self): self.assertRaises( diff --git a/tests/boundaries/test_boundary_integrals.py b/tests/boundaries/test_boundary_integrals.py index dbe224a76..6574b55dd 100644 --- a/tests/boundaries/test_boundary_integrals.py +++ b/tests/boundaries/test_boundary_integrals.py @@ -3,6 +3,8 @@ import discretize from discretize.utils import cart2cyl, cyl2cart +gen = np.random.default_rng(2552) + def u(*args): if len(args) == 1: @@ -86,6 +88,7 @@ class Test1DBoundaryIntegral(discretize.tests.OrderTest): meshDimension = 1 expectedOrders = 2 meshSizes = [4, 8, 16, 32, 64, 128] + rng = gen def getError(self): mesh = self.M @@ -136,6 +139,7 @@ class Test2DBoundaryIntegral(discretize.tests.OrderTest): meshDimension = 2 expectedOrders = [2, 2, 2, 2, 1] meshSizes = [4, 8, 16, 32, 64, 128] + rng = gen def getError(self): mesh = self.M @@ -216,6 +220,7 @@ class Test3DBoundaryIntegral(discretize.tests.OrderTest): meshDimension = 3 expectedOrders = [2, 1, 2, 2, 2, 2] meshSizes = [4, 8, 16, 32] + rng = gen def getError(self): mesh = self.M diff --git a/tests/boundaries/test_boundary_poisson.py b/tests/boundaries/test_boundary_poisson.py index a88e9924c..25f4f17bd 100644 --- a/tests/boundaries/test_boundary_poisson.py +++ b/tests/boundaries/test_boundary_poisson.py @@ -6,6 +6,8 @@ from discretize import utils from pymatsolver import Solver, Pardiso +gen = np.random.default_rng(42) + class TestCC1D_InhomogeneousDirichlet(discretize.tests.OrderTest): name = "1D - Dirichlet" @@ -13,6 +15,7 @@ class TestCC1D_InhomogeneousDirichlet(discretize.tests.OrderTest): meshDimension = 1 expectedOrders = 2 meshSizes = [4, 8, 16, 32, 64, 128] + rng = gen def getError(self): # Test function @@ -56,6 +59,7 @@ class TestCC2D_InhomogeneousDirichlet(discretize.tests.OrderTest): meshDimension = 2 expectedOrders = [2, 2, 1] meshSizes = [4, 8, 16, 32, 64] + rng = gen def getError(self): # Test function @@ -99,6 +103,7 @@ class TestCC1D_InhomogeneousNeumann(discretize.tests.OrderTest): meshDimension = 1 expectedOrders = 2 meshSizes = [4, 8, 16, 32, 64, 128] + rng = gen def getError(self): # Test function @@ -166,6 +171,7 @@ class TestCC2D_InhomogeneousNeumann(discretize.tests.OrderTest): expectedOrders = [2, 2, 1] meshSizes = [4, 8, 16, 32] # meshSizes = [4] + rng = gen def getError(self): # Test function @@ -228,6 +234,7 @@ class TestCC1D_InhomogeneousMixed(discretize.tests.OrderTest): meshDimension = 1 expectedOrders = 2 meshSizes = [4, 8, 16, 32, 64, 128] + rng = gen def getError(self): # Test function @@ -286,6 +293,7 @@ class TestCC2D_InhomogeneousMixed(discretize.tests.OrderTest): meshDimension = 2 expectedOrders = [2, 2, 1] meshSizes = [2, 4, 8, 16] + rng = gen # meshSizes = [4] def getError(self): @@ -353,6 +361,7 @@ class TestCC3D_InhomogeneousMixed(discretize.tests.OrderTest): meshDimension = 3 expectedOrders = [2, 2, 2] meshSizes = [2, 4, 8, 16, 32] + rng = gen def getError(self): # Test function @@ -441,6 +450,7 @@ class TestN1D_boundaries(discretize.tests.OrderTest): meshDimension = 1 expectedOrders = 2 meshSizes = [2, 4, 8, 16, 32, 64, 128] + rng = gen # meshSizes = [4] def getError(self): @@ -520,6 +530,7 @@ class TestN2D_boundaries(discretize.tests.OrderTest): expectedOrders = 2 tolerance = [0.8, 0.8, 0.6] meshSizes = [8, 16, 32, 64] + rng = gen # meshSizes = [4] def getError(self): @@ -621,6 +632,7 @@ class TestN3D_boundaries(discretize.tests.OrderTest): expectedOrders = 2 tolerance = 0.6 meshSizes = [2, 4, 8, 16, 32] + rng = gen # meshSizes = [4] def getError(self): diff --git a/tests/cyl/test_cyl3D.py b/tests/cyl/test_cyl3D.py index 1c98b60ce..94d9ed6c6 100644 --- a/tests/cyl/test_cyl3D.py +++ b/tests/cyl/test_cyl3D.py @@ -4,7 +4,7 @@ import discretize from discretize import utils -np.random.seed(16) +rng = np.random.default_rng(16) TOL = 1e-1 diff --git a/tests/tree/test_tree_operators.py b/tests/tree/test_tree_operators.py index 22b2e2c8b..8d6d94f4b 100644 --- a/tests/tree/test_tree_operators.py +++ b/tests/tree/test_tree_operators.py @@ -30,9 +30,6 @@ ) ) -# np.random.seed(None) -# np.random.seed(7) - class TestCellGrad2D(discretize.tests.OrderTest): name = "Cell Gradient 2D, using cellGradx and cellGrady" From b719523a48bae590be69eb573910aa193e218bc5 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Thu, 9 May 2024 16:25:44 -0600 Subject: [PATCH 02/16] pass through rng for OrderTest class --- discretize/tests.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/discretize/tests.py b/discretize/tests.py index f5ebbe2b5..057f1a79a 100644 --- a/discretize/tests.py +++ b/discretize/tests.py @@ -305,7 +305,7 @@ class for the given operator. Within the test class, the user sets the parameter meshDimension = 3 rng = np.random.default_rng() - def setupMesh(self, nC): + def setupMesh(self, nC, rng=None): """Generate mesh and set as current mesh for testing. Parameters @@ -318,7 +318,9 @@ def setupMesh(self, nC): Float Maximum cell width for the mesh """ - mesh, max_h = setup_mesh(self._meshType, nC, self.meshDimension, rng=self.rng) + if rng is None: + rng = self.rng + mesh, max_h = setup_mesh(self._meshType, nC, self.meshDimension, rng=rng) self.M = mesh return max_h @@ -363,7 +365,7 @@ def orderTest(self, rng=None): ) def test_func(n_cells): - max_h = self.setupMesh(n_cells) + max_h = self.setupMesh(n_cells, rng=rng) err = self.getError() return err, max_h From 15acf3e314d5af3567a2ab64a2d8c2a1cfa96411 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Thu, 9 May 2024 16:30:21 -0600 Subject: [PATCH 03/16] add bit of docs --- discretize/tests.py | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/discretize/tests.py b/discretize/tests.py index 057f1a79a..378b30e6f 100644 --- a/discretize/tests.py +++ b/discretize/tests.py @@ -99,6 +99,9 @@ def setup_mesh(mesh_type, nC, nDim, rng=None): number of base mesh cells and must be a power of 2. nDim : int The dimension of the mesh. Must be 1, 2 or 3. + rng : numpy.random.Generator, int, optional + The random number generator to use for the adjoint test, if an integer or None + it used to seed a new `numpy.random.default_rng`. Only used if mesh_type is 'random' Returns ------- @@ -237,6 +240,9 @@ class for the given operator. Within the test class, the user sets the parameter for the meshes used in the convergence test; e.g. [4, 8, 16, 32] meshDimension : int Mesh dimension. Must be 1, 2 or 3 + rng : numpy.random.Generator, int, optional + The random number generator to use for the adjoint test, if an integer or None + it used to seed a new `numpy.random.default_rng`. Only used if mesh_type is 'random' Notes ----- @@ -312,6 +318,9 @@ def setupMesh(self, nC, rng=None): ---------- nC : int Number of cells along each axis. + rng : numpy.random.Generator, int, optional + The random number generator to use for the adjoint test, if an integer or None + it used to seed a new `numpy.random.default_rng`. Only used if mesh_type is 'random' Returns ------- From ab74637b3ca004f7590a89e5f0c88cc6ce9ed692 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Sun, 1 Sep 2024 00:19:58 -0600 Subject: [PATCH 04/16] update tests to use random generator simplex passes tree tests passing cyl passing boundaries passing, also removes pymatsolver imports base passes --- discretize/base/base_mesh.py | 24 +++---- discretize/mixins/mpl_mod.py | 4 +- discretize/tests.py | 47 +++++++------ discretize/tree_mesh.py | 3 +- discretize/utils/interpolation_utils.py | 7 +- discretize/utils/matrix_utils.py | 20 +++--- discretize/utils/mesh_utils.py | 3 +- examples/plot_cahn_hilliard.py | 4 +- examples/plot_dc_resistivity.py | 8 +-- tests/base/test_coordutils.py | 19 +++--- tests/base/test_interpolation.py | 32 ++++----- tests/base/test_operators.py | 28 ++------ tests/base/test_tensor.py | 4 +- tests/base/test_tensor_innerproduct.py | 9 --- tests/base/test_tensor_innerproduct_derivs.py | 24 +++---- tests/base/test_tensor_omf.py | 3 +- tests/base/test_tests.py | 9 ++- tests/base/test_utils.py | 45 +++++++------ tests/base/test_view.py | 4 +- tests/base/test_volume_avg.py | 66 +++++++++++-------- tests/boundaries/test_boundary_integrals.py | 6 +- tests/boundaries/test_boundary_maxwell.py | 8 +-- tests/boundaries/test_boundary_poisson.py | 36 ++++------ tests/boundaries/test_errors.py | 9 ++- tests/boundaries/test_tensor_boundary.py | 11 ++-- .../test_tensor_boundary_poisson.py | 11 ++-- tests/cyl/test_cyl.py | 10 +-- tests/cyl/test_cyl3D.py | 2 - tests/cyl/test_cylOperators.py | 4 +- tests/cyl/test_cyl_innerproducts.py | 14 ++-- tests/cyl/test_cyl_operators.py | 6 +- tests/simplex/test_inner_products.py | 14 ++-- tests/simplex/test_utils.py | 20 +++--- tests/tree/test_refine.py | 23 +++---- tests/tree/test_tree.py | 44 +++++++------ tests/tree/test_tree_innerproduct_derivs.py | 32 ++++----- tests/tree/test_tree_interpolation.py | 2 - tests/tree/test_tree_operators.py | 23 +++---- tests/tree/test_tree_plotting.py | 26 ++++---- tests/tree/test_tree_utils.py | 1 - tutorials/inner_products/1_basic.py | 10 +-- .../inner_products/2_physical_properties.py | 14 ++-- tutorials/inner_products/4_advanced.py | 6 +- tutorials/pde/1_poisson.py | 5 +- tutorials/pde/2_advection_diffusion.py | 6 +- 45 files changed, 350 insertions(+), 356 deletions(-) diff --git a/discretize/base/base_mesh.py b/discretize/base/base_mesh.py index 617c270c6..90fb2e331 100644 --- a/discretize/base/base_mesh.py +++ b/discretize/base/base_mesh.py @@ -2419,15 +2419,15 @@ def get_face_inner_product_deriv( >>> import numpy as np >>> import matplotlib as mpl >>> mpl.rcParams.update({'font.size': 14}) - >>> np.random.seed(45) + >>> rng = np.random.default_rng(45) >>> mesh = TensorMesh([[(1, 4)], [(1, 4)]]) Define a model, and a random vector to multiply the derivative with, then we grab the respective derivative function and calculate the sparse matrix, - >>> m = np.random.rand(mesh.nC) # physical property parameters - >>> u = np.random.rand(mesh.nF) # vector of shape (n_faces) + >>> m = rng.random(mesh.nC) # physical property parameters + >>> u = rng.random(mesh.nF) # vector of shape (n_faces) >>> Mf = mesh.get_face_inner_product(m) >>> F = mesh.get_face_inner_product_deriv(m) # Function handle >>> dFdm_u = F(u) @@ -2458,8 +2458,8 @@ def get_face_inner_product_deriv( function handle :math:`\mathbf{F}(\mathbf{u})` and plot the evaluation of this function on a spy plot. - >>> m = np.random.rand(mesh.nC, 3) # anisotropic physical property parameters - >>> u = np.random.rand(mesh.nF) # vector of shape (n_faces) + >>> m = rng.random((mesh.nC, 3)) # anisotropic physical property parameters + >>> u = rng.random(mesh.nF) # vector of shape (n_faces) >>> Mf = mesh.get_face_inner_product(m) >>> F = mesh.get_face_inner_product_deriv(m) # Function handle >>> dFdm_u = F(u) @@ -2602,14 +2602,14 @@ def get_edge_inner_product_deriv( >>> import numpy as np >>> import matplotlib as mpl >>> mpl.rcParams.update({'font.size': 14}) - >>> np.random.seed(45) + >>> rng = np.random.default_rng(45) >>> mesh = TensorMesh([[(1, 4)], [(1, 4)]]) Next we create a random isotropic model vector, and a random vector to multiply the derivative with (for illustration purposes). - >>> m = np.random.rand(mesh.nC) # physical property parameters - >>> u = np.random.rand(mesh.nF) # vector of shape (n_edges) + >>> m = rng.random(mesh.nC) # physical property parameters + >>> u = rng.random(mesh.nF) # vector of shape (n_edges) >>> Me = mesh.get_edge_inner_product(m) >>> F = mesh.get_edge_inner_product_deriv(m) # Function handle >>> dFdm_u = F(u) @@ -2640,8 +2640,8 @@ def get_edge_inner_product_deriv( function handle :math:`\mathbf{F}(\mathbf{u})` and plot the evaluation of this function on a spy plot. - >>> m = np.random.rand(mesh.nC, 3) # physical property parameters - >>> u = np.random.rand(mesh.nF) # vector of shape (n_edges) + >>> m = rng.random((mesh.nC, 3)) # physical property parameters + >>> u = rng.random(mesh.nF) # vector of shape (n_edges) >>> Me = mesh.get_edge_inner_product(m) >>> F = mesh.get_edge_inner_product_deriv(m) # Function handle >>> dFdm_u = F(u) @@ -4128,9 +4128,9 @@ def get_interpolation_matrix( >>> from discretize import TensorMesh >>> import numpy as np >>> import matplotlib.pyplot as plt - >>> np.random.seed(14) + >>> rng = np.random.default_rng(14) - >>> locs = np.random.rand(50)*0.8+0.1 + >>> locs = rng.random(50)*0.8+0.1 >>> dense = np.linspace(0, 1, 200) >>> fun = lambda x: np.cos(2*np.pi*x) diff --git a/discretize/mixins/mpl_mod.py b/discretize/mixins/mpl_mod.py index 331439020..be417829a 100644 --- a/discretize/mixins/mpl_mod.py +++ b/discretize/mixins/mpl_mod.py @@ -470,7 +470,7 @@ def plot_slice( >>> from matplotlib import pyplot as plt >>> import discretize - >>> from pymatsolver import Solver + >>> from scipy.sparse.linalg import spsolve >>> hx = [(5, 2, -1.3), (2, 4), (5, 2, 1.3)] >>> hy = [(2, 2, -1.3), (2, 6), (2, 2, 1.3)] >>> hz = [(2, 2, -1.3), (2, 6), (2, 2, 1.3)] @@ -482,7 +482,7 @@ def plot_slice( >>> q[[4, 4], [4, 4], [2, 6]]=[-1, 1] >>> q = discretize.utils.mkvc(q) >>> A = M.face_divergence * M.cell_gradient - >>> b = Solver(A) * (q) + >>> b = spsolve(A, q) and finaly, plot the vector values of the result, which are defined on faces diff --git a/discretize/tests.py b/discretize/tests.py index adb62a4a4..4848371e4 100644 --- a/discretize/tests.py +++ b/discretize/tests.py @@ -100,8 +100,9 @@ def setup_mesh(mesh_type, nC, nDim, rng=None): nDim : int The dimension of the mesh. Must be 1, 2 or 3. rng : numpy.random.Generator, int, optional - The random number generator to use for the adjoint test, if an integer or None - it used to seed a new `numpy.random.default_rng`. Only used if mesh_type is 'random' + If ``random`` is in `mesh_type`, this is the random number generator to use for + creating a random mesh. If an integer or None it is used to seed a new + `numpy.random.default_rng`. Returns ------- @@ -309,9 +310,9 @@ class for the given operator. Within the test class, the user sets the parameter meshTypes = ["uniformTensorMesh"] _meshType = meshTypes[0] meshDimension = 3 - rng = np.random.default_rng() + rng = None - def setupMesh(self, nC, rng=None): + def setupMesh(self, nC): """Generate mesh and set as current mesh for testing. Parameters @@ -327,9 +328,7 @@ def setupMesh(self, nC, rng=None): Float Maximum cell width for the mesh """ - if rng is None: - rng = self.rng - mesh, max_h = setup_mesh(self._meshType, nC, self.meshDimension, rng=rng) + mesh, max_h = setup_mesh(self._meshType, nC, self.meshDimension, rng=self.rng) self.M = mesh return max_h @@ -373,8 +372,11 @@ def orderTest(self, rng=None): "expectedOrders must have the same length as the meshTypes" ) + if rng is not None: + self.rng = rng + def test_func(n_cells): - max_h = self.setupMesh(n_cells, rng=rng) + max_h = self.setupMesh(n_cells) err = self.getError() return err, max_h @@ -574,8 +576,8 @@ def check_derivative( Parameters ---------- - fctn : function - Function handle + fctn : callable + The function to test. x0 : numpy.ndarray Point at which to check derivative num : int, optional @@ -584,7 +586,7 @@ def check_derivative( If *True*, plot the convergence of the approximation of the derivative dx : numpy.ndarray, optional Step direction. By default, this parameter is set to *None* and a random - step direction is chosen. + step direction is chosen using `rng`. expectedOrder : int, optional The expected order of convergence for the numerical derivative tolerance : float, optional @@ -595,8 +597,9 @@ def check_derivative( An axis object for the convergence plot if *plotIt = True*. Otherwise, the function will create a new axis. rng : numpy.random.Generator, int, optional - The random number generator to use for the adjoint test, if an integer or None - it used to seed a new `numpy.random.default_rng`. + If `dx` is ``None``, This is the random number generator to use for + generating a step direction. If an integer or None, it is used to seed + a new `numpy.random.default_rng`. Returns ------- @@ -608,10 +611,11 @@ def check_derivative( >>> from discretize import tests, utils >>> import numpy as np >>> import matplotlib.pyplot as plt + >>> rng = np.random.default_rng(786412) >>> def simplePass(x): ... return np.sin(x), utils.sdiag(np.cos(x)) - >>> passed = tests.check_derivative(simplePass, np.random.randn(5)) + >>> passed = tests.check_derivative(simplePass, rng.standard_normal(5), rng=rng) ==================== check_derivative ==================== iter h |ft-f0| |ft-f0-h*J0*dx| Order --------------------------------------------------------- @@ -628,7 +632,6 @@ def check_derivative( __tracebackhide__ = True # matplotlib is a soft dependencies for discretize, # lazy-loaded to decrease load time of discretize. - rng = np.random.default_rng(rng) try: import matplotlib @@ -646,6 +649,7 @@ def check_derivative( x0 = mkvc(x0) if dx is None: + rng = np.random.default_rng(rng) dx = rng.standard_normal(len(x0)) h = np.logspace(-1, -num, num) @@ -712,11 +716,14 @@ def _plot_it(axes, passed): # Thus it has no higher order derivatives. pass else: - test = np.mean(order1) > tolerance * expectedOrder + order_mean = np.mean(order1) + expected = tolerance * expectedOrder + test = order_mean > expected if not test: raise AssertionError( - f"\n Order mean {np.mean(order1)} is not greater than" - f" {tolerance} of the expected order {expectedOrder}." + f"\n Order mean {order_mean} is not greater than" + f" {expected} = tolerance: {tolerance} " + f"* expected order: {expectedOrder}." ) print("{0!s} PASS! {1!s}".format("=" * 25, "=" * 25)) print(_happiness_rng.choice(happiness) + "\n") @@ -844,8 +851,8 @@ def assert_isadjoint( returned as boolean and a message is printed. rng : numpy.random.Generator, int, optional - The random number generator to use for the adjoint test, if an integer or None - it used to seed a new `numpy.random.default_rng`. + The random number generator to use for the adjoint test. If an integer or None + it is used to seed a new `numpy.random.default_rng`. Returns ------- diff --git a/discretize/tree_mesh.py b/discretize/tree_mesh.py index fe2db06cf..0478ad79e 100644 --- a/discretize/tree_mesh.py +++ b/discretize/tree_mesh.py @@ -436,7 +436,8 @@ def refine_bounding_box( >>> import matplotlib.pyplot as plt >>> import matplotlib.patches as patches >>> mesh = discretize.TreeMesh([32, 32]) - >>> points = np.random.rand(20, 2) * 0.25 + 3/8 + >>> rng = np.random.default_rng(852) + >>> points = rng.random((20, 2)) * 0.25 + 3/8 Now we want to refine to the maximum level, with no padding the in `x` direction and `2` cells in `y`. At the second highest level we want 2 padding diff --git a/discretize/utils/interpolation_utils.py b/discretize/utils/interpolation_utils.py index 9c6ba159a..2c3c866e2 100644 --- a/discretize/utils/interpolation_utils.py +++ b/discretize/utils/interpolation_utils.py @@ -84,11 +84,11 @@ def interpolation_matrix(locs, x, y=None, z=None): >>> from discretize import TensorMesh >>> import numpy as np >>> import matplotlib.pyplot as plt - >>> np.random.seed(14) + >>> rng = np.random.default_rng(14) Create an interpolation matrix - >>> locs = np.random.rand(50)*0.8+0.1 + >>> locs = rng.random(50)*0.8+0.1 >>> x = np.linspace(0, 1, 7) >>> dense = np.linspace(0, 1, 200) >>> fun = lambda x: np.cos(2*np.pi*x) @@ -217,6 +217,7 @@ def volume_average(mesh_in, mesh_out, values=None, output=None): >>> import numpy as np >>> from discretize import TensorMesh + >>> rng = np.random.default_rng(853) >>> h1 = np.ones(32) >>> h2 = np.ones(16)*2 >>> mesh_in = TensorMesh([h1, h1]) @@ -226,7 +227,7 @@ def volume_average(mesh_in, mesh_out, values=None, output=None): interpolate it to the output mesh. >>> from discretize.utils import volume_average - >>> model1 = np.random.rand(mesh_in.nC) + >>> model1 = rng.random(mesh_in.nC) >>> model2 = volume_average(mesh_in, mesh_out, model1) Because these two meshes' cells are perfectly aligned, but the output mesh diff --git a/discretize/utils/matrix_utils.py b/discretize/utils/matrix_utils.py index 70fd73ef9..75d922723 100644 --- a/discretize/utils/matrix_utils.py +++ b/discretize/utils/matrix_utils.py @@ -35,8 +35,9 @@ def mkvc(x, n_dims=1, **kwargs): >>> from discretize.utils import mkvc >>> import numpy as np + >>> rng = np.random.default_rng(856) - >>> a = np.random.rand(3, 2) + >>> a = rng.random(3, 2) >>> a array([[0.33534155, 0.25334363], [0.07147884, 0.81080958], @@ -570,7 +571,8 @@ def get_subarray(A, ind): >>> from discretize.utils import get_subarray >>> import numpy as np - >>> A = np.random.rand(3, 3) + >>> rng = np.random.default_rng(421) + >>> A = rng.random((3, 3)) >>> A array([[1.07969034e-04, 9.78613931e-01, 6.62123429e-01], [8.80722877e-01, 7.61035691e-01, 7.42546796e-01], @@ -1167,6 +1169,7 @@ def make_property_tensor(mesh, tensor): >>> import numpy as np >>> import matplotlib.pyplot as plt >>> import matplotlib as mpl + >>> rng = np.random.default_rng(421) Define a 2D tensor mesh @@ -1176,9 +1179,9 @@ def make_property_tensor(mesh, tensor): Define a physical property for all cases (2D) >>> sigma_scalar = 4. - >>> sigma_isotropic = np.random.randint(1, 10, mesh.nC) - >>> sigma_anisotropic = np.random.randint(1, 10, (mesh.nC, 2)) - >>> sigma_tensor = np.random.randint(1, 10, (mesh.nC, 3)) + >>> sigma_isotropic = rng.integers(1, 10, mesh.nC) + >>> sigma_anisotropic = rng.integers(1, 10, (mesh.nC, 2)) + >>> sigma_tensor = rng.integers(1, 10, (mesh.nC, 3)) Construct the property tensor in each case @@ -1319,6 +1322,7 @@ def inverse_property_tensor(mesh, tensor, return_matrix=False, **kwargs): >>> import numpy as np >>> import matplotlib.pyplot as plt >>> import matplotlib as mpl + >>> rng = np.random.default_rng(421) Define a 2D tensor mesh @@ -1328,9 +1332,9 @@ def inverse_property_tensor(mesh, tensor, return_matrix=False, **kwargs): Define a physical property for all cases (2D) >>> sigma_scalar = 4. - >>> sigma_isotropic = np.random.randint(1, 10, mesh.nC) - >>> sigma_anisotropic = np.random.randint(1, 10, (mesh.nC, 2)) - >>> sigma_tensor = np.random.randint(1, 10, (mesh.nC, 3)) + >>> sigma_isotropic = rng.integers(1, 10, mesh.nC) + >>> sigma_anisotropic = rng.integers(1, 10, (mesh.nC, 2)) + >>> sigma_tensor = rng.integers(1, 10, (mesh.nC, 3)) Construct the property tensor in each case diff --git a/discretize/utils/mesh_utils.py b/discretize/utils/mesh_utils.py index d5a7b500d..22ea4fac7 100644 --- a/discretize/utils/mesh_utils.py +++ b/discretize/utils/mesh_utils.py @@ -428,8 +428,9 @@ def mesh_builder_xyz( >>> import discretize >>> import matplotlib.pyplot as plt >>> import numpy as np + >>> rng = np.random.default_rng(87142) - >>> xy_loc = np.random.randn(8,2) + >>> xy_loc = rng.standard_normal((8,2)) >>> mesh = discretize.utils.mesh_builder_xyz( ... xy_loc, [0.1, 0.1], depth_core=0.5, ... padding_distance=[[1,2], [1,0]], diff --git a/examples/plot_cahn_hilliard.py b/examples/plot_cahn_hilliard.py index 6ccfd2c2f..2784df405 100644 --- a/examples/plot_cahn_hilliard.py +++ b/examples/plot_cahn_hilliard.py @@ -44,9 +44,9 @@ """ import discretize -from pymatsolver import Solver import numpy as np import matplotlib.pyplot as plt +from scipy.sparse.linalg import spsolve def run(plotIt=True, n=60): @@ -92,7 +92,7 @@ def run(plotIt=True, n=60): MAT = dt * A * d2fdphi2 - I - dt * A * L rhs = (dt * A * d2fdphi2 - I) * phi - dt * A * dfdphi - phi = Solver(MAT) * rhs + phi = spsolve(MAT, rhs) if elapsed > capture[jj]: PHIS += [(elapsed, phi.copy())] diff --git a/examples/plot_dc_resistivity.py b/examples/plot_dc_resistivity.py index d8ee4abda..4f1cd8eac 100644 --- a/examples/plot_dc_resistivity.py +++ b/examples/plot_dc_resistivity.py @@ -6,9 +6,9 @@ """ import discretize -from pymatsolver import SolverLU import numpy as np import matplotlib.pyplot as plt +from scipy.solver.linalg import spsolve def run(plotIt=True): @@ -37,12 +37,10 @@ def DCfun(mesh, pts): # Step3: Solve DC problem (LU solver) AtM, rhstM = DCfun(tM, pts) - AinvtM = SolverLU(AtM) - phitM = AinvtM * rhstM + phitM = spsolve(AtM, rhstM) ArM, rhsrM = DCfun(rM, pts) - AinvrM = SolverLU(ArM) - phirM = AinvrM * rhsrM + phirM = spsolve(ArM, rhsrM) if not plotIt: return diff --git a/tests/base/test_coordutils.py b/tests/base/test_coordutils.py index 97d565df3..efa15ba37 100644 --- a/tests/base/test_coordutils.py +++ b/tests/base/test_coordutils.py @@ -4,15 +4,16 @@ tol = 1e-15 +rng = np.random.default_rng(523) + class coorutilsTest(unittest.TestCase): def test_rotation_matrix_from_normals(self): - np.random.seed(0) - v0 = np.random.rand(3) + v0 = rng.random(3) v0 *= 1.0 / np.linalg.norm(v0) np.random.seed(5) - v1 = np.random.rand(3) + v1 = rng.random(3) v1 *= 1.0 / np.linalg.norm(v1) Rf = utils.rotation_matrix_from_normals(v0, v1) @@ -22,12 +23,11 @@ def test_rotation_matrix_from_normals(self): self.assertTrue(np.linalg.norm(utils.mkvc(Ri.dot(v1) - v0)) < tol) def test_rotate_points_from_normals(self): - np.random.seed(10) - v0 = np.random.rand(3) + v0 = rng.random(3) v0 *= 1.0 / np.linalg.norm(v0) np.random.seed(15) - v1 = np.random.rand(3) + v1 = rng.random(3) v1 *= 1.0 / np.linalg.norm(v1) v2 = utils.mkvc(utils.rotate_points_from_normals(utils.mkvc(v0, 2).T, v0, v1)) @@ -35,16 +35,15 @@ def test_rotate_points_from_normals(self): self.assertTrue(np.linalg.norm(v2 - v1) < tol) def test_rotateMatrixFromNormals(self): - np.random.seed(20) - n0 = np.random.rand(3) + n0 = rng.random(3) n0 *= 1.0 / np.linalg.norm(n0) np.random.seed(25) - n1 = np.random.rand(3) + n1 = rng.random(3) n1 *= 1.0 / np.linalg.norm(n1) np.random.seed(30) - scale = np.random.rand(100, 1) + scale = rng.random((100, 1)) XYZ0 = scale * n0 XYZ1 = scale * n1 diff --git a/tests/base/test_interpolation.py b/tests/base/test_interpolation.py index cea6abd02..bd46a2c93 100644 --- a/tests/base/test_interpolation.py +++ b/tests/base/test_interpolation.py @@ -3,8 +3,6 @@ import discretize -gen = np.random.default_rng(182) - MESHTYPES = ["uniformTensorMesh", "randomTensorMesh"] TOLERANCES = [0.9, 0.5, 0.5] call1 = lambda fun, xyz: fun(xyz) @@ -43,13 +41,13 @@ class TestInterpolation1D(discretize.tests.OrderTest): - LOCS = np.random.rand(50) * 0.6 + 0.2 name = "Interpolation 1D" meshTypes = MESHTYPES tolerance = TOLERANCES meshDimension = 1 meshSizes = [8, 16, 32, 64, 128] - rng = gen + rng = np.random.default_rng(5136) + LOCS = rng.random(50) * 0.6 + 0.2 def getError(self): funX = lambda x: np.cos(2 * np.pi * x) @@ -96,12 +94,12 @@ def test_outliers(self): class TestInterpolation2d(discretize.tests.OrderTest): name = "Interpolation 2D" - LOCS = np.random.rand(50, 2) * 0.6 + 0.2 meshTypes = MESHTYPES tolerance = TOLERANCES meshDimension = 2 meshSizes = [8, 16, 32, 64] - rng = gen + rng = np.random.default_rng(2457) + LOCS = rng.random((50, 2)) * 0.6 + 0.2 def getError(self): funX = lambda x, y: np.cos(2 * np.pi * y) @@ -180,14 +178,12 @@ def test_exceptions(self): class TestInterpolationSymCyl(discretize.tests.OrderTest): name = "Interpolation Symmetric 3D" - LOCS = np.c_[ - np.random.rand(4) * 0.6 + 0.2, np.zeros(4), np.random.rand(4) * 0.6 + 0.2 - ] meshTypes = ["uniform_symmetric_CylMesh"] # MESHTYPES + tolerance = 0.6 meshDimension = 3 meshSizes = [32, 64, 128, 256] - rng = gen + rng = np.random.default_rng(81756234) + LOCS = np.c_[rng.random(4) * 0.6 + 0.2, np.zeros(4), rng.random(4) * 0.6 + 0.2] def getError(self): funX = lambda x, y: np.cos(2 * np.pi * y) @@ -247,15 +243,15 @@ def test_orderEy(self): class TestInterpolationCyl(discretize.tests.OrderTest): name = "Interpolation Cylindrical 3D" - LOCS = np.c_[ - np.random.rand(20) * 0.6 + 0.2, - 2 * np.pi * (np.random.rand(20) * 0.6 + 0.2), - np.random.rand(20) * 0.6 + 0.2, - ] meshTypes = ["uniformCylMesh", "randomCylMesh"] # MESHTYPES + meshDimension = 3 meshSizes = [8, 16, 32, 64] - rng = gen + rng = np.random.default_rng(876234) + LOCS = np.c_[ + rng.random(20) * 0.6 + 0.2, + 2 * np.pi * (rng.random(20) * 0.6 + 0.2), + rng.random(20) * 0.6 + 0.2, + ] def getError(self): func = lambda x, y, z: np.cos(2 * np.pi * x) + np.cos(y) + np.cos(2 * np.pi * z) @@ -318,13 +314,13 @@ def test_orderEz(self): class TestInterpolation3D(discretize.tests.OrderTest): + rng = np.random.default_rng(234) name = "Interpolation" - LOCS = np.random.rand(50, 3) * 0.6 + 0.2 + LOCS = rng.random((50, 3)) * 0.6 + 0.2 meshTypes = MESHTYPES tolerance = TOLERANCES meshDimension = 3 meshSizes = [8, 16, 32, 64] - rng = gen def getError(self): funX = lambda x, y, z: np.cos(2 * np.pi * y) diff --git a/tests/base/test_operators.py b/tests/base/test_operators.py index 1679c2a25..b5fc225a3 100644 --- a/tests/base/test_operators.py +++ b/tests/base/test_operators.py @@ -45,7 +45,6 @@ class TestCurl(discretize.tests.OrderTest): name = "Curl" meshTypes = MESHTYPES - rng = gen def getError(self): # fun: i (cos(y)) + j (cos(z)) + k (cos(x)) @@ -83,7 +82,6 @@ class TestCurl2D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 2 meshSizes = [8, 16, 32, 64] - rng = gen def getError(self): # Test function @@ -110,7 +108,6 @@ def test_order(self): # 1 # because of the averaging involved in the ghost point. u_b = (u_n + u_g)/2 # ) # meshSizes = [8, 16, 32, 64] -# rng = gen # # def getError(self): # # Test function @@ -138,7 +135,6 @@ class TestCellGrad2D_Dirichlet(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 2 meshSizes = [8, 16, 32, 64] - rng = gen def getError(self): # Test function @@ -167,7 +163,6 @@ class TestCellGrad3D_Dirichlet(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 3 meshSizes = [8, 16, 32] - rng = gen def getError(self): # Test function @@ -219,7 +214,6 @@ class TestCellGrad2D_Neumann(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 2 meshSizes = [8, 16, 32, 64] - rng = gen def getError(self): # Test function @@ -248,7 +242,6 @@ class TestCellGrad3D_Neumann(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 3 meshSizes = [8, 16, 32] - rng = gen def getError(self): # Test function @@ -299,7 +292,6 @@ class TestFaceDiv3D(discretize.tests.OrderTest): name = "Face Divergence 3D" meshTypes = MESHTYPES meshSizes = [8, 16, 32, 64] - rng = gen def getError(self): # Test function @@ -335,7 +327,6 @@ class TestFaceDiv2D(discretize.tests.OrderTest): meshTypes = MESHTYPES meshDimension = 2 meshSizes = [8, 16, 32, 64] - rng = gen def getError(self): # Test function @@ -360,7 +351,6 @@ def test_order(self): class TestNodalGrad(discretize.tests.OrderTest): name = "Nodal Gradient" meshTypes = MESHTYPES - rng = gen def getError(self): # Test function @@ -388,7 +378,6 @@ class TestNodalGrad2D(discretize.tests.OrderTest): name = "Nodal Gradient 2D" meshTypes = MESHTYPES meshDimension = 2 - rng = gen def getError(self): # Test function @@ -416,7 +405,6 @@ class TestAveraging1D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 1 meshSizes = [16, 32, 64] - rng = gen def getError(self): num = self.getAve(self.M) * self.getHere(self.M) @@ -518,8 +506,8 @@ def test_orderE2FV(self): class TestAverating2DSimple(unittest.TestCase): def setUp(self): - hx = np.random.rand(10) - hy = np.random.rand(10) + hx = gen.random(10) + hy = gen.random(10) self.mesh = discretize.TensorMesh([hx, hy]) def test_constantEdges(self): @@ -538,7 +526,6 @@ class TestAveraging2D(discretize.tests.OrderTest): meshTypes = MESHTYPES meshDimension = 2 meshSizes = [16, 32, 64] - rng = gen def getError(self): num = self.getAve(self.M) * self.getHere(self.M) @@ -637,9 +624,9 @@ def test_orderCC2FV(self): class TestAverating3DSimple(unittest.TestCase): def setUp(self): - hx = np.random.rand(10) - hy = np.random.rand(10) - hz = np.random.rand(10) + hx = gen.random(10) + hy = gen.random(10) + hz = gen.random(10) self.mesh = discretize.TensorMesh([hx, hy, hz]) def test_constantEdges(self): @@ -658,7 +645,6 @@ class TestAveraging3D(discretize.tests.OrderTest): meshTypes = MESHTYPES meshDimension = 3 meshSizes = [16, 32, 64] - rng = gen def getError(self): num = self.getAve(self.M) * self.getHere(self.M) @@ -792,7 +778,7 @@ def test_DivCurl(self): mesh, _ = discretize.tests.setup_mesh( meshType, self.meshSize, self.meshDimension ) - v = np.random.rand(mesh.nE) + v = gen.random(mesh.nE) divcurlv = mesh.face_divergence * (mesh.edge_curl * v) rel_err = np.linalg.norm(divcurlv) / np.linalg.norm(v) passed = rel_err < self.tol @@ -806,7 +792,7 @@ def test_CurlGrad(self): mesh, _ = discretize.tests.setup_mesh( meshType, self.meshSize, self.meshDimension ) - v = np.random.rand(mesh.nN) + v = gen.random(mesh.nN) curlgradv = mesh.edge_curl * (mesh.nodal_gradient * v) rel_err = np.linalg.norm(curlgradv) / np.linalg.norm(v) passed = rel_err < self.tol diff --git a/tests/base/test_tensor.py b/tests/base/test_tensor.py index cf42f226f..5d880f5bc 100644 --- a/tests/base/test_tensor.py +++ b/tests/base/test_tensor.py @@ -2,7 +2,7 @@ import numpy as np import unittest import discretize -from pymatsolver import Solver +from scipy.sparse.linalg import spsolve TOL = 1e-10 @@ -299,7 +299,7 @@ def getError(self): err = np.linalg.norm((sA - sN), np.inf) else: fA = fun(self.M.gridCC) - fN = Solver(D * G) * (sol(self.M.gridCC)) + fN = spsolve(D * G, sol(self.M.gridCC)) err = np.linalg.norm((fA - fN), np.inf) return err diff --git a/tests/base/test_tensor_innerproduct.py b/tests/base/test_tensor_innerproduct.py index 8ca27543d..7f2ebd9db 100644 --- a/tests/base/test_tensor_innerproduct.py +++ b/tests/base/test_tensor_innerproduct.py @@ -4,8 +4,6 @@ from discretize import TensorMesh from discretize.utils import sdinv -gen = np.random.default_rng(50) - class TestInnerProducts(discretize.tests.OrderTest): """Integrate an function over a unit cube domain @@ -14,7 +12,6 @@ class TestInnerProducts(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh", "uniformCurv", "rotateCurv"] meshDimension = 3 meshSizes = [16, 32] - rng = gen def getError(self): call = lambda fun, xyz: fun(xyz[:, 0], xyz[:, 1], xyz[:, 2]) @@ -174,7 +171,6 @@ class TestInnerProductsFaceProperties3D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 3 meshSizes = [16, 32] - rng = gen def getError(self): call = lambda fun, xyz: fun(xyz[:, 0], xyz[:, 1], xyz[:, 2]) @@ -321,7 +317,6 @@ class TestInnerProductsEdgeProperties3D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 3 meshSizes = [16, 32] - rng = gen def getError(self): call = lambda fun, xyz: fun(xyz[:, 0], xyz[:, 1], xyz[:, 2]) @@ -410,7 +405,6 @@ class TestInnerProducts2D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh", "uniformCurv", "rotateCurv"] meshDimension = 2 meshSizes = [4, 8, 16, 32, 64, 128] - rng = gen def getError(self): z = 5 # Because 5 is just such a great number. @@ -556,7 +550,6 @@ class TestInnerProductsFaceProperties2D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 2 meshSizes = [8, 16, 32] - rng = gen def getError(self): call = lambda fun, xy: fun(xy[:, 0], xy[:, 1]) @@ -647,7 +640,6 @@ class TestInnerProductsEdgeProperties2D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 2 meshSizes = [8, 16, 32] - rng = gen def getError(self): call = lambda fun, xy: fun(xy[:, 0], xy[:, 1]) @@ -710,7 +702,6 @@ class TestInnerProducts1D(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh"] meshDimension = 1 meshSizes = [4, 8, 16, 32, 64, 128] - rng = gen def getError(self): y = 12 # Because 12 is just such a great number. diff --git a/tests/base/test_tensor_innerproduct_derivs.py b/tests/base/test_tensor_innerproduct_derivs.py index 0273192c6..68d1cf2b8 100644 --- a/tests/base/test_tensor_innerproduct_derivs.py +++ b/tests/base/test_tensor_innerproduct_derivs.py @@ -3,7 +3,7 @@ import discretize from discretize import TensorMesh -np.random.seed(50) +rng = np.random.default_rng(542) class TestInnerProductsDerivsTensor(unittest.TestCase): @@ -19,8 +19,8 @@ def doTestFace( mesh.number(balance=False) elif meshType == "Tensor": mesh = discretize.TensorMesh(h) - v = np.random.rand(mesh.nF) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nC * rep) + v = rng.random(mesh.nF) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nC * rep) def fun(sig): M = mesh.get_face_inner_product( @@ -56,8 +56,8 @@ def doTestEdge( mesh.number(balance=False) elif meshType == "Tensor": mesh = discretize.TensorMesh(h) - v = np.random.rand(mesh.nE) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nC * rep) + v = rng.random(mesh.nE) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nC * rep) def fun(sig): M = mesh.get_edge_inner_product( @@ -341,8 +341,8 @@ def doTestFace(self, h, rep, meshType, invert_model=False, invert_matrix=False): mesh.number(balance=False) elif meshType == "Tensor": mesh = discretize.TensorMesh(h) - v = np.random.rand(mesh.nF) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nF * rep) + v = rng.random(mesh.nF) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nF * rep) def fun(sig): M = mesh.get_face_inner_product_surface( @@ -372,8 +372,8 @@ def doTestEdge(self, h, rep, meshType, invert_model=False, invert_matrix=False): mesh.number(balance=False) elif meshType == "Tensor": mesh = discretize.TensorMesh(h) - v = np.random.rand(mesh.nE) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nF * rep) + v = rng.random(mesh.nE) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nF * rep) def fun(sig): M = mesh.get_edge_inner_product_surface( @@ -477,8 +477,8 @@ def doTestEdge(self, h, rep, meshType, invert_model=False, invert_matrix=False): mesh.number(balance=False) elif meshType == "Tensor": mesh = discretize.TensorMesh(h) - v = np.random.rand(mesh.nE) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nE * rep) + v = rng.random(mesh.nE) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nE * rep) def fun(sig): M = mesh.get_edge_inner_product_line( @@ -544,7 +544,7 @@ class TestTensorSizeErrorRaises(unittest.TestCase): def setUp(self): self.mesh3D = TensorMesh([4, 4, 4]) - self.model = np.random.rand(self.mesh3D.nC) + self.model = rng.random(self.mesh3D.nC) def test_edge_inner_product_surface_deriv(self): self.assertRaises( diff --git a/tests/base/test_tensor_omf.py b/tests/base/test_tensor_omf.py index 467cff20d..d4179d545 100644 --- a/tests/base/test_tensor_omf.py +++ b/tests/base/test_tensor_omf.py @@ -50,6 +50,7 @@ def test_to_omf(self): self.assertTrue(np.allclose(models[name], arr)) def test_from_omf(self): + rng = np.random.default_rng(52134) omf_element = omf.VolumeElement( name="vol_ir", geometry=omf.VolumeGridGeometry( @@ -65,7 +66,7 @@ def test_from_omf(self): omf.ScalarData( name="Random Data", location="cells", - array=np.random.rand(10, 15, 20).flatten(), + array=rng.random((10, 15, 20)).flatten(), ) ], ) diff --git a/tests/base/test_tests.py b/tests/base/test_tests.py index d2f279cf7..e8fc560c3 100644 --- a/tests/base/test_tests.py +++ b/tests/base/test_tests.py @@ -73,20 +73,23 @@ def test_simplePass(self): def simplePass(x): return np.sin(x), sp.diags(np.cos(x)) - check_derivative(simplePass, np.random.randn(5), plotIt=False) + rng = np.random.default_rng(5322) + check_derivative(simplePass, rng.standard_normal(5), plotIt=False, rng=rng) def test_simpleFunction(self): def simpleFunction(x): return np.sin(x), lambda xi: np.cos(x) * xi - check_derivative(simpleFunction, np.random.randn(5), plotIt=False) + rng = np.random.default_rng(5322) + check_derivative(simpleFunction, rng.standard_normal(5), plotIt=False, rng=rng) def test_simpleFail(self): def simpleFail(x): return np.sin(x), -sp.diags(np.cos(x)) + rng = np.random.default_rng(5322) with pytest.raises(AssertionError): - check_derivative(simpleFail, np.random.randn(5), plotIt=False) + check_derivative(simpleFail, rng.standard_normal(5), plotIt=False, rng=rng) @pytest.mark.parametrize("test_type", ["mean", "min", "last", "all", "mean_at_least"]) diff --git a/tests/base/test_utils.py b/tests/base/test_utils.py index d34303f13..bf4ddc5b8 100644 --- a/tests/base/test_utils.py +++ b/tests/base/test_utils.py @@ -123,7 +123,8 @@ def test_index_cube_3D(self): ) def test_invXXXBlockDiagonal(self): - a = [np.random.rand(5, 1) for i in range(4)] + rng = np.random.default_rng(78352) + a = [rng.random((5, 1)) for i in range(4)] B = inverse_2x2_block_diagonal(*a) @@ -137,7 +138,7 @@ def test_invXXXBlockDiagonal(self): Z2 = B * A - sp.identity(10) self.assertTrue(np.linalg.norm(Z2.todense().ravel(), 2) < TOL) - a = [np.random.rand(5, 1) for i in range(9)] + a = [rng.random((5, 1)) for i in range(9)] B = inverse_3x3_block_diagonal(*a) A = sp.vstack( @@ -153,10 +154,11 @@ def test_invXXXBlockDiagonal(self): self.assertTrue(np.linalg.norm(Z3.todense().ravel(), 2) < TOL) def test_inverse_property_tensor2D(self): + rng = np.random.default_rng(763) M = discretize.TensorMesh([6, 6]) - a1 = np.random.rand(M.nC) - a2 = np.random.rand(M.nC) - a3 = np.random.rand(M.nC) + a1 = rng.random(M.nC) + a2 = rng.random(M.nC) + a3 = rng.random(M.nC) prop1 = a1 prop2 = np.c_[a1, a2] prop3 = np.c_[a1, a2, a3] @@ -173,10 +175,11 @@ def test_inverse_property_tensor2D(self): self.assertTrue(np.linalg.norm(Z.todense().ravel(), 2) < TOL) def test_TensorType2D(self): + rng = np.random.default_rng(8546) M = discretize.TensorMesh([6, 6]) - a1 = np.random.rand(M.nC) - a2 = np.random.rand(M.nC) - a3 = np.random.rand(M.nC) + a1 = rng.random(M.nC) + a2 = rng.random(M.nC) + a3 = rng.random(M.nC) prop1 = a1 prop2 = np.c_[a1, a2] prop3 = np.c_[a1, a2, a3] @@ -188,13 +191,14 @@ def test_TensorType2D(self): self.assertTrue(TensorType(M, None) == -1) def test_TensorType3D(self): + rng = np.random.default_rng(78352) M = discretize.TensorMesh([6, 6, 7]) - a1 = np.random.rand(M.nC) - a2 = np.random.rand(M.nC) - a3 = np.random.rand(M.nC) - a4 = np.random.rand(M.nC) - a5 = np.random.rand(M.nC) - a6 = np.random.rand(M.nC) + a1 = rng.random(M.nC) + a2 = rng.random(M.nC) + a3 = rng.random(M.nC) + a4 = rng.random(M.nC) + a5 = rng.random(M.nC) + a6 = rng.random(M.nC) prop1 = a1 prop2 = np.c_[a1, a2, a3] prop3 = np.c_[a1, a2, a3, a4, a5, a6] @@ -206,13 +210,14 @@ def test_TensorType3D(self): self.assertTrue(TensorType(M, None) == -1) def test_inverse_property_tensor3D(self): + rng = np.random.default_rng(78352) M = discretize.TensorMesh([6, 6, 6]) - a1 = np.random.rand(M.nC) - a2 = np.random.rand(M.nC) - a3 = np.random.rand(M.nC) - a4 = np.random.rand(M.nC) - a5 = np.random.rand(M.nC) - a6 = np.random.rand(M.nC) + a1 = rng.random(M.nC) + a2 = rng.random(M.nC) + a3 = rng.random(M.nC) + a4 = rng.random(M.nC) + a5 = rng.random(M.nC) + a6 = rng.random(M.nC) prop1 = a1 prop2 = np.c_[a1, a2, a3] prop3 = np.c_[a1, a2, a3, a4, a5, a6] diff --git a/tests/base/test_view.py b/tests/base/test_view.py index bcab25849..64baf98df 100644 --- a/tests/base/test_view.py +++ b/tests/base/test_view.py @@ -6,8 +6,6 @@ import pytest -np.random.seed(16) - TOL = 1e-1 @@ -51,7 +49,7 @@ def test_incorrectAxesWarnings(self): def test_plot_image(self): with self.assertRaises(NotImplementedError): - self.mesh.plot_image(np.random.rand(self.mesh.nC)) + self.mesh.plot_image(np.empty(self.mesh.nC)) if __name__ == "__main__": diff --git a/tests/base/test_volume_avg.py b/tests/base/test_volume_avg.py index ed319bf9e..652c6a521 100644 --- a/tests/base/test_volume_avg.py +++ b/tests/base/test_volume_avg.py @@ -7,9 +7,10 @@ class TestVolumeAverage(unittest.TestCase): def test_tensor_to_tensor(self): - h1 = np.random.rand(16) + rng = np.random.default_rng(68723) + h1 = rng.random(16) h1 /= h1.sum() - h2 = np.random.rand(16) + h2 = rng.random(16) h2 /= h2.sum() h1s = [] @@ -21,7 +22,7 @@ def test_tensor_to_tensor(self): mesh1 = discretize.TensorMesh(h1s) mesh2 = discretize.TensorMesh(h2s) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) @@ -40,9 +41,10 @@ def test_tensor_to_tensor(self): self.assertAlmostEqual(vol1, vol2) def test_tree_to_tree(self): - h1 = np.random.rand(16) + rng = np.random.default_rng(68723) + h1 = rng.random(16) h1 /= h1.sum() - h2 = np.random.rand(16) + h2 = rng.random(16) h2 /= h2.sum() h1s = [h1] @@ -60,7 +62,7 @@ def test_tree_to_tree(self): mesh2 = discretize.TreeMesh(h2s) mesh2.insert_cells([insert_2], [4]) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) @@ -79,9 +81,10 @@ def test_tree_to_tree(self): self.assertAlmostEqual(vol1, vol2) def test_tree_to_tensor(self): - h1 = np.random.rand(16) + rng = np.random.default_rng(68723) + h1 = rng.random(16) h1 /= h1.sum() - h2 = np.random.rand(16) + h2 = rng.random(16) h2 /= h2.sum() h1s = [h1] @@ -96,7 +99,7 @@ def test_tree_to_tensor(self): mesh1.insert_cells([insert_1], [4]) mesh2 = discretize.TensorMesh(h2s) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) @@ -115,9 +118,10 @@ def test_tree_to_tensor(self): self.assertAlmostEqual(vol1, vol2) def test_tensor_to_tree(self): - h1 = np.random.rand(16) + rng = np.random.default_rng(68723) + h1 = rng.random(16) h1 /= h1.sum() - h2 = np.random.rand(16) + h2 = rng.random(16) h2 /= h2.sum() h1s = [h1] @@ -132,7 +136,7 @@ def test_tensor_to_tree(self): mesh2 = discretize.TreeMesh(h2s) mesh2.insert_cells([insert_2], [4]) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) @@ -151,9 +155,10 @@ def test_tensor_to_tree(self): self.assertAlmostEqual(vol1, vol2) def test_errors(self): - h1 = np.random.rand(16) + rng = np.random.default_rng(68723) + h1 = rng.random(16) h1 /= h1.sum() - h2 = np.random.rand(16) + h2 = rng.random(16) h2 /= h2.sum() mesh1D = discretize.TensorMesh([h1]) mesh2D = discretize.TensorMesh([h1, h1]) @@ -175,9 +180,9 @@ def test_errors(self): # Gives mismatching mesh dimensions volume_average(mesh2D, mesh3D) - model1 = np.random.randn(mesh2D.nC) - bad_model1 = np.random.randn(3) - bad_model2 = np.random.rand(1) + model1 = rng.standard_normal(mesh2D.nC) + bad_model1 = rng.standard_normal(3) + bad_model2 = rng.random(1) # gives input values with incorrect lengths with self.assertRaises(ValueError): volume_average(mesh2D, mesh2, bad_model1) @@ -185,7 +190,8 @@ def test_errors(self): volume_average(mesh2D, mesh2, model1, bad_model2) def test_tree_to_tree_same_base(self): - h1 = np.random.rand(16) + rng = np.random.default_rng(68723) + h1 = rng.random(16) h1 /= h1.sum() h1s = [h1] @@ -201,7 +207,7 @@ def test_tree_to_tree_same_base(self): mesh2 = discretize.TreeMesh(h1s) mesh2.insert_cells([insert_2], [4]) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) @@ -220,7 +226,8 @@ def test_tree_to_tree_same_base(self): self.assertAlmostEqual(vol1, vol2) def test_tree_to_tensor_same_base(self): - h1 = np.random.rand(16) + rng = np.random.default_rng(867532) + h1 = rng.random(16) h1 /= h1.sum() h1s = [h1] @@ -233,7 +240,7 @@ def test_tree_to_tensor_same_base(self): mesh1.insert_cells([insert_1], [4]) mesh2 = discretize.TensorMesh(h1s) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) @@ -252,7 +259,8 @@ def test_tree_to_tensor_same_base(self): self.assertAlmostEqual(vol1, vol2) def test_tensor_to_tree_same_base(self): - h1 = np.random.rand(16) + rng = np.random.default_rng(91) + h1 = rng.random(16) h1 /= h1.sum() h1s = [h1] @@ -265,7 +273,7 @@ def test_tensor_to_tree_same_base(self): mesh2 = discretize.TreeMesh(h1s) mesh2.insert_cells([insert_2], [4]) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) @@ -284,6 +292,7 @@ def test_tensor_to_tree_same_base(self): self.assertAlmostEqual(vol1, vol2) def test_tensor_to_tensor_sub(self): + rng = np.random.default_rng(9867153) h1 = np.ones(32) h2 = np.ones(16) @@ -296,7 +305,7 @@ def test_tensor_to_tensor_sub(self): mesh1 = discretize.TensorMesh(h1s) mesh2 = discretize.TensorMesh(h2s) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) @@ -320,6 +329,7 @@ def test_tensor_to_tensor_sub(self): self.assertAlmostEqual(vol1, vol2) def test_tree_to_tree_sub(self): + rng = np.random.default_rng(987263) h1 = np.ones(32) h2 = np.ones(16) @@ -338,7 +348,7 @@ def test_tree_to_tree_sub(self): mesh2 = discretize.TreeMesh(h2s) mesh2.insert_cells([insert_2], [4]) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) @@ -362,6 +372,7 @@ def test_tree_to_tree_sub(self): self.assertAlmostEqual(vol1, vol2) def test_tree_to_tensor_sub(self): + rng = np.random.default_rng(5) h1 = np.ones(32) h2 = np.ones(16) @@ -377,7 +388,7 @@ def test_tree_to_tensor_sub(self): mesh1.insert_cells([insert_1], [4]) mesh2 = discretize.TensorMesh(h2s) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) @@ -401,6 +412,7 @@ def test_tree_to_tensor_sub(self): self.assertAlmostEqual(vol1, vol2) def test_tensor_to_tree_sub(self): + rng = np.random.default_rng(1) h1 = np.ones(32) h2 = np.ones(16) @@ -416,7 +428,7 @@ def test_tensor_to_tree_sub(self): mesh2 = discretize.TreeMesh(h2s) mesh2.insert_cells([insert_2], [4]) - in_put = np.random.rand(mesh1.nC) + in_put = rng.random(mesh1.nC) out_put = np.empty(mesh2.nC) # test the three ways of calling... out1 = volume_average(mesh1, mesh2, in_put, out_put) diff --git a/tests/boundaries/test_boundary_integrals.py b/tests/boundaries/test_boundary_integrals.py index 6574b55dd..a9a47e7bf 100644 --- a/tests/boundaries/test_boundary_integrals.py +++ b/tests/boundaries/test_boundary_integrals.py @@ -3,8 +3,6 @@ import discretize from discretize.utils import cart2cyl, cyl2cart -gen = np.random.default_rng(2552) - def u(*args): if len(args) == 1: @@ -88,7 +86,6 @@ class Test1DBoundaryIntegral(discretize.tests.OrderTest): meshDimension = 1 expectedOrders = 2 meshSizes = [4, 8, 16, 32, 64, 128] - rng = gen def getError(self): mesh = self.M @@ -139,7 +136,6 @@ class Test2DBoundaryIntegral(discretize.tests.OrderTest): meshDimension = 2 expectedOrders = [2, 2, 2, 2, 1] meshSizes = [4, 8, 16, 32, 64, 128] - rng = gen def getError(self): mesh = self.M @@ -220,7 +216,7 @@ class Test3DBoundaryIntegral(discretize.tests.OrderTest): meshDimension = 3 expectedOrders = [2, 1, 2, 2, 2, 2] meshSizes = [4, 8, 16, 32] - rng = gen + rng = np.random.default_rng(57681234) def getError(self): mesh = self.M diff --git a/tests/boundaries/test_boundary_maxwell.py b/tests/boundaries/test_boundary_maxwell.py index 1ffd35888..ea19ec063 100644 --- a/tests/boundaries/test_boundary_maxwell.py +++ b/tests/boundaries/test_boundary_maxwell.py @@ -1,7 +1,7 @@ import numpy as np import discretize from discretize import utils -from pymatsolver import Pardiso +from scipy.sparse.linalg import spsolve class TestFz2D_InhomogeneousDirichlet(discretize.tests.OrderTest): @@ -33,7 +33,7 @@ def getError(self): A = V @ C @ MeI @ C.T @ V + V rhs = V @ q_ana + V @ C @ MeI @ M_be @ ez_bc - ez_test = Pardiso(A.tocsr()) * rhs + ez_test = spsolve(A, rhs) if self._meshType == "rotateCurv": err = np.linalg.norm(mesh.cell_volumes * (ez_test - ez_ana)) else: @@ -52,7 +52,7 @@ class TestE3D_Inhomogeneous(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh", "uniformTree", "rotateCurv"] meshDimension = 3 expectedOrders = [2, 2, 1] - meshSizes = [8, 16, 32] + meshSizes = [4, 8, 16] def getError(self): # Test function @@ -99,7 +99,7 @@ def q_fun(x): A = C.T @ Mf @ C + Me rhs = Me @ q_ana + M_be * h_bc - e_test = Pardiso(A.tocsr(), is_symmetric=True, is_positive_definite=True) * rhs + e_test = spsolve(A, rhs) diff = e_test - e_ana if "Face" in self.myTest: diff --git a/tests/boundaries/test_boundary_poisson.py b/tests/boundaries/test_boundary_poisson.py index 0a199d1d6..3c8a0eaa7 100644 --- a/tests/boundaries/test_boundary_poisson.py +++ b/tests/boundaries/test_boundary_poisson.py @@ -4,9 +4,7 @@ import unittest import discretize from discretize import utils -from pymatsolver import Solver, Pardiso - -gen = np.random.default_rng(42) +from scipy.sparse.linalg import spsolve class TestCC1D_InhomogeneousDirichlet(discretize.tests.OrderTest): @@ -15,7 +13,6 @@ class TestCC1D_InhomogeneousDirichlet(discretize.tests.OrderTest): meshDimension = 1 expectedOrders = 2 meshSizes = [4, 8, 16, 32, 64, 128] - rng = gen def getError(self): # Test function @@ -42,7 +39,7 @@ def getError(self): A = V @ D @ MfI @ G rhs = V @ q_ana - V @ D @ MfI @ M_bf @ phi_bc - phi_test = Solver(A) * rhs + phi_test = spsolve(A, rhs) err = np.linalg.norm((phi_test - phi_ana)) / np.sqrt(mesh.n_cells) return err @@ -59,7 +56,6 @@ class TestCC2D_InhomogeneousDirichlet(discretize.tests.OrderTest): meshDimension = 2 expectedOrders = [2, 2, 1] meshSizes = [4, 8, 16, 32, 64] - rng = gen def getError(self): # Test function @@ -83,7 +79,7 @@ def getError(self): A = V @ D @ MfI @ G rhs = V @ q_ana - V @ D @ MfI @ M_bf @ phi_bc - phi_test = Solver(A) * rhs + phi_test = spsolve(A, rhs) if self._meshType == "rotateCurv": err = np.linalg.norm(mesh.cell_volumes * (phi_test - phi_ana)) else: @@ -103,7 +99,6 @@ class TestCC1D_InhomogeneousNeumann(discretize.tests.OrderTest): meshDimension = 1 expectedOrders = 2 meshSizes = [4, 8, 16, 32, 64, 128] - rng = gen def getError(self): # Test function @@ -171,7 +166,6 @@ class TestCC2D_InhomogeneousNeumann(discretize.tests.OrderTest): expectedOrders = [2, 2, 1] meshSizes = [4, 8, 16, 32] # meshSizes = [4] - rng = gen def getError(self): # Test function @@ -234,7 +228,6 @@ class TestCC1D_InhomogeneousMixed(discretize.tests.OrderTest): meshDimension = 1 expectedOrders = 2 meshSizes = [4, 8, 16, 32, 64, 128] - rng = gen def getError(self): # Test function @@ -268,10 +261,10 @@ def getError(self): rhs = V @ q_ana - V @ D @ MfI @ b_bc if self.myTest == "xc": - xc = Solver(A) * rhs + xc = spsolve(A, rhs) err = np.linalg.norm(xc - xc_ana) / np.sqrt(mesh.n_cells) elif self.myTest == "xcJ": - xc = Solver(A) * rhs + xc = spsolve(A, rhs) j = MfI @ ((-G + B_bc) @ xc + b_bc) err = np.linalg.norm(j - j_ana, np.inf) return err @@ -293,7 +286,6 @@ class TestCC2D_InhomogeneousMixed(discretize.tests.OrderTest): meshDimension = 2 expectedOrders = [2, 2, 1] meshSizes = [2, 4, 8, 16] - rng = gen # meshSizes = [4] def getError(self): @@ -340,7 +332,7 @@ def getError(self): A = V @ D @ MfI @ (-G + B_bc) rhs = V @ q_ana - V @ D @ MfI @ b_bc - phi_test = Solver(A) * rhs + phi_test = spsolve(A, rhs) if self._meshType == "rotateCurv": err = np.linalg.norm(mesh.cell_volumes * (phi_test - phi_ana)) @@ -360,8 +352,7 @@ class TestCC3D_InhomogeneousMixed(discretize.tests.OrderTest): meshTypes = ["uniformTensorMesh", "uniformTree", "rotateCurv"] meshDimension = 3 expectedOrders = [2, 2, 2] - meshSizes = [2, 4, 8, 16, 32] - rng = gen + meshSizes = [2, 4, 8, 16] def getError(self): # Test function @@ -430,7 +421,7 @@ def getError(self): A = V @ D @ MfI @ (-G + B_bc) rhs = V @ q_ana - V @ D @ MfI @ b_bc - phi_test = Pardiso(A.tocsr()) * rhs + phi_test = spsolve(A, rhs) if self._meshType == "rotateCurv": err = np.linalg.norm(mesh.cell_volumes * (phi_test - phi_ana)) @@ -450,7 +441,6 @@ class TestN1D_boundaries(discretize.tests.OrderTest): meshDimension = 1 expectedOrders = 2 meshSizes = [2, 4, 8, 16, 32, 64, 128] - rng = gen # meshSizes = [4] def getError(self): @@ -498,7 +488,7 @@ def getError(self): rhs = P_f.T @ rhs - (P_f.T @ A @ P_b) @ phi_ana[[0]] A = P_f.T @ A @ P_f - phi_test = Solver(A) * rhs + phi_test = spsolve(A, rhs) if self.boundary_type == "Nuemann": phi_test = P_f @ phi_test + P_b @ phi_ana[[0]] @@ -530,7 +520,6 @@ class TestN2D_boundaries(discretize.tests.OrderTest): expectedOrders = 2 tolerance = [0.8, 0.8, 0.6] meshSizes = [8, 16, 32, 64] - rng = gen # meshSizes = [4] def getError(self): @@ -601,7 +590,7 @@ def getError(self): rhs = P_f.T @ rhs - (P_f.T @ A @ P_b) @ phi_ana[[0]] A = P_f.T @ A @ P_f - phi_test = Solver(A) * rhs + phi_test = spsolve(A, rhs) if self.boundary_type == "Nuemann": phi_test = P_f @ phi_test + P_b @ phi_ana[[0]] @@ -631,8 +620,7 @@ class TestN3D_boundaries(discretize.tests.OrderTest): meshDimension = 3 expectedOrders = 2 tolerance = 0.6 - meshSizes = [2, 4, 8, 16, 32] - rng = gen + meshSizes = [2, 4, 8, 16] # meshSizes = [4] def getError(self): @@ -727,7 +715,7 @@ def getError(self): rhs = P_f.T @ rhs - (P_f.T @ A @ P_b) @ phi_ana[[0]] A = P_f.T @ A @ P_f - phi_test = Pardiso(A.tocsr()) * rhs + phi_test = spsolve(A, rhs) if self.boundary_type == "Nuemann": phi_test = P_f @ phi_test + P_b @ phi_ana[[0]] diff --git a/tests/boundaries/test_errors.py b/tests/boundaries/test_errors.py index b8475acda..41430b0a1 100644 --- a/tests/boundaries/test_errors.py +++ b/tests/boundaries/test_errors.py @@ -3,6 +3,9 @@ import discretize +rng = np.random.default_rng(53679) + + class RobinOperatorTest(unittest.TestCase): def setUp(self): self.mesh = discretize.TensorMesh([18, 20, 32]) @@ -28,7 +31,7 @@ def testCellGradBroadcasting(self): np.testing.assert_equal(bt, b1) self.assertEqual((B1 - B_t).nnz, 0) - gamma = np.random.rand(n_boundary_faces, 2) + gamma = rng.random((n_boundary_faces, 2)) B1, b1 = mesh.cell_gradient_weak_form_robin( alpha=0.5, beta=1.5, gamma=gamma[:, 0] ) @@ -78,7 +81,7 @@ def testEdgeDivBroadcasting(self): np.testing.assert_allclose(bt, b1) np.testing.assert_allclose(B1.data, B_t.data) - gamma = np.random.rand(n_boundary_faces, 2) + gamma = rng.random((n_boundary_faces, 2)) B1, b1 = mesh.edge_divergence_weak_form_robin( alpha=0.5, beta=1.5, gamma=gamma[:, 0] ) @@ -90,7 +93,7 @@ def testEdgeDivBroadcasting(self): np.testing.assert_allclose(B1.data, B3.data) np.testing.assert_allclose(np.c_[b1, b2], b3) - gamma = np.random.rand(n_boundary_nodes, 2) + gamma = rng.random((n_boundary_nodes, 2)) B1, b1 = mesh.edge_divergence_weak_form_robin( alpha=0.5, beta=1.5, gamma=gamma[:, 0] ) diff --git a/tests/boundaries/test_tensor_boundary.py b/tests/boundaries/test_tensor_boundary.py index c32a55bac..94be4aa0a 100644 --- a/tests/boundaries/test_tensor_boundary.py +++ b/tests/boundaries/test_tensor_boundary.py @@ -1,7 +1,7 @@ import numpy as np import unittest import discretize -from pymatsolver import Solver +from scipy.sparse.linalg import spsolve MESHTYPES = ["uniformTensorMesh"] @@ -225,8 +225,7 @@ def q_fun(x): if self.myTest == "xc": # TODO: fix the null space - Ainv = Solver(A) - xc = Ainv * rhs + xc = spsolve(A, rhs) err = np.linalg.norm((xc - xc_ana), np.inf) else: NotImplementedError @@ -320,8 +319,7 @@ def gamma_fun(alpha, beta, phi, phi_deriv): A = Div * MfrhoI * G if self.myTest == "xc": - Ainv = Solver(A) - xc = Ainv * rhs + xc = spsolve(A, rhs) err = np.linalg.norm((xc - xc_ana), np.inf) else: NotImplementedError @@ -452,8 +450,7 @@ def gamma_fun(alpha, beta, phi, phi_deriv): if self.myTest == "xc": # TODO: fix the null space - Ainv = Solver(A) - xc = Ainv * rhs + xc = spsolve(A, rhs) err = np.linalg.norm((xc - xc_ana), np.inf) else: NotImplementedError diff --git a/tests/boundaries/test_tensor_boundary_poisson.py b/tests/boundaries/test_tensor_boundary_poisson.py index df3763117..5f5b191d6 100644 --- a/tests/boundaries/test_tensor_boundary_poisson.py +++ b/tests/boundaries/test_tensor_boundary_poisson.py @@ -3,7 +3,7 @@ import unittest import discretize from discretize import utils -from pymatsolver import Solver, SolverCG +from scipy.sparse.linalg import spsolve MESHTYPES = ["uniformTensorMesh"] @@ -55,13 +55,12 @@ def getError(self): err = np.linalg.norm((q - q_ana), np.inf) elif self.myTest == "xc": # TODO: fix the null space - solver = SolverCG(A, maxiter=1000) - xc = solver * rhs + xc = spsolve(A, rhs) print("ACCURACY", np.linalg.norm(utils.mkvc(A * xc) - rhs)) err = np.linalg.norm((xc - xc_ana), np.inf) elif self.myTest == "xcJ": # TODO: fix the null space - xc = Solver(A) * rhs + xc = spsolve(A, rhs) print(np.linalg.norm(utils.mkvc(A * xc) - rhs)) j = McI * (G * xc + P * phi_bc) err = np.linalg.norm((j - j_ana), np.inf) @@ -141,10 +140,10 @@ def getError(self): elif self.myTest == "q": err = np.linalg.norm((q - q_ana), np.inf) elif self.myTest == "xc": - xc = Solver(A) * (rhs) + xc = spsolve(A, rhs) err = np.linalg.norm((xc - xc_ana), np.inf) elif self.myTest == "xcJ": - xc = Solver(A) * (rhs) + xc = spsolve(A, rhs) j = McI * (G * xc + P * bc) err = np.linalg.norm((j - j_ana), np.inf) diff --git a/tests/cyl/test_cyl.py b/tests/cyl/test_cyl.py index c3bdd22af..1b239bf7a 100644 --- a/tests/cyl/test_cyl.py +++ b/tests/cyl/test_cyl.py @@ -5,7 +5,7 @@ import discretize from discretize import tests, utils -np.random.seed(13) +rng = np.random.default_rng(87564123) class TestCylSymmetricMesh(unittest.TestCase): @@ -388,8 +388,8 @@ class TestCellGrad2D_Dirichlet(unittest.TestCase): # self.orderTest() def setUp(self): - hx = np.random.rand(10) - hz = np.random.rand(10) + hx = rng.random(10) + hz = rng.random(10) self.mesh = discretize.CylindricalMesh([hx, 1, hz]) def test_NotImplementedError(self): @@ -399,8 +399,8 @@ def test_NotImplementedError(self): class TestAveragingSimple(unittest.TestCase): def setUp(self): - hx = np.random.rand(10) - hz = np.random.rand(10) + hx = rng.random(10) + hz = rng.random(10) self.mesh = discretize.CylindricalMesh([hx, 1, hz]) def test_simpleEdges(self): diff --git a/tests/cyl/test_cyl3D.py b/tests/cyl/test_cyl3D.py index 94d9ed6c6..4a7110046 100644 --- a/tests/cyl/test_cyl3D.py +++ b/tests/cyl/test_cyl3D.py @@ -4,8 +4,6 @@ import discretize from discretize import utils -rng = np.random.default_rng(16) - TOL = 1e-1 diff --git a/tests/cyl/test_cylOperators.py b/tests/cyl/test_cylOperators.py index 1dca64497..4486fca19 100644 --- a/tests/cyl/test_cylOperators.py +++ b/tests/cyl/test_cylOperators.py @@ -6,15 +6,13 @@ import discretize from discretize import tests -np.random.seed(16) - TOL = 1e-1 # ----------------------------- Test Operators ------------------------------ # -MESHTYPES = ["uniformCylMesh", "randomCylMesh"] +MESHTYPES = ["uniformCylMesh"] call2 = lambda fun, xyz: fun(xyz[:, 0], xyz[:, 2]) call3 = lambda fun, xyz: fun(xyz[:, 0], xyz[:, 1], xyz[:, 2]) cyl_row2 = lambda g, xfun, yfun: np.c_[call2(xfun, g), call2(yfun, g)] diff --git a/tests/cyl/test_cyl_innerproducts.py b/tests/cyl/test_cyl_innerproducts.py index b9752bfcd..d74859485 100644 --- a/tests/cyl/test_cyl_innerproducts.py +++ b/tests/cyl/test_cyl_innerproducts.py @@ -9,7 +9,7 @@ TOL = 1e-1 TOLD = 0.7 # tolerance on deriv checks -np.random.seed(99) +rng = np.random.default_rng(99) class FaceInnerProductFctsIsotropic(object): @@ -465,8 +465,8 @@ class TestCylInnerProducts_Deriv(unittest.TestCase): def setUp(self): n = 2 self.mesh = discretize.CylindricalMesh([n, 1, n]) - self.face_vec = np.random.rand(self.mesh.nF) - self.edge_vec = np.random.rand(self.mesh.nE) + self.face_vec = rng.random(self.mesh.nF) + self.edge_vec = rng.random(self.mesh.nE) # make up a smooth function self.x0 = 2 * self.mesh.gridCC[:, 0] ** 2 + self.mesh.gridCC[:, 2] ** 4 @@ -579,8 +579,8 @@ class TestCylInnerProductsAnisotropic_Deriv(unittest.TestCase): def setUp(self): n = 60 self.mesh = discretize.CylindricalMesh([n, 1, n]) - self.face_vec = np.random.rand(self.mesh.nF) - self.edge_vec = np.random.rand(self.mesh.nE) + self.face_vec = rng.random(self.mesh.nF) + self.edge_vec = rng.random(self.mesh.nE) # make up a smooth function self.x0 = np.array( [2 * self.mesh.gridCC[:, 0] ** 2 + self.mesh.gridCC[:, 2] ** 4] @@ -745,8 +745,8 @@ class TestCylInnerProductsFaceProperties_Deriv(unittest.TestCase): def setUp(self): n = 2 self.mesh = discretize.CylindricalMesh([n, 1, n]) - self.face_vec = np.random.rand(self.mesh.nF) - self.edge_vec = np.random.rand(self.mesh.nE) + self.face_vec = rng.random(self.mesh.nF) + self.edge_vec = rng.random(self.mesh.nE) # make up a smooth function self.x0 = np.r_[ 2 * self.mesh.gridFx[:, 0] ** 2 + self.mesh.gridFx[:, 2] ** 4, diff --git a/tests/cyl/test_cyl_operators.py b/tests/cyl/test_cyl_operators.py index e7e78c807..fa43d6b33 100644 --- a/tests/cyl/test_cyl_operators.py +++ b/tests/cyl/test_cyl_operators.py @@ -13,6 +13,8 @@ variable_names=list("RTZ"), ) +rng = np.random.default_rng(563279) + def lambdify_vector(variabls, u_vecs, func): funcs = [sp.lambdify(variabls, func.coeff(u_hat), "numpy") for u_hat in u_vecs] @@ -346,7 +348,7 @@ def get_error(n_cells): def test_mimetic_div_curl(mesh_type): mesh, _ = setup_mesh(mesh_type, 10) - v = np.random.rand(mesh.n_edges) + v = rng.random(mesh.n_edges) divcurlv = mesh.face_divergence @ (mesh.edge_curl @ v) np.testing.assert_allclose(divcurlv, 0, atol=1e-11) @@ -355,7 +357,7 @@ def test_mimetic_div_curl(mesh_type): def test_mimetic_curl_grad(mesh_type): mesh, _ = setup_mesh(mesh_type, 10) - v = np.random.rand(mesh.n_nodes) + v = rng.random(mesh.n_nodes) divcurlv = mesh.edge_curl @ (mesh.nodal_gradient @ v) np.testing.assert_allclose(divcurlv, 0, atol=1e-11) diff --git a/tests/simplex/test_inner_products.py b/tests/simplex/test_inner_products.py index 88d454f8a..83f6875a4 100644 --- a/tests/simplex/test_inner_products.py +++ b/tests/simplex/test_inner_products.py @@ -4,6 +4,8 @@ import scipy.sparse as sp from discretize.utils import example_simplex_mesh +rng = np.random.default_rng(4421) + def u(*args): if len(args) == 1: @@ -332,8 +334,8 @@ class TestInnerProductsDerivs(unittest.TestCase): def doTestFace(self, h, rep): nodes, simplices = example_simplex_mesh(h) mesh = discretize.SimplexMesh(nodes, simplices) - v = np.random.rand(mesh.n_faces) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nC * rep) + v = rng.random(mesh.n_faces) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nC * rep) def fun(sig): M = mesh.get_face_inner_product(sig) @@ -346,8 +348,8 @@ def fun(sig): def doTestEdge(self, h, rep): nodes, simplices = example_simplex_mesh(h) mesh = discretize.SimplexMesh(nodes, simplices) - v = np.random.rand(mesh.n_edges) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nC * rep) + v = rng.random(mesh.n_edges) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nC * rep) def fun(sig): M = mesh.get_edge_inner_product(sig) @@ -544,7 +546,7 @@ def setUp(self): def test_bad_model_size(self): mesh = self.mesh - bad_model = np.random.rand(mesh.n_cells, 5) + bad_model = rng.random((mesh.n_cells, 5)) with self.assertRaises(ValueError): mesh.get_face_inner_product(bad_model) with self.assertRaises(ValueError): @@ -552,7 +554,7 @@ def test_bad_model_size(self): def test_cant_invert(self): mesh = self.mesh - good_model = np.random.rand(mesh.n_cells) + good_model = rng.random(mesh.n_cells) with self.assertRaises(NotImplementedError): mesh.get_face_inner_product(good_model, invert_matrix=True) with self.assertRaises(NotImplementedError): diff --git a/tests/simplex/test_utils.py b/tests/simplex/test_utils.py index 87138c1b9..e0c4fc1a4 100644 --- a/tests/simplex/test_utils.py +++ b/tests/simplex/test_utils.py @@ -13,6 +13,8 @@ except ImportError: has_vtk = False +rng = np.random.default_rng(87916253) + class SimplexTests(unittest.TestCase): def test_init_errors(self): @@ -42,7 +44,7 @@ def test_init_errors(self): with self.assertRaises(ValueError): # pass bad dimensionality - discretize.SimplexMesh(np.random.rand(10, 4), simplices[:, :-1]) + discretize.SimplexMesh(np.ones((10, 4)), simplices[:, :-1]) def test_find_containing(self): n = 4 @@ -81,11 +83,11 @@ def test_image_plotting(self): points, simplices = discretize.utils.example_simplex_mesh((n, n)) mesh = discretize.SimplexMesh(points, simplices) - cc_dat = np.random.rand(mesh.n_cells) - n_dat = np.random.rand(mesh.n_nodes) - f_dat = np.random.rand(mesh.n_faces) - e_dat = np.random.rand(mesh.n_edges) - ccv_dat = np.random.rand(mesh.n_cells, 2) + cc_dat = rng.random(mesh.n_cells) + n_dat = rng.random(mesh.n_nodes) + f_dat = rng.random(mesh.n_faces) + e_dat = rng.random(mesh.n_edges) + ccv_dat = rng.random((mesh.n_cells, 2)) mesh.plot_image(cc_dat) mesh.plot_image(ccv_dat, v_type="CCv", view="vec") @@ -106,7 +108,7 @@ def test_image_plotting(self): points, simplices = discretize.utils.example_simplex_mesh((n, n, n)) mesh = discretize.SimplexMesh(points, simplices) - cc_dat = np.random.rand(mesh.n_cells) + cc_dat = rng.random(mesh.n_cells) mesh.plot_image(cc_dat) plt.close("all") @@ -128,7 +130,7 @@ def test_2D_vtk(self): n = 5 points, simplices = discretize.utils.example_simplex_mesh((n, n)) mesh = discretize.SimplexMesh(points, simplices) - cc_dat = np.random.rand(mesh.n_cells) + cc_dat = rng.random(mesh.n_cells) vtk_obj = mesh.to_vtk(models={"info": cc_dat}) @@ -149,7 +151,7 @@ def test_3D_vtk(self): n = 5 points, simplices = discretize.utils.example_simplex_mesh((n, n, n)) mesh = discretize.SimplexMesh(points, simplices) - cc_dat = np.random.rand(mesh.n_cells) + cc_dat = rng.random(mesh.n_cells) vtk_obj = mesh.to_vtk(models={"info": cc_dat}) diff --git a/tests/tree/test_refine.py b/tests/tree/test_refine.py index e2ec53da4..6b9ae81cb 100644 --- a/tests/tree/test_refine.py +++ b/tests/tree/test_refine.py @@ -160,9 +160,9 @@ def refine_triangle(cell): def test_triangle_errors(): - not_triangles_array = np.random.rand(4, 2, 2) - triangle3 = np.random.rand(3, 3) - triangles2 = np.random.rand(10, 3, 2) + not_triangles_array = np.empty((4, 2, 2)) + triangle3 = np.empty((3, 3)) + triangles2 = np.empty((10, 3, 2)) levels = np.full(8, -1) mesh1 = discretize.TreeMesh([64, 64]) @@ -324,9 +324,9 @@ def refine_simplex(cell): def test_tetra_errors(): - not_simplex_array = np.random.rand(4, 3, 3) - simplex = np.random.rand(4, 2) - simplices = np.random.rand(10, 4, 3) + not_simplex_array = np.empty((4, 3, 3)) + simplex = np.empty((4, 2)) + simplices = np.empty((10, 4, 3)) levels = np.full(8, -1) mesh1 = discretize.TreeMesh([32, 32, 32]) @@ -461,24 +461,25 @@ def test_refine_triang_prism_errors(): mesh.refine_vertical_trianglular_prism(xyz[:, :-1], h, -1) # incorrect levels and triangles - ps = np.random.rand(10, 3, 3) + ps = np.empty((10, 3, 3)) with pytest.raises(ValueError): mesh.refine_vertical_trianglular_prism(ps, h, [-1, -2]) # incorrect heights and triangles - ps = np.random.rand(10, 3, 3) + ps = np.empty((10, 3, 3)) with pytest.raises(ValueError): mesh.refine_vertical_trianglular_prism(ps, [h, h], -1) # negative heights - ps = np.random.rand(10, 3, 3) + ps = np.empty((10, 3, 3)) with pytest.raises(ValueError): mesh.refine_vertical_trianglular_prism(ps, -h, -1) def test_bounding_box(): # No padding - xyz = np.random.rand(20, 2) * 0.25 + 3 / 8 + rng = np.random.default_rng(51623978) + xyz = rng.random((20, 2)) * 0.25 + 3 / 8 mesh1 = discretize.TreeMesh([32, 32]) mesh1.refine_bounding_box(xyz, -1, None) @@ -512,7 +513,7 @@ def test_bounding_box(): def test_bounding_box_errors(): mesh1 = discretize.TreeMesh([32, 32]) - xyz = np.random.rand(20, 3) + xyz = np.empty((20, 3)) # incorrect padding shape with pytest.raises(ValueError): mesh1.refine_bounding_box(xyz, -1, [[2, 3, 4]]) diff --git a/tests/tree/test_tree.py b/tests/tree/test_tree.py index d63310a3d..878a5a73f 100644 --- a/tests/tree/test_tree.py +++ b/tests/tree/test_tree.py @@ -5,12 +5,14 @@ TOL = 1e-8 +rng = np.random.default_rng(6234) + class TestSimpleQuadTree(unittest.TestCase): def test_counts(self): nc = 8 - h1 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 - h2 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 + h1 = rng.random(nc) * nc * 0.5 + nc * 0.5 + h2 = rng.random(nc) * nc * 0.5 + nc * 0.5 h = [hi / np.sum(hi) for hi in [h1, h2]] # normalize M = discretize.TreeMesh(h) points = np.array([[0.1, 0.1]]) @@ -132,9 +134,9 @@ def test_serialization(self): class TestOcTree(unittest.TestCase): def test_counts(self): nc = 8 - h1 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 - h2 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 - h3 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 + h1 = rng.random(nc) * nc * 0.5 + nc * 0.5 + h2 = rng.random(nc) * nc * 0.5 + nc * 0.5 + h3 = rng.random(nc) * nc * 0.5 + nc * 0.5 h = [hi / np.sum(hi) for hi in [h1, h2, h3]] # normalize M = discretize.TreeMesh(h, levels=3) points = np.array([[0.2, 0.1, 0.7], [0.8, 0.4, 0.2]]) @@ -311,8 +313,8 @@ def refinefcn(cell): def test_cell_nodes(self): # 2D nc = 8 - h1 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 - h2 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 + h1 = rng.random(nc) * nc * 0.5 + nc * 0.5 + h2 = rng.random(nc) * nc * 0.5 + nc * 0.5 h = [hi / np.sum(hi) for hi in [h1, h2]] # normalize M = discretize.TreeMesh(h) points = np.array([[0.2, 0.1], [0.8, 0.4]]) @@ -326,9 +328,9 @@ def test_cell_nodes(self): # 3D nc = 8 - h1 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 - h2 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 - h3 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 + h1 = rng.random(nc) * nc * 0.5 + nc * 0.5 + h2 = rng.random(nc) * nc * 0.5 + nc * 0.5 + h3 = rng.random(nc) * nc * 0.5 + nc * 0.5 h = [hi / np.sum(hi) for hi in [h1, h2, h3]] # normalize M = discretize.TreeMesh(h, levels=3) points = np.array([[0.2, 0.1, 0.7], [0.8, 0.4, 0.2]]) @@ -346,8 +348,8 @@ class TestTreeMeshNodes: def sample_mesh(self, request): """Return a sample TreeMesh""" nc = 8 - h1 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 - h2 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 + h1 = rng.random(nc) * nc * 0.5 + nc * 0.5 + h2 = rng.random(nc) * nc * 0.5 + nc * 0.5 if request.param == "2D": h = [hi / np.sum(hi) for hi in [h1, h2]] # normalize mesh = discretize.TreeMesh(h) @@ -355,7 +357,7 @@ def sample_mesh(self, request): levels = np.array([1, 2]) mesh.insert_cells(points, levels, finalize=True) else: - h3 = np.random.rand(nc) * nc * 0.5 + nc * 0.5 + h3 = rng.random(nc) * nc * 0.5 + nc * 0.5 h = [hi / np.sum(hi) for hi in [h1, h2, h3]] # normalize mesh = discretize.TreeMesh(h, levels=3) points = np.array([[0.2, 0.1, 0.7], [0.8, 0.4, 0.2]]) @@ -407,12 +409,12 @@ def function(cell): self.M = M def test_fx(self): - r = np.random.rand(self.M.nFx) + r = rng.random(self.M.nFx) P = self.M.get_interpolation_matrix(self.M.gridFx, "Fx") self.assertLess(np.abs(P[:, : self.M.nFx] * r - r).max(), TOL) def test_fy(self): - r = np.random.rand(self.M.nFy) + r = rng.random(self.M.nFy) P = self.M.get_interpolation_matrix(self.M.gridFy, "Fy") self.assertLess(np.abs(P[:, self.M.nFx :] * r - r).max(), TOL) @@ -437,36 +439,36 @@ def function(cell): self.M = M def test_Fx(self): - r = np.random.rand(self.M.nFx) + r = rng.random(self.M.nFx) P = self.M.get_interpolation_matrix(self.M.gridFx, "Fx") self.assertLess(np.abs(P[:, : self.M.nFx] * r - r).max(), TOL) def test_Fy(self): - r = np.random.rand(self.M.nFy) + r = rng.random(self.M.nFy) P = self.M.get_interpolation_matrix(self.M.gridFy, "Fy") self.assertLess( np.abs(P[:, self.M.nFx : (self.M.nFx + self.M.nFy)] * r - r).max(), TOL ) def test_Fz(self): - r = np.random.rand(self.M.nFz) + r = rng.random(self.M.nFz) P = self.M.get_interpolation_matrix(self.M.gridFz, "Fz") self.assertLess(np.abs(P[:, (self.M.nFx + self.M.nFy) :] * r - r).max(), TOL) def test_Ex(self): - r = np.random.rand(self.M.nEx) + r = rng.random(self.M.nEx) P = self.M.get_interpolation_matrix(self.M.gridEx, "Ex") self.assertLess(np.abs(P[:, : self.M.nEx] * r - r).max(), TOL) def test_Ey(self): - r = np.random.rand(self.M.nEy) + r = rng.random(self.M.nEy) P = self.M.get_interpolation_matrix(self.M.gridEy, "Ey") self.assertLess( np.abs(P[:, self.M.nEx : (self.M.nEx + self.M.nEy)] * r - r).max(), TOL ) def test_Ez(self): - r = np.random.rand(self.M.nEz) + r = rng.random(self.M.nEz) P = self.M.get_interpolation_matrix(self.M.gridEz, "Ez") self.assertLess(np.abs(P[:, (self.M.nEx + self.M.nEy) :] * r - r).max(), TOL) diff --git a/tests/tree/test_tree_innerproduct_derivs.py b/tests/tree/test_tree_innerproduct_derivs.py index 96885a749..4bd060cea 100644 --- a/tests/tree/test_tree_innerproduct_derivs.py +++ b/tests/tree/test_tree_innerproduct_derivs.py @@ -2,6 +2,8 @@ import unittest import discretize +rng = np.random.default_rng(678423) + class TestInnerProductsDerivsTensor(unittest.TestCase): def doTestFace( @@ -15,8 +17,8 @@ def doTestFace( mesh.refine(lambda xc: 3) elif meshType == "Tensor": mesh = discretize.TensorMesh(h) - v = np.random.rand(mesh.nF) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nC * rep) + v = rng.random(mesh.nF) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nC * rep) def fun(sig): M = mesh.get_face_inner_product( @@ -38,7 +40,7 @@ def fun(sig): fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=rng) def doTestEdge( self, h, rep, fast, meshType, invert_model=False, invert_matrix=False @@ -51,8 +53,8 @@ def doTestEdge( mesh.refine(lambda xc: 3) elif meshType == "Tensor": mesh = discretize.TensorMesh(h) - v = np.random.rand(mesh.nE) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nC * rep) + v = rng.random(mesh.nE) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nC * rep) def fun(sig): M = mesh.get_edge_inner_product( @@ -74,7 +76,7 @@ def fun(sig): fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=rng) def test_FaceIP_2D_float_Tree(self): self.assertTrue(self.doTestFace([8, 8], 0, False, "Tree")) @@ -169,8 +171,8 @@ def doTestFace(self, h, rep, meshType, invert_model=False, invert_matrix=False): mesh.refine(lambda xc: 3) elif meshType == "Tensor": mesh = discretize.TensorMesh(h) - v = np.random.rand(mesh.nF) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nF * rep) + v = rng.random(mesh.nF) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nF * rep) def fun(sig): M = mesh.get_face_inner_product_surface( @@ -191,7 +193,7 @@ def fun(sig): # fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=rng) def doTestEdge(self, h, rep, meshType, invert_model=False, invert_matrix=False): if meshType == "Curv": @@ -202,8 +204,8 @@ def doTestEdge(self, h, rep, meshType, invert_model=False, invert_matrix=False): mesh.refine(lambda xc: 3) elif meshType == "Tensor": mesh = discretize.TensorMesh(h) - v = np.random.rand(mesh.nE) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nF * rep) + v = rng.random(mesh.nE) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nF * rep) def fun(sig): M = mesh.get_edge_inner_product_surface( @@ -224,7 +226,7 @@ def fun(sig): # fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=rng) def test_FaceIP_2D_float_fast_Tree(self): self.assertTrue(self.doTestFace([8, 8], 0, "Tree")) @@ -261,8 +263,8 @@ def doTestEdge(self, h, rep, meshType, invert_model=False, invert_matrix=False): mesh.refine(lambda xc: 3) elif meshType == "Tensor": mesh = discretize.TensorMesh(h) - v = np.random.rand(mesh.nE) - sig = np.random.rand(1) if rep == 0 else np.random.rand(mesh.nE * rep) + v = rng.random(mesh.nE) + sig = rng.random(1) if rep == 0 else rng.random(mesh.nE * rep) def fun(sig): M = mesh.get_edge_inner_product_line( @@ -283,7 +285,7 @@ def fun(sig): # fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=rng) def test_EdgeIP_2D_float_fast_Tree(self): self.assertTrue(self.doTestEdge([8, 8], 0, "Tree")) diff --git a/tests/tree/test_tree_interpolation.py b/tests/tree/test_tree_interpolation.py index 2886a2970..d4688364b 100644 --- a/tests/tree/test_tree_interpolation.py +++ b/tests/tree/test_tree_interpolation.py @@ -42,7 +42,6 @@ class TestInterpolation2d(discretize.tests.OrderTest): """ name = "Interpolation 2D" - # LOCS = np.random.rand(50, 2)*0.6+0.2 # location_type = 'Ex' X, Y = np.mgrid[0:1:250j, 0:1:250j] LOCS = np.c_[X.reshape(-1), Y.reshape(-1)] @@ -129,7 +128,6 @@ def test_orderEy(self): class TestInterpolation3D(discretize.tests.OrderTest): name = "Interpolation" - # LOCS = np.random.rand(50, 3)*0.6+0.2 X, Y, Z = np.mgrid[0:1:50j, 0:1:50j, 0:1:50j] LOCS = np.c_[X.reshape(-1), Y.reshape(-1), Z.reshape(-1)] meshTypes = MESHTYPES diff --git a/tests/tree/test_tree_operators.py b/tests/tree/test_tree_operators.py index 8d6d94f4b..483385d70 100644 --- a/tests/tree/test_tree_operators.py +++ b/tests/tree/test_tree_operators.py @@ -38,6 +38,7 @@ class TestCellGrad2D(discretize.tests.OrderTest): meshSizes = [8, 16] # because of the averaging involved in the ghost point. u_b = (u_n + u_g)/2 expectedOrders = 1 + rng = np.random.default_rng(87964213) def getError(self): # Test function @@ -55,7 +56,6 @@ def getError(self): return err def test_order(self): - np.random.seed(7) self.orderTest() @@ -66,6 +66,7 @@ class TestCellGrad3D(discretize.tests.OrderTest): meshSizes = [8, 16] # because of the averaging involved in the ghost point. u_b = (u_n + u_g)/2 expectedOrders = 1 + rng = np.random.default_rng(6957832) def getError(self): # Test function @@ -105,7 +106,6 @@ def getError(self): return err def test_order(self): - np.random.seed(6) self.orderTest() @@ -114,6 +114,7 @@ class TestFaceDivxy2D(discretize.tests.OrderTest): meshTypes = MESHTYPES meshDimension = 2 meshSizes = [16, 32] + rng = np.random.default_rng(19647823) def getError(self): # Test function @@ -136,7 +137,6 @@ def getError(self): return err def test_order(self): - np.random.seed(4) self.orderTest() @@ -144,6 +144,7 @@ class TestFaceDiv3D(discretize.tests.OrderTest): name = "Face Divergence 3D" meshTypes = MESHTYPES meshSizes = [8, 16, 32] + rng = np.random.default_rng(81725364) def getError(self): fx = lambda x, y, z: np.sin(2 * np.pi * x) @@ -164,7 +165,6 @@ def getError(self): return np.linalg.norm((divF - divF_ana), np.inf) def test_order(self): - np.random.seed(7) self.orderTest() @@ -173,6 +173,7 @@ class TestFaceDivxyz3D(discretize.tests.OrderTest): meshTypes = MESHTYPES meshDimension = 3 meshSizes = [8, 16, 32] + rng = np.random.default_rng(6172824) def getError(self): # Test function @@ -202,13 +203,12 @@ def getError(self): return err def test_order(self): - np.random.seed(7) self.orderTest() class TestCurl(discretize.tests.OrderTest): name = "Curl" - meshTypes = ["notatreeTree", "uniformTree"] # , 'randomTree']#, 'uniformTree'] + meshTypes = ["notatreeTree", "uniformTree"] meshSizes = [8, 16] # , 32] expectedOrders = [2, 1] # This is due to linear interpolation in the Re projection @@ -238,13 +238,12 @@ def getError(self): return err def test_order(self): - np.random.seed(7) self.orderTest() class TestNodalGrad(discretize.tests.OrderTest): name = "Nodal Gradient" - meshTypes = ["notatreeTree", "uniformTree"] # ['randomTree', 'uniformTree'] + meshTypes = ["notatreeTree", "uniformTree"] meshSizes = [8, 16] # , 32] expectedOrders = [2, 1] @@ -267,13 +266,12 @@ def getError(self): return err def test_order(self): - np.random.seed(7) self.orderTest() class TestNodalGrad2D(discretize.tests.OrderTest): name = "Nodal Gradient 2D" - meshTypes = ["notatreeTree", "uniformTree"] # ['randomTree', 'uniformTree'] + meshTypes = ["notatreeTree", "uniformTree"] meshSizes = [8, 16] # , 32] expectedOrders = [2, 1] meshDimension = 2 @@ -296,7 +294,6 @@ def getError(self): return err def test_order(self): - np.random.seed(7) self.orderTest() @@ -918,7 +915,7 @@ def test_order1_edges_invert_model(self): class TestTreeAveraging2D(discretize.tests.OrderTest): """Integrate an function over a unit cube domain using edgeInnerProducts and faceInnerProducts.""" - meshTypes = ["notatreeTree", "uniformTree"] # 'randomTree'] + meshTypes = ["notatreeTree", "uniformTree"] meshDimension = 2 meshSizes = [4, 8, 16] expectedOrders = [2, 1] @@ -998,7 +995,7 @@ def test_orderCC2F(self): class TestAveraging3D(discretize.tests.OrderTest): name = "Averaging 3D" - meshTypes = ["notatreeTree", "uniformTree"] # , 'randomTree'] + meshTypes = ["notatreeTree", "uniformTree"] meshDimension = 3 meshSizes = [8, 16] expectedOrders = [2, 1] diff --git a/tests/tree/test_tree_plotting.py b/tests/tree/test_tree_plotting.py index a932ed037..8ed4623be 100644 --- a/tests/tree/test_tree_plotting.py +++ b/tests/tree/test_tree_plotting.py @@ -6,6 +6,8 @@ matplotlib.use("Agg") +rng = np.random.default_rng(4213678) + class TestOcTreePlotting(unittest.TestCase): def setUp(self): @@ -22,8 +24,8 @@ def test_plot_slice(self): mesh.plot_grid(faces=True, edges=True, nodes=True) # CC plot - mod_cc = np.random.rand(len(mesh)) + 1j * np.random.rand(len(mesh)) - mod_cc[np.random.rand(len(mesh)) < 0.2] = np.nan + mod_cc = rng.random(len(mesh)) + 1j * rng.random(len(mesh)) + mod_cc[rng.random(len(mesh)) < 0.2] = np.nan mesh.plot_slice(mod_cc, normal="X", grid=True) mesh.plot_slice(mod_cc, normal="Y", ax=ax) @@ -31,11 +33,11 @@ def test_plot_slice(self): mesh.plot_slice(mod_cc, view="imag", ax=ax) mesh.plot_slice(mod_cc, view="abs", ax=ax) - mod_ccv = np.random.rand(len(mesh), 3) + mod_ccv = rng.random((len(mesh), 3)) mesh.plot_slice(mod_ccv, v_type="CCv", view="vec", ax=ax) # F plot tests - mod_f = np.random.rand(mesh.n_faces) + mod_f = rng.random(mesh.n_faces) mesh.plot_slice(mod_f, v_type="Fx", ax=ax) mesh.plot_slice(mod_f, v_type="Fy", ax=ax) mesh.plot_slice(mod_f, v_type="Fz", ax=ax) @@ -43,7 +45,7 @@ def test_plot_slice(self): mesh.plot_slice(mod_f, v_type="F", view="vec", ax=ax) # E plot tests - mod_e = np.random.rand(mesh.n_edges) + mod_e = rng.random(mesh.n_edges) mesh.plot_slice(mod_e, v_type="Ex", ax=ax) mesh.plot_slice(mod_e, v_type="Ey", ax=ax) mesh.plot_slice(mod_e, v_type="Ez", ax=ax) @@ -51,7 +53,7 @@ def test_plot_slice(self): mesh.plot_slice(mod_e, v_type="E", view="vec", ax=ax) # Nodes - mod_n = np.random.rand(mesh.n_nodes) + mod_n = rng.random(mesh.n_nodes) mesh.plot_slice(mod_n, v_type="N") plt.close("all") @@ -70,32 +72,32 @@ def test_plot_slice(self): mesh.plot_grid(faces=True, edges=True, nodes=True) # CC plot - mod_cc = np.random.rand(len(mesh)) + 1j * np.random.rand(len(mesh)) - mod_cc[np.random.rand(len(mesh)) < 0.2] = np.nan + mod_cc = rng.random(len(mesh)) + 1j * rng.random(len(mesh)) + mod_cc[rng.random(len(mesh)) < 0.2] = np.nan mesh.plot_image(mod_cc) mesh.plot_image(mod_cc, ax=ax) mesh.plot_image(mod_cc, view="imag", ax=ax) mesh.plot_image(mod_cc, view="abs", ax=ax) - mod_ccv = np.random.rand(len(mesh), 2) + mod_ccv = rng.random((len(mesh), 2)) mesh.plot_image(mod_ccv, v_type="CCv", view="vec", ax=ax) # F plot tests - mod_f = np.random.rand(mesh.n_faces) + mod_f = rng.random(mesh.n_faces) mesh.plot_image(mod_f, v_type="Fx", ax=ax) mesh.plot_image(mod_f, v_type="Fy", ax=ax) mesh.plot_image(mod_f, v_type="F", ax=ax) mesh.plot_image(mod_f, v_type="F", view="vec", ax=ax) # E plot tests - mod_e = np.random.rand(mesh.n_edges) + mod_e = rng.random(mesh.n_edges) mesh.plot_image(mod_e, v_type="Ex", ax=ax) mesh.plot_image(mod_e, v_type="Ey", ax=ax) mesh.plot_image(mod_e, v_type="E", ax=ax) mesh.plot_image(mod_e, v_type="E", view="vec", ax=ax) # Nodes - mod_n = np.random.rand(mesh.n_nodes) + mod_n = rng.random(mesh.n_nodes) mesh.plot_image(mod_n, v_type="N", ax=ax) plt.close("all") diff --git a/tests/tree/test_tree_utils.py b/tests/tree/test_tree_utils.py index e868058a5..692d27f56 100644 --- a/tests/tree/test_tree_utils.py +++ b/tests/tree/test_tree_utils.py @@ -3,7 +3,6 @@ from discretize.utils import mesh_builder_xyz, refine_tree_xyz TOL = 1e-8 -np.random.seed(12) class TestRefineOcTree(unittest.TestCase): diff --git a/tutorials/inner_products/1_basic.py b/tutorials/inner_products/1_basic.py index bd86eec09..4d317251d 100644 --- a/tutorials/inner_products/1_basic.py +++ b/tutorials/inner_products/1_basic.py @@ -74,6 +74,8 @@ import matplotlib.pyplot as plt import numpy as np +rng = np.random.default_rng(8572) + # sphinx_gallery_thumbnail_number = 2 @@ -257,10 +259,10 @@ def fcn_y(xy, sig): Mf_inv = mesh.get_face_inner_product(invert_matrix=True) # Generate some random vectors -phi_c = np.random.rand(mesh.nC) -# phi_n = np.random.rand(mesh.nN) -vec_e = np.random.rand(mesh.nE) -vec_f = np.random.rand(mesh.nF) +phi_c = rng.random(mesh.nC) +# phi_n = rng.random(mesh.nN) +vec_e = rng.random(mesh.nE) +vec_f = rng.random(mesh.nF) # Generate some random vectors norm_c = np.linalg.norm(phi_c - Mc_inv.dot(Mc.dot(phi_c))) diff --git a/tutorials/inner_products/2_physical_properties.py b/tutorials/inner_products/2_physical_properties.py index 531e8b077..d5845da5d 100644 --- a/tutorials/inner_products/2_physical_properties.py +++ b/tutorials/inner_products/2_physical_properties.py @@ -75,6 +75,8 @@ import numpy as np import matplotlib.pyplot as plt +rng = np.random.default_rng(87236) + # sphinx_gallery_thumbnail_number = 1 ##################################################### @@ -181,17 +183,17 @@ mesh = TensorMesh([h, h, h]) # Isotropic case: (nC, ) numpy array -sig = np.random.rand(mesh.nC) # sig for each cell +sig = rng.random(mesh.nC) # sig for each cell Me1 = mesh.get_edge_inner_product(sig) # Edges inner product matrix Mf1 = mesh.get_face_inner_product(sig) # Faces inner product matrix # Linear case: (nC, dim) numpy array -sig = np.random.rand(mesh.nC, mesh.dim) +sig = rng.random((mesh.nC, mesh.dim)) Me2 = mesh.get_edge_inner_product(sig) Mf2 = mesh.get_face_inner_product(sig) # Anisotropic case: (nC, 3) for 2D and (nC, 6) for 3D -sig = np.random.rand(mesh.nC, 6) +sig = rng.random((mesh.nC, 6)) Me3 = mesh.get_edge_inner_product(sig) Mf3 = mesh.get_face_inner_product(sig) @@ -235,17 +237,17 @@ mesh = TensorMesh([h, h, h]) # Isotropic case: (nC, ) numpy array -sig = np.random.rand(mesh.nC) +sig = rng.random(mesh.nC) Me1_inv = mesh.get_edge_inner_product(sig, invert_matrix=True) Mf1_inv = mesh.get_face_inner_product(sig, invert_matrix=True) # Diagonal anisotropic: (nC, dim) numpy array -sig = np.random.rand(mesh.nC, mesh.dim) +sig = rng.random((mesh.nC, mesh.dim)) Me2_inv = mesh.get_edge_inner_product(sig, invert_matrix=True) Mf2_inv = mesh.get_face_inner_product(sig, invert_matrix=True) # Full anisotropic: (nC, 3) for 2D and (nC, 6) for 3D -sig = np.random.rand(mesh.nC, 6) +sig = rng.random((mesh.nC, 6)) Me3 = mesh.get_edge_inner_product(sig) Mf3 = mesh.get_face_inner_product(sig) diff --git a/tutorials/inner_products/4_advanced.py b/tutorials/inner_products/4_advanced.py index 8f97d1641..e384f5f24 100644 --- a/tutorials/inner_products/4_advanced.py +++ b/tutorials/inner_products/4_advanced.py @@ -21,6 +21,8 @@ import numpy as np import matplotlib.pyplot as plt +rng = np.random.default_rng(4321) + ##################################################### # Constitive Relations and Differential Operators @@ -74,8 +76,8 @@ # Make basic mesh h = np.ones(10) mesh = TensorMesh([h, h, h]) -sig = np.random.rand(mesh.nC) # isotropic -Sig = np.random.rand(mesh.nC, 6) # anisotropic +sig = rng.random(mesh.nC) # isotropic +Sig = rng.random((mesh.nC, 6)) # anisotropic # Inner product matricies Mc = sdiag(mesh.cell_volumes * sig) # Inner product matrix (centers) diff --git a/tutorials/pde/1_poisson.py b/tutorials/pde/1_poisson.py index 17dfc3095..d8665d1b7 100644 --- a/tutorials/pde/1_poisson.py +++ b/tutorials/pde/1_poisson.py @@ -92,7 +92,7 @@ from discretize import TensorMesh -from pymatsolver import SolverLU +from scipy.sparse.linalg import spsolve import matplotlib.pyplot as plt import numpy as np from discretize.utils import sdiag @@ -124,8 +124,7 @@ rho[kpos] = 1 # LU factorization and solve -AinvM = SolverLU(A) -phi = AinvM * rho +phi = spsolve(A, rho) # Compute electric fields E = Mf_inv * DIV.T * Mc * phi diff --git a/tutorials/pde/2_advection_diffusion.py b/tutorials/pde/2_advection_diffusion.py index 8783360db..7cc7bf326 100644 --- a/tutorials/pde/2_advection_diffusion.py +++ b/tutorials/pde/2_advection_diffusion.py @@ -164,7 +164,7 @@ # from discretize import TensorMesh -from pymatsolver import SolverLU +from scipy.sparse.linalg import splu import matplotlib.pyplot as plt import matplotlib as mpl import numpy as np @@ -224,7 +224,7 @@ B = I + dt * M s = Mc_inv * q -Binv = SolverLU(B) +Binv = splu(B) # Plot the vector field @@ -252,7 +252,7 @@ n = 3 for ii in range(300): - p = Binv * (p + s) + p = Binv.solve(p + s) if ii + 1 in (1, 25, 50, 100, 200, 300): ax[n] = fig.add_subplot(3, 3, n + 1) From 7a997c19a943e2cd2086e51c5f241363bcf2132f Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Sun, 1 Sep 2024 00:25:26 -0600 Subject: [PATCH 05/16] use seed in check derivative test calls --- tests/base/test_tensor_innerproduct_derivs.py | 12 +-- tests/base/test_tests.py | 6 +- tests/cyl/test_cyl_innerproducts.py | 88 ++++++++++++++----- tests/simplex/test_inner_products.py | 6 +- tests/tree/test_tree_innerproduct_derivs.py | 12 +-- 5 files changed, 87 insertions(+), 37 deletions(-) diff --git a/tests/base/test_tensor_innerproduct_derivs.py b/tests/base/test_tensor_innerproduct_derivs.py index 68d1cf2b8..e6e1ff21a 100644 --- a/tests/base/test_tensor_innerproduct_derivs.py +++ b/tests/base/test_tensor_innerproduct_derivs.py @@ -42,7 +42,7 @@ def fun(sig): fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=452) def doTestEdge( self, h, rep, fast, meshType, invert_model=False, invert_matrix=False @@ -79,7 +79,9 @@ def fun(sig): fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, rng=4567 + ) def test_FaceIP_1D_float(self): self.assertTrue(self.doTestFace([10], 0, False, "Tensor")) @@ -360,7 +362,7 @@ def fun(sig): rep, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=421) def doTestEdge(self, h, rep, meshType, invert_model=False, invert_matrix=False): if meshType == "Curv": @@ -391,7 +393,7 @@ def fun(sig): rep, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=31) def test_FaceIP_2D_float(self): self.assertTrue(self.doTestFace([10, 4], 0, "Tensor")) @@ -500,7 +502,7 @@ def fun(sig): rep, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=64) def test_EdgeIP_2D_float(self): self.assertTrue(self.doTestEdge([10, 4], 0, "Tensor")) diff --git a/tests/base/test_tests.py b/tests/base/test_tests.py index e8fc560c3..f14b0130b 100644 --- a/tests/base/test_tests.py +++ b/tests/base/test_tests.py @@ -74,14 +74,14 @@ def simplePass(x): return np.sin(x), sp.diags(np.cos(x)) rng = np.random.default_rng(5322) - check_derivative(simplePass, rng.standard_normal(5), plotIt=False, rng=rng) + check_derivative(simplePass, rng.standard_normal(5), plotIt=False, rng=42) def test_simpleFunction(self): def simpleFunction(x): return np.sin(x), lambda xi: np.cos(x) * xi rng = np.random.default_rng(5322) - check_derivative(simpleFunction, rng.standard_normal(5), plotIt=False, rng=rng) + check_derivative(simpleFunction, rng.standard_normal(5), plotIt=False, rng=23) def test_simpleFail(self): def simpleFail(x): @@ -89,7 +89,7 @@ def simpleFail(x): rng = np.random.default_rng(5322) with pytest.raises(AssertionError): - check_derivative(simpleFail, rng.standard_normal(5), plotIt=False, rng=rng) + check_derivative(simpleFail, rng.standard_normal(5), plotIt=False, rng=64) @pytest.mark.parametrize("test_type", ["mean", "min", "last", "all", "mean_at_least"]) diff --git a/tests/cyl/test_cyl_innerproducts.py b/tests/cyl/test_cyl_innerproducts.py index d74859485..b3e96f18e 100644 --- a/tests/cyl/test_cyl_innerproducts.py +++ b/tests/cyl/test_cyl_innerproducts.py @@ -478,7 +478,9 @@ def fun(x): print("Testing FaceInnerProduct Isotropic") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=532 + ) ) def test_FaceInnerProductIsotropicDerivInvProp(self): @@ -491,7 +493,9 @@ def fun(x): print("Testing FaceInnerProduct Isotropic InvProp") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=75 + ) ) def test_FaceInnerProductIsotropicDerivInvMat(self): @@ -504,7 +508,9 @@ def fun(x): print("Testing FaceInnerProduct Isotropic InvMat") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=1 + ) ) def test_FaceInnerProductIsotropicDerivInvPropInvMat(self): @@ -519,7 +525,9 @@ def fun(x): print("Testing FaceInnerProduct Isotropic InvProp InvMat") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=74 + ) ) def test_EdgeInnerProductIsotropicDeriv(self): @@ -530,7 +538,9 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=345 + ) ) def test_EdgeInnerProductIsotropicDerivInvProp(self): @@ -543,7 +553,9 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic InvProp") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=643 + ) ) def test_EdgeInnerProductIsotropicDerivInvMat(self): @@ -556,7 +568,9 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic InvMat") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=363 + ) ) def test_EdgeInnerProductIsotropicDerivInvPropInvMat(self): @@ -571,7 +585,9 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic InvProp InvMat") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=773 + ) ) @@ -603,7 +619,9 @@ def fun(x): print("Testing FaceInnerProduct Anisotropic") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=2436 + ) ) def test_FaceInnerProductAnisotropicDerivInvProp(self): @@ -621,7 +639,9 @@ def fun(x): print("Testing FaceInnerProduct Anisotropic InvProp") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=634 + ) ) def test_FaceInnerProductAnisotropicDerivInvMat(self): @@ -639,7 +659,9 @@ def fun(x): print("Testing FaceInnerProduct Anisotropic InvMat") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=222 + ) ) def test_FaceInnerProductAnisotropicDerivInvPropInvMat(self): @@ -661,7 +683,9 @@ def fun(x): print("Testing FaceInnerProduct Anisotropic InvProp InvMat") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=654 + ) ) def test_EdgeInnerProductAnisotropicDeriv(self): @@ -679,7 +703,9 @@ def fun(x): print("Testing EdgeInnerProduct Anisotropic") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=7754 + ) ) def test_EdgeInnerProductAnisotropicDerivInvProp(self): @@ -697,7 +723,9 @@ def fun(x): print("Testing EdgeInnerProduct Anisotropic InvProp") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=1164 + ) ) def test_EdgeInnerProductAnisotropicDerivInvMat(self): @@ -715,7 +743,9 @@ def fun(x): print("Testing EdgeInnerProduct Anisotropic InvMat") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=643 + ) ) def test_EdgeInnerProductAnisotropicDerivInvPropInvMat(self): @@ -737,7 +767,9 @@ def fun(x): print("Testing EdgeInnerProduct Anisotropic InvProp InvMat") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=8654 + ) ) @@ -761,7 +793,9 @@ def fun(x): print("Testing FaceInnerProduct Isotropic (Face Properties)") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=234 + ) ) def test_FaceInnerProductIsotropicDerivInvProp(self): @@ -774,7 +808,9 @@ def fun(x): print("Testing FaceInnerProduct Isotropic InvProp (Face Properties)") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=7543 + ) ) def test_FaceInnerProductIsotropicDerivInvMat(self): @@ -787,7 +823,9 @@ def fun(x): print("Testing FaceInnerProduct Isotropic InvMat (Face Properties)") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=2745725 + ) ) def test_EdgeInnerProductIsotropicDeriv(self): @@ -798,7 +836,9 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic (Face Properties)") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=6654 + ) ) def test_EdgeInnerProductIsotropicDerivInvProp(self): @@ -811,7 +851,9 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic InvProp (Face Properties)") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=4564 + ) ) def test_EdgeInnerProductIsotropicDerivInvMat(self): @@ -824,5 +866,7 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic InvMat (Face Properties)") return self.assertTrue( - tests.check_derivative(fun, self.x0, num=7, tolerance=TOLD, plotIt=False) + tests.check_derivative( + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=2355 + ) ) diff --git a/tests/simplex/test_inner_products.py b/tests/simplex/test_inner_products.py index 83f6875a4..676faa027 100644 --- a/tests/simplex/test_inner_products.py +++ b/tests/simplex/test_inner_products.py @@ -343,7 +343,9 @@ def fun(sig): return M * v, Md(v) print("Face", rep) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, rng=5352 + ) def doTestEdge(self, h, rep): nodes, simplices = example_simplex_mesh(h) @@ -357,7 +359,7 @@ def fun(sig): return M * v, Md(v) print("Edge", rep) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=532) def test_FaceIP_2D_float(self): self.assertTrue(self.doTestFace([10, 4], 0)) diff --git a/tests/tree/test_tree_innerproduct_derivs.py b/tests/tree/test_tree_innerproduct_derivs.py index 4bd060cea..708680c7f 100644 --- a/tests/tree/test_tree_innerproduct_derivs.py +++ b/tests/tree/test_tree_innerproduct_derivs.py @@ -40,7 +40,9 @@ def fun(sig): fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=rng) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, rng=4421 + ) def doTestEdge( self, h, rep, fast, meshType, invert_model=False, invert_matrix=False @@ -76,7 +78,7 @@ def fun(sig): fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=rng) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=643) def test_FaceIP_2D_float_Tree(self): self.assertTrue(self.doTestFace([8, 8], 0, False, "Tree")) @@ -193,7 +195,7 @@ def fun(sig): # fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=rng) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=677) def doTestEdge(self, h, rep, meshType, invert_model=False, invert_matrix=False): if meshType == "Curv": @@ -226,7 +228,7 @@ def fun(sig): # fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=rng) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=543) def test_FaceIP_2D_float_fast_Tree(self): self.assertTrue(self.doTestFace([8, 8], 0, "Tree")) @@ -285,7 +287,7 @@ def fun(sig): # fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=rng) + return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=23) def test_EdgeIP_2D_float_fast_Tree(self): self.assertTrue(self.doTestEdge([8, 8], 0, "Tree")) From 40df4b73f5a922703029d167809d5f36d5ec35a6 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Sun, 1 Sep 2024 00:30:44 -0600 Subject: [PATCH 06/16] use seed in adjoint tests --- tests/base/test_tests.py | 11 ++++++++++- 1 file changed, 10 insertions(+), 1 deletion(-) diff --git a/tests/base/test_tests.py b/tests/base/test_tests.py index f14b0130b..9a4475ca2 100644 --- a/tests/base/test_tests.py +++ b/tests/base/test_tests.py @@ -26,6 +26,7 @@ def test_defaults(self, capsys): mesh1.n_cells, mesh2.n_cells, assert_error=False, + rng=41, ) out2, _ = capsys.readouterr() assert out1 @@ -38,6 +39,7 @@ def test_defaults(self, capsys): lambda v: P.T * v, mesh1.n_cells, mesh2.n_cells, + rng=42, ) def test_different_shape(self): @@ -53,7 +55,13 @@ def adj(inp): out = np.expand_dims(inp, 1) return np.tile(out, nt) - assert_isadjoint(fwd, adj, shape_u=(4, nt), shape_v=(4,)) + assert_isadjoint( + fwd, + adj, + shape_u=(4, nt), + shape_v=(4,), + rng=42, + ) def test_complex_clinear(self): # The complex conjugate is self-adjoint, real-linear. @@ -65,6 +73,7 @@ def test_complex_clinear(self): complex_u=True, complex_v=True, clinear=False, + rng=112, ) From 2702e2d55526f0472bb08f150a80c5567bc73065 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Sun, 1 Sep 2024 00:35:50 -0600 Subject: [PATCH 07/16] incorrect expansion --- discretize/utils/mesh_utils.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/discretize/utils/mesh_utils.py b/discretize/utils/mesh_utils.py index 22ea4fac7..7aeeaf0d8 100644 --- a/discretize/utils/mesh_utils.py +++ b/discretize/utils/mesh_utils.py @@ -75,7 +75,7 @@ def random_model(shape, seed=None, anisotropy=None, its=100, bounds=None): if type(shape) in num_types: shape = (shape,) # make it a tuple for consistency - mr = rng.random(*shape) + mr = rng.random(shape) if anisotropy is None: if len(shape) == 1: smth = np.array([1, 10.0, 1], dtype=float) From 6c862848c2c889388d9cfa2bf1ccf1075b54ca1f Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Sun, 1 Sep 2024 00:35:57 -0600 Subject: [PATCH 08/16] correct import --- examples/plot_dc_resistivity.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/plot_dc_resistivity.py b/examples/plot_dc_resistivity.py index 4f1cd8eac..8bfdf76f0 100644 --- a/examples/plot_dc_resistivity.py +++ b/examples/plot_dc_resistivity.py @@ -8,7 +8,7 @@ import discretize import numpy as np import matplotlib.pyplot as plt -from scipy.solver.linalg import spsolve +from scipy.sparse.linalg import spsolve def run(plotIt=True): From 567bce0247b45c6b65980bbf615b21a29c900094 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Sun, 1 Sep 2024 00:39:24 -0600 Subject: [PATCH 09/16] remove pymatsolver from testing and doc requirements --- .ci/environment_test.yml | 1 - pyproject.toml | 1 - 2 files changed, 2 deletions(-) diff --git a/.ci/environment_test.yml b/.ci/environment_test.yml index 91748c856..bdff1b8ce 100644 --- a/.ci/environment_test.yml +++ b/.ci/environment_test.yml @@ -16,7 +16,6 @@ dependencies: - numpydoc>=1.5 - jupyter - graphviz - - pymatsolver>=0.1.2 - pillow - pooch # testing diff --git a/pyproject.toml b/pyproject.toml index ca54b91cd..3d777ef86 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -70,7 +70,6 @@ doc = [ "numpydoc>=1.5", "jupyter", "graphviz", - "pymatsolver>=0.1.2", "pillow", "pooch", "discretize[all]", From 7ffd4e4a5050ce7cc8c136d8868aaaaf5df20706 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Sun, 1 Sep 2024 12:30:28 -0600 Subject: [PATCH 10/16] add some spaces to environment_test --- .ci/environment_test.yml | 3 +++ 1 file changed, 3 insertions(+) diff --git a/.ci/environment_test.yml b/.ci/environment_test.yml index bdff1b8ce..90568a70b 100644 --- a/.ci/environment_test.yml +++ b/.ci/environment_test.yml @@ -4,11 +4,13 @@ channels: dependencies: - numpy>=1.22.4 - scipy>=1.8 + # optionals - vtk>=6 - pyvista - omf - matplotlib + # documentation - sphinx - pydata-sphinx-theme==0.13.3 @@ -18,6 +20,7 @@ dependencies: - graphviz - pillow - pooch + # testing - sympy - pytest From 1fe51228ae4799972839154018be477fcb790686 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Sun, 1 Sep 2024 12:38:45 -0600 Subject: [PATCH 11/16] minor syntax errors in old docs --- docs/release/0.7.1-notes.rst | 2 +- docs/release/0.8.1-notes.rst | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/docs/release/0.7.1-notes.rst b/docs/release/0.7.1-notes.rst index 1eeab2052..292e52d77 100644 --- a/docs/release/0.7.1-notes.rst +++ b/docs/release/0.7.1-notes.rst @@ -42,4 +42,4 @@ Pull requests * `#256 `__: Update mpl_mod.py * `#258 `__: Numpy docstrings api review * `#262 `__: Fix wrong colour for fullspaces - again / -`#264 `__: patch for fullspace slicer colorscales +* `#264 `__: patch for fullspace slicer colorscales diff --git a/docs/release/0.8.1-notes.rst b/docs/release/0.8.1-notes.rst index e3536e2f3..a2ac19517 100644 --- a/docs/release/0.8.1-notes.rst +++ b/docs/release/0.8.1-notes.rst @@ -54,5 +54,5 @@ Pull requests * `#284 `__: Improve load time * `#285 `__: zeros_outside * `#286 `__: Cell node tree -* `#288 `__: Allow np.int_ +* `#288 `__: Allow ``np.int_`` * `#289 `__: 0.8.1 Release From 88b5d0a12c3cdb74895bd5307a91d50356007069 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Sun, 1 Sep 2024 12:39:43 -0600 Subject: [PATCH 12/16] fix invalid inline literal --- docs/release/0.10.0-notes.rst | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/docs/release/0.10.0-notes.rst b/docs/release/0.10.0-notes.rst index 7cf5c73f7..05e0e71ab 100644 --- a/docs/release/0.10.0-notes.rst +++ b/docs/release/0.10.0-notes.rst @@ -22,7 +22,7 @@ Build system ``discretize`` now uses a ``pyproject.toml`` file with a ``meson-python`` backend to build the compiled external modules (used for the ``TreeMesh``, ``SimplexMesh``, and interpolation functions. Moving away from a ``setup.py`` file allows us to reliably control the build environment -seperate from the install environment, as the build requirements are not the same as the runtime +separate from the install environment, as the build requirements are not the same as the runtime requirements. We will also begin distributing many more pre-compiled wheels on pypi for Windows, MacOS, and Linux @@ -30,11 +30,11 @@ systems from python 3.8 to 3.12. Our goal is to provide a pre-compiled wheel for scipy provides wheels for. Tensor Mesh ------------- -You can now directly index a ``TensorMesh``, and it will then return a ``TensorCell`` object. -This functionality mimics what is currently available in ``TreeMesh``s. +----------- +You can now directly index a ``TensorMesh`` and return a ``TensorCell`` object. +This functionality mimics what is currently available in ``TreeMesh``. -``Tensor Mesh`` also has a ``cell_nodes`` property that list the indices of each node of every +``TensorMesh`` also has a ``cell_nodes`` property that list the indices of each node of every cell (again similar to the ``TreeMesh``). Face Properties From 9ffbf260f76bb2a27357b00310f294a5c71cc3e6 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Sun, 1 Sep 2024 13:01:51 -0600 Subject: [PATCH 13/16] update docstrings --- discretize/tests.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/discretize/tests.py b/discretize/tests.py index 4848371e4..ea4340a29 100644 --- a/discretize/tests.py +++ b/discretize/tests.py @@ -227,7 +227,7 @@ class for the given operator. Within the test class, the user sets the parameter Name the convergence test meshTypes : list of str List denoting the mesh types on which the convergence will be tested. - List entries must of the list {'uniformTensorMesh', 'randomTensorMesh', + List entries must be of the list {'uniformTensorMesh', 'randomTensorMesh', 'uniformCylindricalMesh', 'randomCylindricalMesh', 'uniformTree', 'randomTree', 'uniformCurv', 'rotateCurv', 'sphereCurv'} expectedOrders : float or list of float (default = 2.0) @@ -242,8 +242,9 @@ class for the given operator. Within the test class, the user sets the parameter meshDimension : int Mesh dimension. Must be 1, 2 or 3 rng : numpy.random.Generator, int, optional - The random number generator to use for the adjoint test, if an integer or None - it used to seed a new `numpy.random.default_rng`. Only used if mesh_type is 'random' + If ``random`` is in `mesh_type`, this is the random number generator + used generate the random meshes, if an ``int`` or ``None``, it used to seed + a new `numpy.random.default_rng`. Notes ----- From 7c84ba282051675801d86acf60d605bd2c711223 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Tue, 3 Sep 2024 21:56:42 -0600 Subject: [PATCH 14/16] change rng argument to random_seed --- discretize/tests.py | 43 +++++++++--------- tests/base/test_tensor_innerproduct_derivs.py | 18 +++++--- tests/base/test_tests.py | 20 ++++++--- tests/cyl/test_cyl_innerproducts.py | 44 +++++++++---------- tests/simplex/test_inner_products.py | 6 ++- tests/tree/test_tree_innerproduct_derivs.py | 18 +++++--- 6 files changed, 86 insertions(+), 63 deletions(-) diff --git a/discretize/tests.py b/discretize/tests.py index ea4340a29..eb122aa3a 100644 --- a/discretize/tests.py +++ b/discretize/tests.py @@ -81,7 +81,7 @@ _happiness_rng = np.random.default_rng() -def setup_mesh(mesh_type, nC, nDim, rng=None): +def setup_mesh(mesh_type, nC, nDim, random_seed=None): """Generate arbitrary mesh for testing. For the mesh type, number of cells along each axis and dimension specified, @@ -99,9 +99,9 @@ def setup_mesh(mesh_type, nC, nDim, rng=None): number of base mesh cells and must be a power of 2. nDim : int The dimension of the mesh. Must be 1, 2 or 3. - rng : numpy.random.Generator, int, optional + random_seed : numpy.random.Generator, int, optional If ``random`` is in `mesh_type`, this is the random number generator to use for - creating a random mesh. If an integer or None it is used to seed a new + creating that random mesh. If an integer or None it is used to seed a new `numpy.random.default_rng`. Returns @@ -110,7 +110,7 @@ def setup_mesh(mesh_type, nC, nDim, rng=None): A discretize mesh of class specified by the input argument *mesh_type* """ if "random" in mesh_type: - rng = np.random.default_rng(rng) + rng = np.random.default_rng(random_seed) if "TensorMesh" in mesh_type: if "uniform" in mesh_type: h = [nC, nC, nC] @@ -221,7 +221,7 @@ class for the given operator. Within the test class, the user sets the parameter OrderTest inherits from :class:`unittest.TestCase`. - Parameters + Attributes ---------- name : str Name the convergence test @@ -241,7 +241,7 @@ class for the given operator. Within the test class, the user sets the parameter for the meshes used in the convergence test; e.g. [4, 8, 16, 32] meshDimension : int Mesh dimension. Must be 1, 2 or 3 - rng : numpy.random.Generator, int, optional + random_seed : numpy.random.Generator, int, optional If ``random`` is in `mesh_type`, this is the random number generator used generate the random meshes, if an ``int`` or ``None``, it used to seed a new `numpy.random.default_rng`. @@ -311,7 +311,7 @@ class for the given operator. Within the test class, the user sets the parameter meshTypes = ["uniformTensorMesh"] _meshType = meshTypes[0] meshDimension = 3 - rng = None + random_seed = None def setupMesh(self, nC): """Generate mesh and set as current mesh for testing. @@ -320,16 +320,15 @@ def setupMesh(self, nC): ---------- nC : int Number of cells along each axis. - rng : numpy.random.Generator, int, optional - The random number generator to use for the adjoint test, if an integer or None - it used to seed a new `numpy.random.default_rng`. Only used if mesh_type is 'random' Returns ------- Float Maximum cell width for the mesh """ - mesh, max_h = setup_mesh(self._meshType, nC, self.meshDimension, rng=self.rng) + mesh, max_h = setup_mesh( + self._meshType, nC, self.meshDimension, random_seed=self.random_seed + ) self.M = mesh return max_h @@ -348,7 +347,7 @@ def getError(self): """ return 1.0 - def orderTest(self, rng=None): + def orderTest(self, random_seed=None): """Perform an order test. For number of cells specified in meshSizes setup mesh, call getError @@ -373,8 +372,8 @@ def orderTest(self, rng=None): "expectedOrders must have the same length as the meshTypes" ) - if rng is not None: - self.rng = rng + if random_seed is not None: + self.random_seed = random_seed def test_func(n_cells): max_h = self.setupMesh(n_cells) @@ -568,7 +567,7 @@ def check_derivative( tolerance=0.85, eps=1e-10, ax=None, - rng=None, + random_seed=None, ): """Perform a basic derivative check. @@ -597,8 +596,8 @@ def check_derivative( ax : matplotlib.pyplot.Axes, optional An axis object for the convergence plot if *plotIt = True*. Otherwise, the function will create a new axis. - rng : numpy.random.Generator, int, optional - If `dx` is ``None``, This is the random number generator to use for + random_seed : numpy.random.Generator, int, optional + If `dx` is ``None``, this is the random number generator to use for generating a step direction. If an integer or None, it is used to seed a new `numpy.random.default_rng`. @@ -616,7 +615,7 @@ def check_derivative( >>> def simplePass(x): ... return np.sin(x), utils.sdiag(np.cos(x)) - >>> passed = tests.check_derivative(simplePass, rng.standard_normal(5), rng=rng) + >>> passed = tests.check_derivative(simplePass, rng.standard_normal(5), random_seed=rng) ==================== check_derivative ==================== iter h |ft-f0| |ft-f0-h*J0*dx| Order --------------------------------------------------------- @@ -650,7 +649,7 @@ def check_derivative( x0 = mkvc(x0) if dx is None: - rng = np.random.default_rng(rng) + rng = np.random.default_rng(random_seed) dx = rng.standard_normal(len(x0)) h = np.logspace(-1, -num, num) @@ -800,7 +799,7 @@ def assert_isadjoint( rtol=1e-6, atol=0.0, assert_error=True, - rng=None, + random_seed=None, ): r"""Do a dot product test for the forward operator and its adjoint operator. @@ -851,7 +850,7 @@ def assert_isadjoint( assertion error if failed). If set to False, the result of the test is returned as boolean and a message is printed. - rng : numpy.random.Generator, int, optional + random_seed : numpy.random.Generator, int, optional The random number generator to use for the adjoint test. If an integer or None it is used to seed a new `numpy.random.default_rng`. @@ -868,7 +867,7 @@ def assert_isadjoint( """ __tracebackhide__ = True - rng = np.random.default_rng(rng) + rng = np.random.default_rng(random_seed) def random(size, iscomplex): """Create random data of size and dtype of .""" diff --git a/tests/base/test_tensor_innerproduct_derivs.py b/tests/base/test_tensor_innerproduct_derivs.py index e6e1ff21a..af8a94f70 100644 --- a/tests/base/test_tensor_innerproduct_derivs.py +++ b/tests/base/test_tensor_innerproduct_derivs.py @@ -42,7 +42,9 @@ def fun(sig): fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=452) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, random_seed=452 + ) def doTestEdge( self, h, rep, fast, meshType, invert_model=False, invert_matrix=False @@ -80,7 +82,7 @@ def fun(sig): ("harmonic" if invert_model and invert_matrix else "standard"), ) return discretize.tests.check_derivative( - fun, sig, num=5, plotIt=False, rng=4567 + fun, sig, num=5, plotIt=False, random_seed=4567 ) def test_FaceIP_1D_float(self): @@ -362,7 +364,9 @@ def fun(sig): rep, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=421) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, random_seed=421 + ) def doTestEdge(self, h, rep, meshType, invert_model=False, invert_matrix=False): if meshType == "Curv": @@ -393,7 +397,9 @@ def fun(sig): rep, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=31) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, random_seed=31 + ) def test_FaceIP_2D_float(self): self.assertTrue(self.doTestFace([10, 4], 0, "Tensor")) @@ -502,7 +508,9 @@ def fun(sig): rep, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=64) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, random_seed=64 + ) def test_EdgeIP_2D_float(self): self.assertTrue(self.doTestEdge([10, 4], 0, "Tensor")) diff --git a/tests/base/test_tests.py b/tests/base/test_tests.py index 9a4475ca2..fbf789ca2 100644 --- a/tests/base/test_tests.py +++ b/tests/base/test_tests.py @@ -26,7 +26,7 @@ def test_defaults(self, capsys): mesh1.n_cells, mesh2.n_cells, assert_error=False, - rng=41, + random_seed=41, ) out2, _ = capsys.readouterr() assert out1 @@ -39,7 +39,7 @@ def test_defaults(self, capsys): lambda v: P.T * v, mesh1.n_cells, mesh2.n_cells, - rng=42, + random_seed=42, ) def test_different_shape(self): @@ -60,7 +60,7 @@ def adj(inp): adj, shape_u=(4, nt), shape_v=(4,), - rng=42, + random_seed=42, ) def test_complex_clinear(self): @@ -73,7 +73,7 @@ def test_complex_clinear(self): complex_u=True, complex_v=True, clinear=False, - rng=112, + random_seed=112, ) @@ -83,14 +83,18 @@ def simplePass(x): return np.sin(x), sp.diags(np.cos(x)) rng = np.random.default_rng(5322) - check_derivative(simplePass, rng.standard_normal(5), plotIt=False, rng=42) + check_derivative( + simplePass, rng.standard_normal(5), plotIt=False, random_seed=42 + ) def test_simpleFunction(self): def simpleFunction(x): return np.sin(x), lambda xi: np.cos(x) * xi rng = np.random.default_rng(5322) - check_derivative(simpleFunction, rng.standard_normal(5), plotIt=False, rng=23) + check_derivative( + simpleFunction, rng.standard_normal(5), plotIt=False, random_seed=23 + ) def test_simpleFail(self): def simpleFail(x): @@ -98,7 +102,9 @@ def simpleFail(x): rng = np.random.default_rng(5322) with pytest.raises(AssertionError): - check_derivative(simpleFail, rng.standard_normal(5), plotIt=False, rng=64) + check_derivative( + simpleFail, rng.standard_normal(5), plotIt=False, random_seed=64 + ) @pytest.mark.parametrize("test_type", ["mean", "min", "last", "all", "mean_at_least"]) diff --git a/tests/cyl/test_cyl_innerproducts.py b/tests/cyl/test_cyl_innerproducts.py index b3e96f18e..875b54140 100644 --- a/tests/cyl/test_cyl_innerproducts.py +++ b/tests/cyl/test_cyl_innerproducts.py @@ -479,7 +479,7 @@ def fun(x): print("Testing FaceInnerProduct Isotropic") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=532 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=532 ) ) @@ -494,7 +494,7 @@ def fun(x): print("Testing FaceInnerProduct Isotropic InvProp") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=75 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=75 ) ) @@ -509,7 +509,7 @@ def fun(x): print("Testing FaceInnerProduct Isotropic InvMat") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=1 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=1 ) ) @@ -526,7 +526,7 @@ def fun(x): print("Testing FaceInnerProduct Isotropic InvProp InvMat") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=74 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=74 ) ) @@ -539,7 +539,7 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=345 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=345 ) ) @@ -554,7 +554,7 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic InvProp") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=643 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=643 ) ) @@ -569,7 +569,7 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic InvMat") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=363 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=363 ) ) @@ -586,7 +586,7 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic InvProp InvMat") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=773 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=773 ) ) @@ -620,7 +620,7 @@ def fun(x): print("Testing FaceInnerProduct Anisotropic") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=2436 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=2436 ) ) @@ -640,7 +640,7 @@ def fun(x): print("Testing FaceInnerProduct Anisotropic InvProp") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=634 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=634 ) ) @@ -660,7 +660,7 @@ def fun(x): print("Testing FaceInnerProduct Anisotropic InvMat") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=222 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=222 ) ) @@ -684,7 +684,7 @@ def fun(x): print("Testing FaceInnerProduct Anisotropic InvProp InvMat") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=654 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=654 ) ) @@ -704,7 +704,7 @@ def fun(x): print("Testing EdgeInnerProduct Anisotropic") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=7754 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=7754 ) ) @@ -724,7 +724,7 @@ def fun(x): print("Testing EdgeInnerProduct Anisotropic InvProp") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=1164 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=1164 ) ) @@ -744,7 +744,7 @@ def fun(x): print("Testing EdgeInnerProduct Anisotropic InvMat") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=643 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=643 ) ) @@ -768,7 +768,7 @@ def fun(x): print("Testing EdgeInnerProduct Anisotropic InvProp InvMat") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=8654 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=8654 ) ) @@ -794,7 +794,7 @@ def fun(x): print("Testing FaceInnerProduct Isotropic (Face Properties)") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=234 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=234 ) ) @@ -809,7 +809,7 @@ def fun(x): print("Testing FaceInnerProduct Isotropic InvProp (Face Properties)") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=7543 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=7543 ) ) @@ -824,7 +824,7 @@ def fun(x): print("Testing FaceInnerProduct Isotropic InvMat (Face Properties)") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=2745725 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=2745725 ) ) @@ -837,7 +837,7 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic (Face Properties)") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=6654 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=6654 ) ) @@ -852,7 +852,7 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic InvProp (Face Properties)") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=4564 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=4564 ) ) @@ -867,6 +867,6 @@ def fun(x): print("Testing EdgeInnerProduct Isotropic InvMat (Face Properties)") return self.assertTrue( tests.check_derivative( - fun, self.x0, num=7, tolerance=TOLD, plotIt=False, rng=2355 + fun, self.x0, num=7, tolerance=TOLD, plotIt=False, random_seed=2355 ) ) diff --git a/tests/simplex/test_inner_products.py b/tests/simplex/test_inner_products.py index 676faa027..a5cae4ce6 100644 --- a/tests/simplex/test_inner_products.py +++ b/tests/simplex/test_inner_products.py @@ -344,7 +344,7 @@ def fun(sig): print("Face", rep) return discretize.tests.check_derivative( - fun, sig, num=5, plotIt=False, rng=5352 + fun, sig, num=5, plotIt=False, random_seed=5352 ) def doTestEdge(self, h, rep): @@ -359,7 +359,9 @@ def fun(sig): return M * v, Md(v) print("Edge", rep) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=532) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, random_seed=532 + ) def test_FaceIP_2D_float(self): self.assertTrue(self.doTestFace([10, 4], 0)) diff --git a/tests/tree/test_tree_innerproduct_derivs.py b/tests/tree/test_tree_innerproduct_derivs.py index 708680c7f..da4ff3026 100644 --- a/tests/tree/test_tree_innerproduct_derivs.py +++ b/tests/tree/test_tree_innerproduct_derivs.py @@ -41,7 +41,7 @@ def fun(sig): ("harmonic" if invert_model and invert_matrix else "standard"), ) return discretize.tests.check_derivative( - fun, sig, num=5, plotIt=False, rng=4421 + fun, sig, num=5, plotIt=False, random_seed=4421 ) def doTestEdge( @@ -78,7 +78,9 @@ def fun(sig): fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=643) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, random_seed=643 + ) def test_FaceIP_2D_float_Tree(self): self.assertTrue(self.doTestFace([8, 8], 0, False, "Tree")) @@ -195,7 +197,9 @@ def fun(sig): # fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=677) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, random_seed=677 + ) def doTestEdge(self, h, rep, meshType, invert_model=False, invert_matrix=False): if meshType == "Curv": @@ -228,7 +232,9 @@ def fun(sig): # fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=543) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, random_seed=543 + ) def test_FaceIP_2D_float_fast_Tree(self): self.assertTrue(self.doTestFace([8, 8], 0, "Tree")) @@ -287,7 +293,9 @@ def fun(sig): # fast, ("harmonic" if invert_model and invert_matrix else "standard"), ) - return discretize.tests.check_derivative(fun, sig, num=5, plotIt=False, rng=23) + return discretize.tests.check_derivative( + fun, sig, num=5, plotIt=False, random_seed=23 + ) def test_EdgeIP_2D_float_fast_Tree(self): self.assertTrue(self.doTestEdge([8, 8], 0, "Tree")) From 7d3617ddb06f4e50040e66e373f2b77ab2b9cbfb Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Wed, 4 Sep 2024 15:14:42 -0600 Subject: [PATCH 15/16] choose working seeds --- tests/tree/test_tree_operators.py | 15 +++++---------- 1 file changed, 5 insertions(+), 10 deletions(-) diff --git a/tests/tree/test_tree_operators.py b/tests/tree/test_tree_operators.py index 483385d70..12a9b1ec2 100644 --- a/tests/tree/test_tree_operators.py +++ b/tests/tree/test_tree_operators.py @@ -38,7 +38,6 @@ class TestCellGrad2D(discretize.tests.OrderTest): meshSizes = [8, 16] # because of the averaging involved in the ghost point. u_b = (u_n + u_g)/2 expectedOrders = 1 - rng = np.random.default_rng(87964213) def getError(self): # Test function @@ -56,7 +55,7 @@ def getError(self): return err def test_order(self): - self.orderTest() + self.orderTest(random_seed=421) class TestCellGrad3D(discretize.tests.OrderTest): @@ -66,7 +65,6 @@ class TestCellGrad3D(discretize.tests.OrderTest): meshSizes = [8, 16] # because of the averaging involved in the ghost point. u_b = (u_n + u_g)/2 expectedOrders = 1 - rng = np.random.default_rng(6957832) def getError(self): # Test function @@ -106,7 +104,7 @@ def getError(self): return err def test_order(self): - self.orderTest() + self.orderTest(5532) class TestFaceDivxy2D(discretize.tests.OrderTest): @@ -114,7 +112,6 @@ class TestFaceDivxy2D(discretize.tests.OrderTest): meshTypes = MESHTYPES meshDimension = 2 meshSizes = [16, 32] - rng = np.random.default_rng(19647823) def getError(self): # Test function @@ -137,14 +134,13 @@ def getError(self): return err def test_order(self): - self.orderTest() + self.orderTest(random_seed=19647823) class TestFaceDiv3D(discretize.tests.OrderTest): name = "Face Divergence 3D" meshTypes = MESHTYPES meshSizes = [8, 16, 32] - rng = np.random.default_rng(81725364) def getError(self): fx = lambda x, y, z: np.sin(2 * np.pi * x) @@ -165,7 +161,7 @@ def getError(self): return np.linalg.norm((divF - divF_ana), np.inf) def test_order(self): - self.orderTest() + self.orderTest(random_seed=81725364) class TestFaceDivxyz3D(discretize.tests.OrderTest): @@ -173,7 +169,6 @@ class TestFaceDivxyz3D(discretize.tests.OrderTest): meshTypes = MESHTYPES meshDimension = 3 meshSizes = [8, 16, 32] - rng = np.random.default_rng(6172824) def getError(self): # Test function @@ -203,7 +198,7 @@ def getError(self): return err def test_order(self): - self.orderTest() + self.orderTest(random_seed=6172824) class TestCurl(discretize.tests.OrderTest): From 9e83796bead117d930c6705363171e48e56c5fc3 Mon Sep 17 00:00:00 2001 From: Joseph Capriotti Date: Wed, 4 Sep 2024 15:28:41 -0600 Subject: [PATCH 16/16] more switches from rng->random_seed --- tests/base/test_interpolation.py | 28 ++++++++++++++++------------ tests/base/test_tensor.py | 1 - 2 files changed, 16 insertions(+), 13 deletions(-) diff --git a/tests/base/test_interpolation.py b/tests/base/test_interpolation.py index bd46a2c93..1b433e9bd 100644 --- a/tests/base/test_interpolation.py +++ b/tests/base/test_interpolation.py @@ -46,8 +46,8 @@ class TestInterpolation1D(discretize.tests.OrderTest): tolerance = TOLERANCES meshDimension = 1 meshSizes = [8, 16, 32, 64, 128] - rng = np.random.default_rng(5136) - LOCS = rng.random(50) * 0.6 + 0.2 + random_seed = np.random.default_rng(55124) + LOCS = random_seed.random(50) * 0.6 + 0.2 def getError(self): funX = lambda x: np.cos(2 * np.pi * x) @@ -98,8 +98,8 @@ class TestInterpolation2d(discretize.tests.OrderTest): tolerance = TOLERANCES meshDimension = 2 meshSizes = [8, 16, 32, 64] - rng = np.random.default_rng(2457) - LOCS = rng.random((50, 2)) * 0.6 + 0.2 + random_seed = np.random.default_rng(2457) + LOCS = random_seed.random((50, 2)) * 0.6 + 0.2 def getError(self): funX = lambda x, y: np.cos(2 * np.pi * y) @@ -182,8 +182,12 @@ class TestInterpolationSymCyl(discretize.tests.OrderTest): tolerance = 0.6 meshDimension = 3 meshSizes = [32, 64, 128, 256] - rng = np.random.default_rng(81756234) - LOCS = np.c_[rng.random(4) * 0.6 + 0.2, np.zeros(4), rng.random(4) * 0.6 + 0.2] + random_seed = np.random.default_rng(81756234) + LOCS = np.c_[ + random_seed.random(4) * 0.6 + 0.2, + np.zeros(4), + random_seed.random(4) * 0.6 + 0.2, + ] def getError(self): funX = lambda x, y: np.cos(2 * np.pi * y) @@ -246,11 +250,11 @@ class TestInterpolationCyl(discretize.tests.OrderTest): meshTypes = ["uniformCylMesh", "randomCylMesh"] # MESHTYPES + meshDimension = 3 meshSizes = [8, 16, 32, 64] - rng = np.random.default_rng(876234) + random_seed = np.random.default_rng(876234) LOCS = np.c_[ - rng.random(20) * 0.6 + 0.2, - 2 * np.pi * (rng.random(20) * 0.6 + 0.2), - rng.random(20) * 0.6 + 0.2, + random_seed.random(20) * 0.6 + 0.2, + 2 * np.pi * (random_seed.random(20) * 0.6 + 0.2), + random_seed.random(20) * 0.6 + 0.2, ] def getError(self): @@ -314,9 +318,9 @@ def test_orderEz(self): class TestInterpolation3D(discretize.tests.OrderTest): - rng = np.random.default_rng(234) + random_seed = np.random.default_rng(234) name = "Interpolation" - LOCS = rng.random((50, 3)) * 0.6 + 0.2 + LOCS = random_seed.random((50, 3)) * 0.6 + 0.2 meshTypes = MESHTYPES tolerance = TOLERANCES meshDimension = 3 diff --git a/tests/base/test_tensor.py b/tests/base/test_tensor.py index 5d880f5bc..6a0fd9078 100644 --- a/tests/base/test_tensor.py +++ b/tests/base/test_tensor.py @@ -278,7 +278,6 @@ def test_cell_nodes(self, mesh): class TestPoissonEqn(discretize.tests.OrderTest): name = "Poisson Equation" meshSizes = [10, 16, 20] - rng = gen def getError(self): # Create some functions to integrate