v0.0.23 Icremental score calculator + examples with it (NQueens, TSP, VRP)

Author: CameleoGrey

examples/object_oriented/nqueens/persistence/CotwinBuilderNQueens.py

@@ -3,13 +3,15 @@
 
 from greyjack.persistence.CotwinBuilderBase import CotwinBuilderBase
 from examples.object_oriented.nqueens.cotwin.NQueensCotwin import NQueensCotwin
-from examples.object_oriented.nqueens.score.NQueensScoreCalculator import NQueensScoreCalculator
+from examples.object_oriented.nqueens.score.PlainScoreCalculatorNQueens import PlainScoreCalculatorNQueens
+from examples.object_oriented.nqueens.score.IncrementalScoreCalculatorNQueens import IncrementalScoreCalculatorNQueens
 from examples.object_oriented.nqueens.cotwin.CotQueen import CotQueen
 from greyjack.variables.GJInteger import GJInteger
 
 
 class CotwinBuilderNQueens(CotwinBuilderBase):
-    def __init__(self):
+    def __init__(self, use_incremental_score_calculator):
+        self.use_incremental_score_calculator = use_incremental_score_calculator
         pass
 
     def build_cotwin(self, domain_model, is_already_initialized):
@@ -27,7 +29,10 @@ def build_cotwin(self, domain_model, is_already_initialized):
 
         nqueens_cotwin = NQueensCotwin()
         nqueens_cotwin.add_planning_entities_list( cot_queens, "queens" )
-        nqueens_cotwin.set_score_calculator( NQueensScoreCalculator() )
+        if self.use_incremental_score_calculator:
+            nqueens_cotwin.set_score_calculator( IncrementalScoreCalculatorNQueens() )
+        else:
+            nqueens_cotwin.set_score_calculator( PlainScoreCalculatorNQueens() )
 
         return nqueens_cotwin
 

examples/object_oriented/nqueens/score/IncrementalScoreCalculatorNQueens.py

