Merge pull request #61 from ed-rhilbert/feature/itinerary_setup

Feature/itinerary setup

docs/model.md

@@ -36,6 +36,7 @@ The following data are retrieved from the results of an OSRD simulation.
 - $train_s \in \{1,..,N_{trains}\}, \forall s \in \{1,..,N_{steps}\}$ : The train associated with the step $s$
 - $zone_s \in \{1,..,N_{zones}\}, \forall s \in \{1,..,N_{steps}\}$ : The zone associated with the step $s$
 - $prev_s \in \{0..N_{steps}\}, \forall s \in \{1,..,N_{steps}\}$ : The step preceding $s$ (0 if None)
+- $next_s \in \{0..N_{steps}\}, \forall s \in \{1,..,N_{steps}\}$ : The step following $s$ (0 if None)
 
 **Arrival and departure times**
 - $min\_arrival_s \in \mathbb{N}, \forall s \in \{1,..,N_{steps}\}$ : The min arrival time of a step
@@ -45,11 +46,27 @@ The following data are retrieved from the results of an OSRD simulation.
 **Constraint depending on the type of zone (lane, block, switch...)**
 - $is\_fixed_s \in \{True, False\}, \forall s \in \{1,..,N_{steps}\}$ : If the duration of a step cannot be changed
 
+**Constants**
+- $T_c$ : The time setup to change between to different itineraries
+
 ## Decision variables
 
 - $arrival_s \in \mathbb{N}, \forall s \in \{1,..,N_{steps}\}$ : The arrival time of a step
 - $departure_s \in \mathbb{N}, \forall s \in \{1,..,N_{steps}\}$ : The departure time of a step
 
+## Intermediate variables
+
+- $prec_i^j \in \{0, 1\}, \forall i \in \{1,..,N_{steps}\}, \forall j \in \{1,..,N_{steps}\}$ : $1$ if step $j$ directly follows step $i$ on $zone_i$. It follows the basic structural constraint :
+$$
+zone_i \ne zone_j \text{ or } train_j = train_i \implies prec_i^j = 0
+$$
+- $first_s \in \{0, 1\}, \forall s \in \{1,..,N_{steps}\}$ : $1$ if the step is the first on $zone_s$.
+- $last_s \in \{0, 1\}, \forall s \in \{1,..,N_{steps}\}$ : $1$ if the step is the last on $zone_s$.
+- $diff\_it_i ^j \in \{0, 1\}, \forall i \in \{1,..,N_{steps}\}, \forall j \in \{1,..,N_{steps}\}$ : $1$ if and only if the steps $i$ and $j$ are followed by seps with different zone. This can be expressed like this :
+$$
+diff\_it_i^j = 1 \Leftrightarrow zone_{next_i} \ne zone_{next_j}
+$$
+
 ## Objective
 
 We minimize the sum of arrival delays
@@ -94,7 +111,43 @@ $$\forall s \in \{1,..,N_{steps}\} \; s.t. \; prev_s \neq 0, \\arrival_s = depar
 
 $$\forall s \in \{1,..,N_{steps}\} \; s.t. \; prev_s = 0, arrival_s = min\_arrival_s$$
 
-7. _OPTIONAL_ A train cannot overtake another train (we use the reference ($min\_arrival$) to determine the order of trains)
+### Enforce order (optional)
+
+7. A train cannot overtake another train (we use the reference ($min\_arrival$) to determine the order of trains)
 
 $$\forall s1, s2 \in \{1,...,N_{steps}\}\;s.t.\;  min\_arrival_{s1} < min\_arrival_{s2},\\
 arrival_{s1} < arrival_{s2}$$
+
+### Precedence constraints
+
+8. There can only be one step being the first in a zone
+$$
+\forall z \in N_{zones}\; \sum_{s \in \{zone_s = z\}}first_s=1
+$$
+
+9. There can only be one step being the last in a zone
+$$
+\forall z \in N_{zones}\; \sum_{s \in \{zone_s = z\}}last_s=1
+$$
+
+10. If a step $i$ follows a step $j$ then the revese cannot be true
+$$
+prec_i^j + prec_j^i \le 1
+$$
+
+11. For every step $i$ the following step $j$ must starts after the end of step $j$. An itinerary setup is added if the
+itineraries are different on the next steps of each $i$ and $j$.
+
+$$
+\forall i \in {1,..., N_{steps}}\;\forall j \in {1,..., N_{steps}}\;\\prec_i^j \implies departure_{i} + T_c \times diff\_it_i ^j \le arrival_{j}
+$$
+
+12. For every step there must be exactly one following step except if it's the last step on its zone
+$$
+\forall i \in {1,..., N_{steps}}\;\sum_{j\in \{j, zone_j = zone_i\}} prec_i^j + last_i = 1
+$$
+
+13. For every step there must be exactly one previous step except if it's the first step on its zone
+$$
+\forall i \in {1,..., N_{steps}}\;\sum_{j\in \{j, zone_j = zone_i\}} prec_j^i + first_i = 1
+$$

