constraint-aware variants for pso, firefly, and sa

Author: benvanwerkhoven

kernel_tuner/strategies/common.py

@@ -3,6 +3,7 @@
 from time import perf_counter
 
 import numpy as np
+from scipy.spatial import distance
 
 from kernel_tuner import util
 from kernel_tuner.searchspace import Searchspace
@@ -88,8 +89,17 @@ def __call__(self, x, check_restrictions=True):
 
         # else check if this is a legal (non-restricted) configuration
         if check_restrictions and self.searchspace.restrictions:
+            legal = self.searchspace.is_param_config_valid(tuple(params))
             params_dict = dict(zip(self.searchspace.tune_params.keys(), params))
-            legal = util.check_restrictions(self.searchspace.restrictions, params_dict, self.tuning_options.verbose)
+
+            if "constraint_aware" in self.tuning_options.strategy_options and self.tuning_options.strategy_options["constraint_aware"]:
+                # attempt to repair
+                new_params = unscale_and_snap_to_nearest_valid(x, params, self.searchspace, self.tuning_options.eps)
+                if new_params:
+                    params = new_params
+                    legal = True
+                    x_int = ",".join([str(i) for i in params])
+
             if not legal:
                 result = params_dict
                 result[self.tuning_options.objective] = util.InvalidConfig()
@@ -243,3 +253,28 @@ def scale_from_params(params, tune_params, eps):
     for i, v in enumerate(tune_params.values()):
         x[i] = 0.5 * eps + v.index(params[i])*eps
     return x
+
+
+
+def unscale_and_snap_to_nearest_valid(x, params, searchspace, eps):
+    """Helper func to snap to the nearest valid configuration"""
+
+    # params is nearest unscaled point, but is not valid
+    neighbors = get_neighbors(params, searchspace)
+
+    if neighbors:
+        # sort on distance to x
+        neighbors.sort(key=lambda y: distance.euclidean(x,scale_from_params(y, searchspace.tune_params, eps)))
+
+        # return closest valid neighbor
+        return neighbors[0]
+
+    return []
+
+
+def get_neighbors(params, searchspace):
+    for neighbor_method in ["strictly-adjacent", "adjacent", "Hamming"]:
+        neighbors = searchspace.get_neighbors_no_cache(tuple(params), neighbor_method=neighbor_method)
+        if len(neighbors) > 0:
+            return neighbors
+    return []

kernel_tuner/strategies/firefly_algorithm.py

@@ -13,7 +13,8 @@
                        maxiter=("Maximum number of iterations", 100),
                        B0=("Maximum attractiveness", 1.0),
                        gamma=("Light absorption coefficient", 1.0),
-                       alpha=("Randomization parameter", 0.2))
+                       alpha=("Randomization parameter", 0.2),
+                       constraint_aware=("constraint-aware optimization (True/False)", True))
 
 def tune(searchspace: Searchspace, runner, tuning_options):
 
@@ -23,7 +24,7 @@ def tune(searchspace: Searchspace, runner, tuning_options):
     # using this instead of get_bounds because scaling is used
     bounds, _, eps = cost_func.get_bounds_x0_eps()
 
-    num_particles, maxiter, B0, gamma, alpha = common.get_options(tuning_options.strategy_options, _options)
+    num_particles, maxiter, B0, gamma, alpha, constraint_aware = common.get_options(tuning_options.strategy_options, _options)
 
     best_score_global = sys.float_info.max
     best_position_global = []
@@ -34,9 +35,10 @@ def tune(searchspace: Searchspace, runner, tuning_options):
         swarm.append(Firefly(bounds))
 
     # ensure particles start from legal points
-    population = list(list(p) for p in searchspace.get_random_sample(num_particles))
-    for i, particle in enumerate(swarm):
-        particle.position = scale_from_params(population[i], searchspace.tune_params, eps)
+    if constraint_aware:
+        population = list(list(p) for p in searchspace.get_random_sample(num_particles))
+        for i, particle in enumerate(swarm):
+            particle.position = scale_from_params(population[i], searchspace.tune_params, eps)
 
     # compute initial intensities
     for j in range(num_particles):

kernel_tuner/strategies/pso.py

@@ -13,7 +13,8 @@
                        maxiter=("Maximum number of iterations", 100),
                        w=("Inertia weight constant", 0.5),
                        c1=("Cognitive constant", 2.0),
-                       c2=("Social constant", 1.0))
+                       c2=("Social constant", 1.0),
+                       constraint_aware=("constraint-aware optimization (True/False)", False))
 
 def tune(searchspace: Searchspace, runner, tuning_options):
 
@@ -24,7 +25,7 @@ def tune(searchspace: Searchspace, runner, tuning_options):
     bounds, _, eps = cost_func.get_bounds_x0_eps()
 
 
-    num_particles, maxiter, w, c1, c2 = common.get_options(tuning_options.strategy_options, _options)
+    num_particles, maxiter, w, c1, c2, constraint_aware = common.get_options(tuning_options.strategy_options, _options)
 
     best_score_global = sys.float_info.max
     best_position_global = []
@@ -35,9 +36,10 @@ def tune(searchspace: Searchspace, runner, tuning_options):
         swarm.append(Particle(bounds))
 
     # ensure particles start from legal points