@@ -0,0 +1,115 @@
+
+from greyjack.score_calculation.score_calculators.IncrementalScoreCalculator import IncrementalScoreCalculator
+from greyjack.score_calculation.scores.SimpleScore import SimpleScore
+from greyjack.score_calculation.scores.ScoreVariants import ScoreVariants
+import polars as pl
+import numpy as np
+import numba
+from numba import jit, int64, float64, vectorize
+
+
+class IncrementalScoreCalculatorNQueens(IncrementalScoreCalculator):
+    def __init__(self):
+        super().__init__()
+        self.score_variant = ScoreVariants.SimpleScore
+        self.add_constraint("all_different", self.all_different)
+        pass
+
+    def all_different(self, planning_entity_dfs, problem_fact_dfs, delta_dfs):
+        # pseudo-incremental manner: just use deltas to reduce planning dfs building cost
+        # the score still computes on the whole vector sample
+        # relatively easy implement and understand (comparable to true incremental or super optimized manner) 
+        # with huge performance boost on large datasets
+
+        queens_df = planning_entity_dfs["queens"]
+        queens_deltas_df = delta_dfs["queens"]
+
+        native_column_ids = queens_df["column_id"].to_numpy()
+        native_row_ids = queens_df["row_id"].to_numpy()
+        target_unique_row_ids_count = native_row_ids.shape[0]
+
+        scores = []
+        splitted_dfs = queens_deltas_df.partition_by(["sample_id"], maintain_order=True, include_key=False)
+        for sample_df in splitted_dfs:
+            delta_df_row_ids = sample_df["candidate_df_row_id"].to_numpy()
+            delta_row_ids = sample_df["row_id"].to_numpy()
+            sample_row_ids = native_row_ids.copy()
+
+            sample_row_ids[delta_df_row_ids] = delta_row_ids
+
+            unique_sample_row_ids = np.unique(sample_row_ids)
+            unique_sample_desc_ids = np.unique(native_column_ids + sample_row_ids)
+            unique_sample_asc_ids = np.unique(native_column_ids - sample_row_ids)
+
+            unique_rows_penalty = target_unique_row_ids_count - unique_sample_row_ids.shape[0]
+            unique_desc_ids_penalty = target_unique_row_ids_count - unique_sample_desc_ids.shape[0]
+            unique_asc_ids_penalty = target_unique_row_ids_count - unique_sample_asc_ids.shape[0]
+
+            total_penalty = unique_rows_penalty + unique_desc_ids_penalty + unique_asc_ids_penalty
+            scores.append(SimpleScore(total_penalty))
+
+        return scores
+    
+    """def all_different(self, planning_entity_dfs, problem_fact_dfs, delta_dfs):
+        # alternative (much more complex and faster) version
+        # after experiments: better works on medium/low dimensions (especially with jit)
+
+        queens_df = planning_entity_dfs["queens"]
+        queens_deltas_df = delta_dfs["queens"]
+
+        native_column_ids = queens_df["column_id"].to_numpy()
+        native_row_ids = queens_df["row_id"].to_numpy()
+        target_unique_row_ids_count = native_row_ids.shape[0]
+
+        n_samples = queens_deltas_df["sample_id"].n_unique() # fast call to Rust function
+        sample_array = queens_deltas_df.select(["sample_id", "candidate_df_row_id", "row_id"]).to_numpy() # guarantee of order
+        scores = compute_score_values(
+            native_column_ids, 
+            native_row_ids, 
+            target_unique_row_ids_count, 
+            sample_array[:, 0],
+            sample_array[:, 1], 
+            sample_array[:, 2], 
+            n_samples
+            ) # jitted code
+        scores = [SimpleScore(score_value) for score_value in scores]
+
+        return scores"""
+    
+@staticmethod
+@jit(nopython=True)
+#@numba.cfunc("float64[:](int64[:], int64[:], int64, int64[:], int64[:], int64[:], int64)")
+def compute_score_values(
+    native_column_ids, 
+    native_row_ids, 
+    target_unique_row_ids_count, 
+    sample_ids,
+    delta_df_row_ids,
+    delta_row_ids,
+    n_samples
+    ):
+
+    score_values = np.zeros((n_samples, ), dtype=np.float64)
+    sample_array_id = 0
+    batch_start = sample_array_id
+    for sample_id in range(n_samples):
+        sample_row_ids = native_row_ids.copy()
+
+        while sample_ids[batch_start] == sample_ids[sample_array_id]:
+            sample_row_ids[delta_df_row_ids[sample_array_id]] = delta_row_ids[sample_array_id]
+            sample_array_id += 1
+            if sample_array_id >= sample_ids.shape[0]:
+                break
+        batch_start = sample_array_id
+        
+        unique_sample_row_ids = np.unique(sample_row_ids)
+        unique_sample_desc_ids = np.unique(native_column_ids + sample_row_ids)
+        unique_sample_asc_ids = np.unique(native_column_ids - sample_row_ids)
+
+        unique_rows_penalty = target_unique_row_ids_count - unique_sample_row_ids.shape[0]
+        unique_desc_ids_penalty = target_unique_row_ids_count - unique_sample_desc_ids.shape[0]
+        unique_asc_ids_penalty = target_unique_row_ids_count - unique_sample_asc_ids.shape[0]
+
+        score_values[sample_id] = unique_rows_penalty + unique_desc_ids_penalty + unique_asc_ids_penalty
+
+    return score_values

examples/object_oriented/nqueens/score/NQueensScoreCalculator.py renamed to examples/object_oriented/nqueens/score/PlainScoreCalculatorNQueens.py

@@ -5,12 +5,8 @@
 from greyjack.score_calculation.scores.ScoreVariants import ScoreVariants
 import polars as pl
 
-#import os
-#os.environ["POLARS_MAX_THREADS"] = "1"
-#os.environ["MIMALLOC_ABANDONED_PAGE_RESET"] = "1"
 
-
-class NQueensScoreCalculator(PlainScoreCalculator):
+class PlainScoreCalculatorNQueens(PlainScoreCalculator):
     def __init__(self):
 
         super().__init__()

examples/object_oriented/nqueens/scripts/solve_nqueens.py