src/cpagent/constraints.py

@@ -15,6 +15,7 @@ def _create_constraints(
     self._add_chaining_constraints(model)
     if not self.allow_change_order:
         self._add_enforce_order_constraints(model)
+    self._add_precedence_constraints(model)
 
 
 def _add_spacing_constraints(
@@ -65,3 +66,55 @@ def _add_enforce_order_constraints(
                 and step['zone'] == other['zone']
             ):
                 model.Add(self.t_in[i] < self.t_in[j])
+
+
+def _add_precedence_constraints(
+        self,
+        model: cp_model.CpModel
+) -> None:
+    """Ensure that the precedence between trains is well respected.
+
+    Multiple constraints are created to ensure we can know which steps
+    follow directly which other steps.
+
+    This corresponds ton constraints from 8 to 13 in the model
+    Parameters
+    ----------
+    model : cp_model.CpModel
+        model to fill
+    """
+    # build a map of step per zone
+    step_per_zone = {}
+    for step in self.steps:
+        if not step['zone'] in step_per_zone:
+            step_per_zone[step['zone']] = []
+        step_per_zone[step['zone']].append(step)
+
+    # Constraints 8 and 9 from the model
+    for steps_of_zone in step_per_zone.values():
+        model.AddExactlyOne(self.firsts[step['idx']] for step in steps_of_zone)
+        model.AddExactlyOne(self.lasts[step['idx']] for step in steps_of_zone)
+
+    for i, step in enumerate(self.steps):
+        all_others_before = [self.lasts[i]]
+        all_others_after = [self.firsts[i]]
+        for j, other in enumerate(self.steps):
+            if (
+                step['train'] != other['train']
+                and step['zone'] == other['zone']
+            ):
+                all_others_before.append(self.precs[i][j])
+                all_others_after.append(self.precs[j][i])
+                # Constraint 10
+                model.AddAtMostOne([self.precs[i][j], self.precs[j][i]])
+                # Constraint 11
+                model.Add(
+                    self.t_out[i]
+                    + self.itinierary_setup * self.diff_itineraries[i][j]
+                    <= self.t_in[j]) \
+                    .OnlyEnforceIf(self.precs[i][j])
+
+        # Constraint 12
+        model.AddExactlyOne(all_others_before)
+        # Constraint 13
+        model.AddExactlyOne(all_others_after)

src/cpagent/cp_agent.py

@@ -23,6 +23,7 @@ class CpAgent(SchedulerAgent):
         _add_spacing_constraints,
         _add_chaining_constraints,
         _add_enforce_order_constraints,
+        _add_precedence_constraints,
         _create_constraints  # must be last because it calls above methods
     )
     from .variables import _create_variables
@@ -42,6 +43,7 @@ class CpAgent(SchedulerAgent):
     extra_delays = None
     max_optimization_time = SOLVER_TIMEOUT
     save_history = False
+    itinierary_setup = 120
 
     # solution
 

src/cpagent/schedule_adapters.py

@@ -18,13 +18,15 @@ class OptimisationStatus(Enum):
     FAILED = 3
 
 
-def build_step(train: str, zone: int, prev_idx: int, min_t_in: int,
+def build_step(idx: int, train: str, zone: int, prev_idx: int, min_t_in: int,
                min_t_out: int, min_duration: int, is_fixed: bool,
-               ponderation: int = 1, overlap: int = 0) -> dict:
+               ponderation: int = 1, overlap: int = 0, next: int = -1) -> dict:
     """Add a step to the regulation problem
 
     Parameters
     ----------
+    idx : int
+        the index of this step
     train : str
         label of the associated train
     zone : int
@@ -41,11 +43,17 @@ def build_step(train: str, zone: int, prev_idx: int, min_t_in: int,
         true if the arrival time must match the min_t_in
     ponderation : float
         The step ponderation in the objective function
+    overlap : int
+        The overlap duration for this step
+    next : int
+        index of the next step
     """
     return {
+        "idx": idx,
         "train": train,
         "zone": zone,
         "prev": prev_idx,
+        "next": next,
         "min_t_in": min_t_in,
         "min_t_out": min_t_out,
         "min_duration": min_duration,
@@ -88,6 +96,7 @@ def steps_from_schedule(
 
     steps = []
 
+    global_idx = 0
     for train_idx, train in enumerate(trains):
         prev_step = -1
         prev_zone = None
@@ -108,7 +117,11 @@ def steps_from_schedule(
                 else weights.loc[zone][train]
             )
 
+            if prev_step >= 0:
+                steps[prev_step]["next"] = global_idx
+
             steps.append(build_step(
+                idx=global_idx,
                 train=train,
                 zone=zones.index(zone),
                 prev_idx=prev_step,
@@ -121,6 +134,7 @@ def steps_from_schedule(
             ))
             prev_zone = zone
             prev_step = len(steps) - 1
+            global_idx = global_idx + 1
 
     return steps
 

src/cpagent/variables.py

@@ -37,3 +37,59 @@ def _create_variables(
             self.t_out[i],
             f"t_out[{i}]")
         for i, step in enumerate(self.steps)]
+
+    # Precedence variables
+
+    # first_si : 1 if the step si is the first step to pass in its zone
+    # last_si : 1 if the step si is the last step to pass in its zone
+    # prec_si_sj : 1 if the step si is just before the step sj
+    #               for this to be true the two steps needs to be the same zone
+    #               and different trains
+    # diff_itinerary_si_sj : 1 if the step si and sj have a different zone
+    #                        AFTER the one they share
+
+    self.firsts = [
+        model.NewIntVar(
+            0,
+            1,
+            f"first_s{i}")
+        for i, _ in enumerate(self.steps)]
+    self.lasts = [
+        model.NewIntVar(
+            0,
+            1,
+            f"last_s{i}")
+        for i, _ in enumerate(self.steps)]
+    self.precs = [
+        [
+            model.NewIntVar(
+                0,
+                1 if (
+                    step_i["train"] != step_j["train"]
+                    and step_i["zone"] == step_j["zone"]
+                )
+                else 0,
+                f"prec_s{i}_s{j}")
+            for i, step_i in enumerate(self.steps)
+        ]
+        for j, step_j in enumerate(self.steps)
+    ]
+
+    self.diff_itineraries = [
+        [
+            model.NewIntVar(
+                1 if (
+                    step_i["next"] >= 0
+                    and step_j["next"] >= 0
+                    and (
+                        self.steps[step_i["next"]]["zone"]
+                        != self.steps[step_j["next"]]["zone"]
+                    )
+                )
+                else 0,
+                1,
+                f"diff_itinierary_s{i}_s{j}")
+            for i, step_i in enumerate(self.steps)
+        ]
+        for j, step_j in enumerate(self.steps)
+    ]

tests/conftest.py

@@ -86,13 +86,13 @@ def use_case_cp_4_zones_switch():
     """
     steps = []
 
-    steps.append(build_step(0, 0, -1, 0, 10, 20, False, 2))
-    steps.append(build_step(0, 2, 0, 10, 20, 10, True, 2))
-    steps.append(build_step(0, 3, 1, 20, 30, 10, False, 2))
+    steps.append(build_step(0, 0, 0, -1, 0, 10, 20, False, 2))
+    steps.append(build_step(1, 0, 2, 0, 10, 20, 10, True, 2))
+    steps.append(build_step(2, 0, 3, 1, 20, 30, 10, False, 2))
 
-    steps.append(build_step(1, 1, -1, 20, 30, 10, False, 2))
-    steps.append(build_step(1, 2, 3, 30, 40, 10, True, 2))
-    steps.append(build_step(1, 3, 4, 40, 50, 10, False, 2))
+    steps.append(build_step(3, 1, 1, -1, 20, 30, 10, False, 2))
+    steps.append(build_step(4, 1, 2, 3, 30, 40, 10, True, 2))
+    steps.append(build_step(5, 1, 3, 4, 40, 50, 10, False, 2))
 
     yield 4, 2, steps
 
@@ -112,18 +112,18 @@ def use_case_delay_conv():
 
     steps = []
 
-    steps.append(build_step(0, 0, -1, 0, 10,
+    steps.append(build_step(0, 0, 0, -1, 0, 10,
                             10 + delay_at_first_departure, False, 1))
-    steps.append(build_step(0, 2, 0, 10, 20, 10, True, 1))
-    steps.append(build_step(0, 3, 1, 20, 30, 10, False, 1))
-    steps.append(build_step(0, 4, 2, 30, 40, 10, True, 1))
-    steps.append(build_step(0, 5, 3, 40, 50, 10, False, 1))
+    steps.append(build_step(1, 0, 2, 0, 10, 20, 10, True, 1))
+    steps.append(build_step(2, 0, 3, 1, 20, 30, 10, False, 1))
+    steps.append(build_step(3, 0, 4, 2, 30, 40, 10, True, 1))
+    steps.append(build_step(4, 0, 5, 3, 40, 50, 10, False, 1))
 
-    steps.append(build_step(1, 1, -1, 20, 30, 10, False, 1))
-    steps.append(build_step(1, 2, 5, 30, 40, 10, True, 1))
-    steps.append(build_step(1, 3, 6, 40, 50, 10, False, 1))
-    steps.append(build_step(1, 4, 6, 50, 60, 10, True, 1))
-    steps.append(build_step(1, 6, 8, 60, 70, 10, False, 1))
+    steps.append(build_step(5, 1, 1, -1, 20, 30, 10, False, 1))
+    steps.append(build_step(6, 1, 2, 5, 30, 40, 10, True, 1))
+    steps.append(build_step(7, 1, 3, 6, 40, 50, 10, False, 1))
+    steps.append(build_step(8, 1, 4, 6, 50, 60, 10, True, 1))
+    steps.append(build_step(9, 1, 6, 8, 60, 70, 10, False, 1))
 
     yield 7, 2, steps
 
@@ -141,11 +141,11 @@ def use_case_infeasible():
     """
     steps = []
 
-    steps.append(build_step(0, 0, -1, 0, 10, 30, False))
-    steps.append(build_step(0, 1, 0, 10, 20, 10, False))
+    steps.append(build_step(0, 0, 0, -1, 0, 10, 30, False))
+    steps.append(build_step(1, 0, 1, 0, 10, 20, 10, False))
 
-    steps.append(build_step(1, 0, -1, 10, 20, 10, False))
-    steps.append(build_step(1, 1, 0, 20, 30, 10, False))
+    steps.append(build_step(2, 1, 0, -1, 10, 20, 10, False))
+    steps.append(build_step(3, 1, 1, 0, 20, 30, 10, False))
 
     yield 2, 2, steps
 
@@ -162,11 +162,11 @@ def use_case_straight_line_2t():
         steps
     """
     steps = []
-    steps.append(build_step(0, 0, -1, 0, 10, 10, False))
-    steps.append(build_step(0, 1, 0, 10, 20, 20, False))
+    steps.append(build_step(0, 0, 0, -1, 0, 10, 10, False))
+    steps.append(build_step(1, 0, 1, 0, 10, 20, 20, False))
 
-    steps.append(build_step(1, 0, -1, 10, 20, 10, False))
-    steps.append(build_step(1, 1, 2, 20, 30, 10, False))
+    steps.append(build_step(2, 1, 0, -1, 10, 20, 10, False))
+    steps.append(build_step(3, 1, 1, 2, 20, 30, 10, False))
 
     yield 2, 2, steps
 
@@ -183,10 +183,10 @@ def use_case_empty_zone():
     """
     steps = []
 
-    steps.append(build_step(0, 0, -1, 0, 10, 10, False))
-    steps.append(build_step(0, 1, 0, 10, 20, 20, False))
+    steps.append(build_step(0, 0, 0, -1, 0, 10, 10, False))
+    steps.append(build_step(1, 0, 1, 0, 10, 20, 20, False))
 
-    steps.append(build_step(1, 0, -1, 10, 20, 10, False))
-    steps.append(build_step(1, 1, 2, 20, 30, 10, False))
+    steps.append(build_step(2, 1, 0, -1, 10, 20, 10, False))
+    steps.append(build_step(3, 1, 1, 2, 20, 30, 10, False))
 
     yield 2, 3, steps

tests/cpagent/test_osrd_adapter.py

@@ -46,12 +46,12 @@ def test_regulation_problem_from_osrd(osrd_point_switch):
 
     oracle_steps = []
 
-    oracle_steps.append(build_step("train0", 0, -1, 400, 642, 242, False, 1, 0)) # noqa
-    oracle_steps.append(build_step("train0", 1, 0, 636, 647, 10, True, 1, 5))
-    oracle_steps.append(build_step("train0", 2, 1, 642, 849, 207, True, 1, 5))
-    oracle_steps.append(build_step("train1", 3, -1, 100, 342, 242, False, 1, 0)) # noqa
-    oracle_steps.append(build_step("train1", 1, 3, 336, 347, 10, True, 1, 5))
-    oracle_steps.append(build_step("train1", 0, 4, 342, 549, 207, True, 1, 5))
+    oracle_steps.append(build_step(0, "train0", 0, -1, 400, 642, 242, False, 1, 0, 1)) # noqa
+    oracle_steps.append(build_step(1, "train0", 1, 0, 636, 647, 10, True, 1, 5, 2)) # noqa
+    oracle_steps.append(build_step(2, "train0", 2, 1, 642, 849, 207, True, 1, 5)) # noqa
+    oracle_steps.append(build_step(3, "train1", 3, -1, 100, 342, 242, False, 1, 0, 4)) # noqa
+    oracle_steps.append(build_step(4, "train1", 1, 3, 336, 347, 10, True, 1, 5, 5)) # noqa
+    oracle_steps.append(build_step(5, "train1", 0, 4, 342, 549, 207, True, 1, 5)) # noqa
 
     assert oracle_steps == steps