Repository: CharJon/GeCO

Issue:
Avoid side effects when using random numbers.
Networkx decorator might be usefull. See e.g. [here](https://networkx.org/documentation/stable/reference/randomness.html).

Discussion:

Patch:
diff --git a/geco/mips/facility_location.py b/geco/mips/facility_location.py
index ec8c2d4..7139f3e 100644
--- a/geco/mips/facility_location.py
+++ b/geco/mips/facility_location.py
@@ -2,8 +2,10 @@ import itertools
 
 import numpy as np
 import pyscipopt as scip
+from networkx.utils import py_random_state
 
 
+@py_random_state(3)
 def capacitated_facility_location(n_customers, n_facilities, ratio, seed=0):
     """
     Generate a Capacited Facility Location problem following
@@ -19,23 +21,21 @@ def capacitated_facility_location(n_customers, n_facilities, ratio, seed=0):
         The desired number of facilities.
     ratio: float
         The desired capacity / demand ratio.
-    seed: int
-        The seed to use for random numbers.
+    seed: integer, random_state, or None
+        Indicator of random number generation state.
     """
-    rng = np.random.RandomState(seed)
-
     # locations for customers
-    c_x = rng.rand(n_customers)
-    c_y = rng.rand(n_customers)
+    c_x = np.array([seed.random() for _ in range(n_customers)])
+    c_y = np.array([seed.random() for _ in range(n_customers)])
 
     # locations for facilities
-    f_x = rng.rand(n_facilities)
-    f_y = rng.rand(n_facilities)
+    f_x = np.array([seed.random() for _ in range(n_facilities)])
+    f_y = np.array([seed.random() for _ in range(n_facilities)])
 
-    demands = rng.randint(5, 35 + 1, size=n_customers)
-    capacities = rng.randint(10, 160 + 1, size=n_facilities)
-    fixed_costs = rng.randint(100, 110 + 1, size=n_facilities) * np.sqrt(capacities) \
-                  + rng.randint(90 + 1, size=n_facilities)
+    demands = np.array(seed.sample(range(5, 35 + 1), k=n_customers))
+    capacities = np.array(seed.sample(range(10, 160 + 1), k=n_facilities))
+    fixed_costs = np.array(seed.sample(range(100, 110 + 1), k=n_facilities) * np.sqrt(capacities)) \
+                  + np.array(seed.sample(range(90 + 1), k=n_facilities))
     fixed_costs = fixed_costs.astype(int)
 
     total_demand = demands.sum()
diff --git a/geco/mips/knapsack.py b/geco/mips/knapsack.py
index 2f889ae..9b221a7 100644
--- a/geco/mips/knapsack.py
+++ b/geco/mips/knapsack.py
@@ -1,20 +1,20 @@
-import random
 import math
 
 import pyscipopt as scip
+from networkx.utils import py_random_state
 
 
-def yang(n, seed):
+@py_random_state(1)
+def yang(n, seed=0):
     """
     Generates knapsack instance parameters according to:
         Yu Yang, Natashia Boland, Bistra Dilkina, Martin Savelsbergh,
         "Learning Generalized Strong Branching for Set Covering,
         Set Packing, and 0-1 Knapsack Problems", 2020.
     """
-    random.seed(seed)
 
     def draw_value():
-        return random.randint(1, 10 * n)
+        return seed.randint(1, 10 * n)
 
     profits = [draw_value() for _ in range(n)]
     weights = [draw_value() for _ in range(n)]
diff --git a/geco/mips/scheduling.py b/geco/mips/scheduling.py
index a12ec87..ce15a0b 100644
--- a/geco/mips/scheduling.py
+++ b/geco/mips/scheduling.py
@@ -3,8 +3,9 @@ This module implements the scheduling problem MIP generation techniques from dif
 """
 
 import itertools
+
 import pyscipopt as scip
-import random
+from networkx.utils import py_random_state
 
 
 def hooker_late_tasks_formulation(
@@ -15,8 +16,7 @@ def hooker_late_tasks_formulation(
         C,
         c,
         r,
-        d,
-        seed=0,
+        d
 ):
     # TODO: use more expressive param names
     """Generates late tasks mip formulation described in section 4 in
@@ -25,13 +25,11 @@ def hooker_late_tasks_formulation(
     number_of_facilities: the number of facilities to schedule on
     number_of_tasks: the number of tasks to assign to facilities
     time_steps: the number of time steps starting from 0 (corresponds to "N" in the paper)
-    seed: used for randomization
     Other parameters follow the same naming used in the paper
 
     Returns:
         model: SCIP model of the late tasks instance
     """
-    random.seed(seed)
     model = scip.Model("Hooker Scheduling Late Tasks Formulation")
     assert min(r) == 0  # TODO: handle the case that timesteps don't start at 0
 
@@ -103,21 +101,21 @@ def generate_hookers_instances():
         )
 
 
-def generate_params(number_of_facilities, number_of_tasks, seed):
-    random.seed(seed)
+@py_random_state(2)
+def generate_params(number_of_facilities, number_of_tasks, seed=0):
     p = {}
 
     for j, i in itertools.product(range(number_of_tasks), range(number_of_facilities)):
         if number_of_tasks < 22:
-            p[j, i] = random.randint(2, 20 + 5 * i)
+            p[j, i] = seed.randint(2, 20 + 5 * i)
         else:
-            p[j, i] = random.randint(5, 20 + 5 * i)
+            p[j, i] = seed.randint(5, 20 + 5 * i)
 
     C = [10] * number_of_facilities
 
     c = {}
     for i in range(number_of_facilities):
-        value = random.randint(1, 10)
+        value = seed.randint(1, 10)
         for j in range(number_of_tasks):
             c[j, i] = value
 
@@ -126,11 +124,11 @@ def generate_params(number_of_facilities, number_of_tasks, seed):
     d = {}
     beta = 20 / 9
     for j in range(number_of_tasks):
-        d[j] = random.uniform(beta * number_of_tasks / 4, beta * number_of_tasks)
+        d[j] = seed.uniform(beta * number_of_tasks / 4, beta * number_of_tasks)
 
     r = {}
     for j, k in itertools.product(range(number_of_tasks), range(number_of_facilities)):
-        r[j, k] = random.randint(1, 9)
+        r[j, k] = seed.randint(1, 9)
 
     return p, C, c, R, d, r
 
@@ -143,8 +141,7 @@ def heinz_formulation(
         c,
         R,
         d,
-        r,
-        seed=0,
+        r
 ):
     """Generates mip formulation according to Model 4 in
     # TODO: Add paper reference
@@ -152,12 +149,10 @@ def heinz_formulation(
     number_of_facilities: the number of facilities to schedule on
     number_of_tasks: the number of tasks to assign to facilities
     time_steps: the number of time steps starting from 0 (corresponds to "N" in the paper)
-    seed: used for randomization
 
     Returns:
         model: SCIP model of the late tasks instance
     """
-    random.seed(seed)
     model = scip.Model("Heinz Scheduling")
     time_steps = range(min(R), int(max(d.values())))
 


Test patch:
diff --git a/geco/mips/tests/test_facility_location.py b/geco/mips/tests/test_facility_location.py
index f404678..d21b2b4 100644
--- a/geco/mips/tests/test_facility_location.py
+++ b/geco/mips/tests/test_facility_location.py
@@ -2,10 +2,10 @@ from geco.mips.facility_location import *
 
 
 def test_capacitated_facility_location():
-    n_customers, n_facilities, ratio = 25, 10, 2
-    m_1 = capacitated_facility_location(n_customers, n_facilities, ratio)
-    assert m_1.getNVars() == n_customers * n_facilities + n_facilities
-    assert m_1.getNConss() == n_customers + n_facilities + 1 + n_customers * n_facilities
-    assert m_1.getObjectiveSense() == "minimize"
-    m_1.optimize()
-    assert 5679 <= m_1.getObjVal() <= 5680
+    n_customers, n_facilities, ratio, seed = 25, 10, 2, 0
+    model = capacitated_facility_location(n_customers, n_facilities, ratio, seed)
+    assert model.getNVars() == n_customers * n_facilities + n_facilities
+    assert model.getNConss() == n_customers + n_facilities + 1 + n_customers * n_facilities
+    assert model.getObjectiveSense() == "minimize"
+    model.optimize()
+    assert 5856 <= model.getObjVal() <= 5857