@@ -22,16 +22,16 @@
 if __name__ == "__main__":
 
     # build domain model
-    domain_builder = DomainBuilderNQueens(256, random_seed=45)
-    cotwin_builder = CotwinBuilderNQueens()
+    domain_builder = DomainBuilderNQueens(1024, random_seed=45)
+    cotwin_builder = CotwinBuilderNQueens(use_incremental_score_calculator=True)
 
     #termination_strategy = StepsLimit(step_count_limit=1000)
-    termination_strategy = TimeSpentLimit(time_seconds_limit=60)
+    #termination_strategy = TimeSpentLimit(time_seconds_limit=60)
     #termination_strategy = ScoreNoImprovement(time_seconds_limit=15)
-    #termination_strategy = ScoreLimit(score_to_compare=[0])
+    termination_strategy = ScoreLimit(score_to_compare=[0])
     agent = TabuSearch(neighbours_count=128, tabu_entity_rate=0.0, 
-                       mutation_rate_multiplier=None, move_probas=None, 
-                       migration_frequency=10, termination_strategy=termination_strategy)
+                       mutation_rate_multiplier=None, move_probas=[0, 1, 0, 0, 0, 0], 
+                       migration_frequency=1, termination_strategy=termination_strategy)
     """agent = GeneticAlgorithm(population_size=128, crossover_probability=0.5, p_best_rate=0.05,
                              tabu_entity_rate=0.0, mutation_rate_multiplier=1.0, move_probas=[0, 1, 0, 0, 0, 0],
                              migration_rate=0.00001, migration_frequency=1, termination_strategy=termination_strategy)"""

examples/object_oriented/nqueens/scripts/solve_nqueens_ortools.py

@@ -0,0 +1,77 @@
+"""OR-Tools solution to the N-queens problem."""
+import sys
+import time
+from ortools.sat.python import cp_model
+
+
+class NQueenSolutionPrinter(cp_model.CpSolverSolutionCallback):
+    """Print intermediate solutions."""
+
+    def __init__(self, queens: list[cp_model.IntVar]):
+        cp_model.CpSolverSolutionCallback.__init__(self)
+        self.__queens = queens
+        self.__solution_count = 0
+        self.__start_time = time.time()
+
+    @property
+    def solution_count(self) -> int:
+        return self.__solution_count
+
+    def on_solution_callback(self):
+        current_time = time.time()
+        print(
+            f"Solution {self.__solution_count}, "
+            f"time = {current_time - self.__start_time} s"
+        )
+        self.__solution_count += 1
+
+        """all_queens = range(len(self.__queens))
+        for i in all_queens:
+            for j in all_queens:
+                if self.value(self.__queens[j]) == i:
+                    # There is a queen in column j, row i.
+                    print("Q", end=" ")
+                else:
+                    print("_", end=" ")
+            print()
+        print()"""
+
+
+
+def main(board_size: int) -> None:
+    # Creates the solver.
+    model = cp_model.CpModel()
+
+    # Creates the variables.
+    # There are `board_size` number of variables, one for a queen in each column
+    # of the board. The value of each variable is the row that the queen is in.
+    queens = [model.new_int_var(0, board_size - 1, f"x_{i}") for i in range(board_size)]
+
+    # Creates the constraints.
+    # All rows must be different.
+    model.add_all_different(queens)
+
+    # No two queens can be on the same diagonal.
+    model.add_all_different(queens[i] + i for i in range(board_size))
+    model.add_all_different(queens[i] - i for i in range(board_size))
+
+    # Solve the model.
+    solver = cp_model.CpSolver()
+    solution_printer = NQueenSolutionPrinter(queens)
+    solver.parameters.enumerate_all_solutions = False
+    solver.solve(model, solution_printer)
+
+    # Statistics.
+    print("\nStatistics")
+    print(f"  conflicts      : {solver.num_conflicts}")
+    print(f"  branches       : {solver.num_branches}")
+    print(f"  wall time      : {solver.wall_time} s")
+    print(f"  solutions found: {solution_printer.solution_count}")
+
+
+if __name__ == "__main__":
+    # By default, solve the 8x8 problem.
+    size = 1024
+    if len(sys.argv) > 1:
+        size = int(sys.argv[1])
+    main(size)

