Add LHS expansion method (SMTorg#303)

* New function in lsh + test_lhs

* Commentaries to lhs new functionality

* Minor changes within the lhs testing file

* Update mfk.py

* Update lhs.py

* Blacked functions

* Errors fixed for the PR

* Errors fixed for the PR

* Changes in the boolean if condition

Co-authored-by: Ruben CONDE-ARENZANA <[email protected]>

smt/applications/mfk.py

@@ -90,7 +90,7 @@ def __call__(self, nb_samples_hifi):
                 xlimits=self.xlimits, criterion="ese", random_state=self.random_state
             )
             doe.append(p(nt[i]))
-
+               
         for i in range(1, self.nlevel)[::-1]:
             ind = []
             d = cdist(doe[i], doe[i - 1], "euclidean")

smt/sampling_methods/lhs.py

@@ -323,3 +323,88 @@ def _ese(self, dim, nt):
             return_hist=True,
         )
         return P
+
+    def expand_lhs(self, x, n_points, method="basic"):
+        """
+        Given a Latin Hypercube Sample (LHS) "x", returns an expanded LHS
+        by adding "n_points" new points.
+
+        Parameters
+        ----------
+        x : array
+            Initial LHS.
+        n_points : integer
+            Number of points that are to be added to the expanded LHS.
+        method : str, optional
+            Methodoly for the construction of the expanded LHS.
+            The default is "basic".
+
+        Returns
+        -------
+        x_new : array
+            Expanded LHS.
+
+        """
+
+        xlimits = self.options["xlimits"]
+
+        new_num = len(x) + n_points
+        if new_num % len(x) != 0:
+            print(
+                "WARNING: The added number of points is not a "
+                "multiple of the initial number of points."
+                "Thus, it cannot be ensured that the output is an LHS."
+            )
+
+        # For future functionalities
+        if method == "basic":
+
+            # Evenly spaced intervals with the final dimension of the LHS
+            intervals = []
+            for i in range(len(xlimits)):
+                intervals.append(np.linspace(xlimits[i][0], xlimits[i][1], new_num + 1))
+
+            # Creates a subspace with the rows and columns that have no points
+            # in the new space
+            subspace_limits = [[]] * len(xlimits)
+            subspace_bool = []
+            for i in range(len(xlimits)):
+                subspace_limits[i] = []
+
+                subspace_bool.append(
+                    [
+                        [
+                            intervals[i][j] < x[kk][i] < intervals[i][j + 1]
+                            for kk in range(len(x))
+                        ]
+                        for j in range(len(intervals[i]) - 1)
+                    ]
+                )
+
+                [
+                    subspace_limits[i].append([intervals[i][ii], intervals[i][ii + 1]])
+                    for ii in range(len(subspace_bool[i]))
+                    if not (True in subspace_bool[i][ii])
+                ]
+
+            # Sampling of the new subspace
+            sampling_new = LHS(xlimits=np.array([[0.0, 1.0]] * len(xlimits)))
+            x_subspace = sampling_new(n_points)
+
+            columnIndex = 0
+            sortedArr = x_subspace[x_subspace[:, columnIndex].argsort()]
+
+            for j in range(len(xlimits)):
+                for i in range(len(sortedArr)):
+                    sortedArr[i, j] = subspace_limits[j][i][0] + sortedArr[i, j] * (
+                        subspace_limits[j][i][1] - subspace_limits[j][i][0]
+                    )
+
+            H = np.zeros_like(sortedArr)
+            for j in range(len(xlimits)):
+                order = np.random.permutation(len(sortedArr))
+                H[:, j] = sortedArr[order, j]
+
+            x_new = np.concatenate((x, H), axis=0)
+
+        return x_new

smt/sampling_methods/tests/test_lhs.py

@@ -30,6 +30,47 @@ def test_random_state(self):
         doe2 = sampling(num)
         self.assertTrue(np.allclose(doe1, doe2))
 
+    def test_expand_lhs(self):
+        import numpy as np
+
+        num = 100
+        new_list = np.linspace(1, 5, 5) * num
+
+        for i in range(len(new_list)):
+            xlimits = np.array([[0.0, 4.0], [0.0, 3.0], [0.0, 3.0], [1.0, 5.0]])
+            sampling = LHS(xlimits=xlimits, criterion="ese")
+
+            x = sampling(num)
+            new = int(new_list[i])
+            new_num = num + new
+
+            x_new = sampling.expand_lhs(x, new)
+
+            intervals = []
+            subspace_bool = []
+            for i in range(len(xlimits)):
+                intervals.append(np.linspace(xlimits[i][0], xlimits[i][1], new_num + 1))
+
+                subspace_bool.append(
+                    [
+                        [
+                            intervals[i][j] < x_new[kk][i] < intervals[i][j + 1]
+                            for kk in range(len(x_new))
+                        ]
+                        for j in range(len(intervals[i]) - 1)
+                    ]
+                )
+
+                self.assertEqual(
+                    True,
+                    all(
+                        [
+                            subspace_bool[i][k].count(True) == 1
+                            for k in range(len(subspace_bool[i]))
+                        ]
+                    ),
+                )
+
 
 if __name__ == "__main__":
     unittest.main()