-    population = list(list(p) for p in searchspace.get_random_sample(num_particles))
-    for i, particle in enumerate(swarm):
-        particle.position = scale_from_params(population[i], searchspace.tune_params, eps)
+    if constraint_aware:
+        population = list(list(p) for p in searchspace.get_random_sample(num_particles))
+        for i, particle in enumerate(swarm):
+            particle.position = scale_from_params(population[i], searchspace.tune_params, eps)
 
     # start optimization
     for i in range(maxiter):

kernel_tuner/strategies/simulated_annealing.py

@@ -10,16 +10,17 @@
 from kernel_tuner.strategies.common import CostFunc
 
 _options = dict(T=("Starting temperature", 1.0),
-                       T_min=("End temperature", 0.001),
-                       alpha=("Alpha parameter", 0.995),
-                       maxiter=("Number of iterations within each annealing step", 1))
+                T_min=("End temperature", 0.001),
+                alpha=("Alpha parameter", 0.995),
+                maxiter=("Number of iterations within each annealing step", 1),
+                constraint_aware=("constraint-aware optimization (True/False)", True))
 
 def tune(searchspace: Searchspace, runner, tuning_options):
     # SA works with real parameter values and does not need scaling
     cost_func = CostFunc(searchspace, tuning_options, runner)
 
     # optimization parameters
-    T, T_min, alpha, niter = common.get_options(tuning_options.strategy_options, _options)
+    T, T_min, alpha, niter, constraint_aware = common.get_options(tuning_options.strategy_options, _options)
     T_start = T
 
     # compute how many iterations would be needed to complete the annealing schedule
@@ -30,7 +31,7 @@ def tune(searchspace: Searchspace, runner, tuning_options):
     max_feval = tuning_options.strategy_options.get("max_fevals", max_iter)
 
     # get random starting point and evaluate cost
-    pos = list(searchspace.get_random_sample(1)[0])
+    pos = generate_starting_point(searchspace, constraint_aware)
     old_cost = cost_func(pos, check_restrictions=False)
 
     # main optimization loop
@@ -46,9 +47,9 @@ def tune(searchspace: Searchspace, runner, tuning_options):
 
         for _ in range(niter):
 
-            new_pos = neighbor(pos, searchspace)
+            new_pos = neighbor(pos, searchspace, constraint_aware)
             try:
-                new_cost = cost_func(new_pos, check_restrictions=False)
+                new_cost = cost_func(new_pos, check_restrictions=not constraint_aware)
             except util.StopCriterionReached as e:
                 if tuning_options.verbose:
                     print(e)
@@ -73,7 +74,7 @@ def tune(searchspace: Searchspace, runner, tuning_options):
             stuck = 0
         c_old = c
         if stuck > 100:
-            pos = list(searchspace.get_random_sample(1)[0])
+            pos = generate_starting_point(searchspace, constraint_aware)
             stuck = 0
 
         # safeguard
@@ -103,11 +104,49 @@ def acceptance_prob(old_cost, new_cost, T, tuning_options):
     return np.exp(((old_cost-new_cost)/old_cost)/T)
 
 
-def neighbor(pos, searchspace: Searchspace):
+def neighbor(pos, searchspace: Searchspace, constraint_aware=True):
     """Return a random neighbor of pos."""
-    # Note: this is not the same as the previous implementation, because it is possible that non-edge parameters remain the same, but suggested configurations will all be within restrictions
-    neighbors = searchspace.get_neighbors(tuple(pos), neighbor_method='Hamming') if random.random() < 0.2 else searchspace.get_neighbors(tuple(pos), neighbor_method='strictly-adjacent')
-    if len(neighbors) > 0:
-        return list(random.choice(neighbors))
-    # if there are no neighbors, return a random configuration
-    return list(searchspace.get_random_sample(1)[0])
+
+    if constraint_aware:
+        # Note: this is not the same as the previous implementation, because it is possible that non-edge parameters remain the same, but suggested configurations will all be within restrictions
+        neighbors = searchspace.get_neighbors(tuple(pos), neighbor_method='Hamming') if random.random() < 0.2 else searchspace.get_neighbors(tuple(pos), neighbor_method='strictly-adjacent')
+        if len(neighbors) > 0:
+            return list(random.choice(neighbors))
+        # if there are no neighbors, return a random configuration
+        return list(searchspace.get_random_sample(1)[0])
+
+    else:
+        tune_params = searchspace.tune_params
+        size = len(pos)
+        pos_out = []
+        # random mutation
+        # expected value is set that values all dimensions attempt to get mutated
+        for i in range(size):
+            key = list(tune_params.keys())[i]
+            values = tune_params[key]
+
+            if random.random() < 0.2:  #replace with random value
+                new_value = random_val(i, tune_params)
+            else: #adjacent value
+                ind = values.index(pos[i])
+                if random.random() > 0.5:
+                    ind += 1
+                else:
+                    ind -= 1
+                ind = min(max(ind, 0), len(values)-1)
+                new_value = values[ind]
+
+            pos_out.append(new_value)
+        return pos_out
+
+def random_val(index, tune_params):
+    """return a random value for a parameter"""
+    key = list(tune_params.keys())[index]
+    return random.choice(tune_params[key])
+
+def generate_starting_point(searchspace: Searchspace, constraint_aware=True):
+    if constraint_aware:
+        return list(searchspace.get_random_sample(1)[0])
+    else:
+        tune_params = searchspace.tune_params
+        return [random_val(i, tune_params) for i in range(len(tune_params))]