examples/object_oriented/tsp/persistence/CotwinBuilder.py

@@ -4,18 +4,21 @@
 from examples.object_oriented.tsp.cotwin.CotStop import CotStop
 from examples.object_oriented.tsp.cotwin.EdgeDistance import EdgeDistance
 from examples.object_oriented.tsp.cotwin.TSPCotwin import TSPCotwin
-from examples.object_oriented.tsp.score.TSPPlainScoreCalculator import TSPPlainScoreCalculator
+from examples.object_oriented.tsp.score.PlainScoreCalculatorTSP import PlainScoreCalculatorTSP
+from examples.object_oriented.tsp.score.IncrementalScoreCalculatorTSP import IncrementalScoreCalculatorTSP
 
 from greyjack.persistence.CotwinBuilderBase import CotwinBuilderBase
 from greyjack.variables.GJInteger import GJInteger
 
 
 class CotwinBuilder(CotwinBuilderBase):
 
-    def __init__(self):
+    def __init__(self, use_incremental_score_calculator):
 
         super().__init__()
 
+        self.use_incremental_score_calculator = use_incremental_score_calculator
+
         pass
 
     def build_cotwin(self, domain_model, is_already_initialized):
@@ -27,7 +30,11 @@ def build_cotwin(self, domain_model, is_already_initialized):
         # (!) better use distance matrix mapping than building and joining large distance dataframe!
         #tsp_cotwin.add_problem_facts_list(self._build_edge_distances(domain_model), "edge_distances")
 
-        score_calculator = TSPPlainScoreCalculator()
+        if self.use_incremental_score_calculator:
+            score_calculator = PlainScoreCalculatorTSP()
+        else:
+            score_calculator = IncrementalScoreCalculatorTSP()
+
         score_calculator.utility_objects["distance_matrix"] = domain_model.distance_matrix
         tsp_cotwin.set_score_calculator( score_calculator )
 

examples/object_oriented/tsp/persistence/DomainBuilder.py

@@ -2,6 +2,7 @@
 
 import re
 import numpy as np
+from numba import jit
 from examples.object_oriented.tsp.domain.TravelSchedule import TravelSchedule
 from examples.object_oriented.tsp.domain.Vehicle import Vehicle
 from examples.object_oriented.tsp.domain.Location import Location
@@ -65,21 +66,26 @@ def build_from_solution(self, solution):
     def build_from_domain(self, domain):
         return super().build_from_domain(domain)
 
-
-
     def _build_distance_matrix(self, locations_list):
 
+        @staticmethod
+        @jit()
+        def compute_distance_matrix(latitudes, longitudes, n_locations):
+            distance_matrix = np.zeros( (n_locations, n_locations), dtype=np.int64 )
+            for i in range(n_locations):
+                for j in range(n_locations):
+                    distance_from_to = np.sqrt((latitudes[i] - latitudes[j])**2 + (longitudes[i] - longitudes[j])**2)
+                    distance_matrix[i][j] = round(1000 * distance_from_to, 0)
+            return distance_matrix
 
         n_locations = len(locations_list)
-        distance_matrix = np.zeros( (n_locations, n_locations), dtype=np.int64 )
-        #distance_matrix = np.zeros((n_locations, n_locations), dtype=np.float32)
-        for i in range(n_locations):
-            for j in range(n_locations):
-                location_from = locations_list[i]
-                location_to = locations_list[j]
-                distance_from_to = location_from.get_distance_to_other_location( location_to )
-                distance_matrix[i][j] = round(1000 * distance_from_to, 0)
-                #distance_matrix[i][j] = distance_from_to
+        latitudes = np.zeros((n_locations, ), dtype=np.float64)
+        longitudes = np.zeros((n_locations, ), dtype=np.float64)
+        for i, location in enumerate(locations_list):
+            latitudes[i] = location.latitude
+            longitudes[i] = location.longitude
+    
+        distance_matrix = compute_distance_matrix(latitudes, longitudes, n_locations)
 
         return distance_matrix