quast_decisiontree.problems.classes.tsp

quast_decisiontree.problems.classes.tsp

TSP problem framework.

logger module-attribute

logger = logging.getLogger('dt_logger')

TSPBadState

Bases: Exception

raised if a TSP instance is found in an invalid state

Source code in src/quast_decisiontree/problems/classes/tsp.py

class TSPBadState(Exception):
    """raised if a TSP instance is found in an invalid state"""

TSP

Bases: OptimizationProblem

class representing an instance of the TSP problem

Source code in src/quast_decisiontree/problems/classes/tsp.py

class TSP(OptimizationProblem):
    """class representing an instance of the TSP problem"""

    direct_encoding_modes = ("QUBO_condensed", "QUBO")
    distance_tol = 1e-5

    def __init__(self, adjacency_matrix: np.ndarray, positions: np.ndarray | None = None) -> None:
        """Constructs a TSP instance from an adjacency matrix, possibly also setting a positions
        array.

        Args:
            adjacency_matrix (np.ndarray): The adjacency matrix containing the distance between
                the cities of the TSP.
            positions (Optional[np.ndarray], optional): An optional list of positions with two
                coordinates for each city. It is not necessary since the TSP is fully defined
                via the adjacency matrix. Positions will be used for plotting. Defaults to None.
        """
        self.optimal_path_length = None
        self.adjacency_matrix = adjacency_matrix
        self.positions = positions
        self.graph = from_numpy_array(self.adjacency_matrix)
        self._check_state()

    @staticmethod
    def tsp_value(city_list: Sequence[int], adjacency_matrix: np.ndarray) -> float:
        """Calculate the total path length for a given city ordering.

        Args:
            city_list: Sequence of city indices representing the tour order.
            adjacency_matrix: Distance matrix between cities.

        Returns:
            Total round-trip path length.
        """
        n = len(city_list)
        total = 0.0
        for i in range(n):
            total += adjacency_matrix[city_list[i]][city_list[(i + 1) % n]]
        return total

    @property
    def adjacency_matrix(self):
        return self._adjacency_matrix

    @adjacency_matrix.setter
    def adjacency_matrix(self, matrix: ArrayLike):
        if matrix is None:
            self._adjacency_matrix = None
        else:
            self._adjacency_matrix = np.array(matrix)

    @property
    def num_nodes(self):
        return np.shape(self.adjacency_matrix)[0]

    @property
    def positions(self):
        return self._positions

    @positions.setter
    def positions(self, poslist: ArrayLike):
        if poslist is None:
            self._positions = None
        else:
            self._positions = np.array(poslist)

    def _check_state(self) -> None:
        """Checks the state of the TSP instance and raises errors if there are problems."""
        shape_adjacency = np.shape(self.adjacency_matrix)
        if len(shape_adjacency) != 2 or shape_adjacency[0] != shape_adjacency[1]:
            raise TSPBadState(f"Invalid shape {shape_adjacency!r} of adjacency matrix.")
        if self.positions is not None:
            shape_positions = np.shape(self.positions)
            if len(shape_positions) != 2 or shape_positions[0] != self.num_nodes:
                raise TSPBadState(f"Invalid shape {shape_positions!r} of position list.")

    def __eq__(self, other: object) -> bool:
        """checks equality of the underlying graphs (understood as isomorphism respecting
        weights)
        """
        if not isinstance(other, TSP):
            return NotImplemented
        return is_isomorphic(
            self.graph,
            other.graph,
            edge_match=lambda x, y: _edge_match(x, y, tol=self.distance_tol),
        )

    def _as_tour(self, result: Sequence | str) -> list:
        """Normalize a solution into a canonical 0-based tour visiting each city once.

        Accepts a one-hot bitstring or an integer sequence in full (``num_nodes``) or
        condensed (``num_nodes - 1``) form, using either 0-based or 1-based city labels.
        ``one_hot_to_integer`` yields 1-based labels; normalization keys off the minimum
        label so both conventions collapse to the same canonical tour.
        """
        if isinstance(result, str):
            result = one_hot_to_integer(result)
        result = list(result)
        num_nodes = self.num_nodes

        if len(result) == num_nodes:
            offset = min(result)
            return [city - offset for city in result]
        if len(result) == num_nodes - 1:
            offset = min(result)
            return [0] + [city - offset + 1 for city in result]

        raise ValueError(f"Result {result!r} does not match a tour of {num_nodes} cities.")

    @classmethod
    def from_coordinate_list(cls, coordinate_list: np.ndarray) -> "TSP":
        """Constructs a TSP instance from a coordinate list by inferring the adjacency matrix
        from the Euclidean distance of the coordinates.

        Args:
            coordinate_list (np.ndarray): A list of points with two coordinates each.

        Returns:
            TSP: The TSP instance defined by the coordinate list.
        """
        coordinates = np.asarray(coordinate_list, dtype=float)
        difference = coordinates[:, None, :] - coordinates[None, :, :]
        adjacency_matrix = np.linalg.norm(difference, axis=-1)

        ins = cls(adjacency_matrix)
        ins.positions = coordinates

        return ins

    @classmethod
    def create_random_instance(cls, size: int, seed: Any = None):
        """create a random TSP instance of the specified size.

        Parameters:
        size       how many nodes (cities) to create
        seed            an optional seed for the random number generator

        Returns:
        the newly created TSP instance
        """
        rng = rd.default_rng(seed)
        coordinate_list = rng.uniform(0, 100, (size, 2))

        ins = cls.from_coordinate_list(coordinate_list)
        return ins

    def display(self) -> go.Figure:
        """draws the problem graph of the TSP instance"""

        if self.positions is None:
            raise AttributeError(
                "The display function cannot be called when the attribute positions is not"
                " specified"
            )

        edge_x = []
        edge_y = []

        for edge in self.graph.edges:
            x1, y1 = self.positions[edge[0]]
            x2, y2 = self.positions[edge[1]]
            edge_x.append(x1)
            edge_x.append(x2)
            edge_x.append(None)
            edge_y.append(y1)
            edge_y.append(y2)
            edge_y.append(None)

        edge_trace = go.Scatter(
            x=edge_x, y=edge_y, line=dict(width=1, color="#000"), hoverinfo="skip", mode="lines"
        )

        node_x = []
        node_y = []
        node_text = []

        for i, pos in enumerate(self.positions):
            x, y = pos
            node_x.append(x)
            node_y.append(y)
            node_text.append(str(i))

        node_trace = go.Scatter(
            x=node_x,
            y=node_y,
            mode="markers+text",
            hoverinfo="x+y",
            marker=dict(color="LightSkyBlue", size=30, line_width=2),
            text=node_text,
        )

        fig = go.Figure(
            data=[edge_trace, node_trace],
            layout=go.Layout(
                plot_bgcolor="#FFF",
                showlegend=False,
                hovermode="closest",
                margin=dict(b=20, l=5, r=5, t=40),
                xaxis=dict(showgrid=False, zeroline=False, showticklabels=False),
                yaxis=dict(showgrid=False, zeroline=False, showticklabels=False),
            ),
        )

        return fig

    def display_solution(self, result: Sequence) -> None:
        """Draws the TSP instance with a given solution.

        Args:
            result (Sequence): The sequence of cities to draw alongside the TSP instance.
        """
        edges = [(result[i], result[(i + 1) % len(result)]) for i in range(len(result))]

        if self.positions is not None:
            draw_networkx(
                self.graph,
                pos={key: value for key, value in enumerate(self.positions)},
                edgelist=edges,
            )
        else:
            draw_networkx(self.graph, edgelist=edges)

    def evaluate_objective(self, result: Sequence | str) -> float:
        """returns the path length of a proposed path.

        Accepts a one-hot bitstring or an integer city sequence in condensed or full form
        (see :meth:`_as_tour` for the accepted conventions).
        """
        tour = self._as_tour(result)
        return self.tsp_value(tour, self.adjacency_matrix)

    def _brute_solve(self) -> None:
        """Tries to find an optimal path by brute force solving."""
        num_nodes = self.num_nodes
        if num_nodes > 8:
            logger.warning(
                "TSP has %d nodes and brute force might be slow. "
                "Consider submitting an optimal path length.",
                num_nodes,
            )
        optimal_path_length = math.inf
        for perm in permutations(range(1, num_nodes), num_nodes - 1):
            optimal_path_length = min(self.evaluate_objective(list(perm)), optimal_path_length)
        self.optimal_path_length = optimal_path_length
        logger.info("Optimal path length computed to %r.", self.optimal_path_length)

    def tsp_length(self, eigenstate: Mapping, num_nodes: int, encoding: str = "one-hot") -> float:
        """Calculates the expectation value of the path length for a given state. Infeasible
        basis states are ignored, for feasible ones the TSP length is calculated and added
        to a weighted average.

        Args:
            eigenstate (Mapping): A dictionary of key : val pairs where the keys correspond to
                the bitstrings (computational basis states) and the vals may be amplitudes,
                probabilities or shot counts associated with the bitstrings.
            num_nodes (int): The number of nodes in the associated TSP.
            encoding (str, optional): The encoding the bitstrings contained in the eigenstate
                come from. Defaults to "one-hot".

        Raises:
            ValueError: If an unknown encoding is given.

        Returns:
            float: The expectation value of the path length in the projection of the eigenstate
                on the feasible subspace.
        """

        numerator = 0
        denominator = 0

        if encoding == "one-hot":
            converter = one_hot_to_integer
        elif encoding == "binary":
            converter = binary_to_integer
        elif encoding == "edge":
            converter = edge_to_integer
        else:
            raise ValueError("Unknown encoding.")

        if math.isclose(sum(np.linalg.norm(x) ** 2 for x in eigenstate.values()), 1):
            for res_string, amplitude in eigenstate.items():
                if self.is_feasible(res_string, num_nodes, encoding):
                    numerator += np.linalg.norm(amplitude) ** 2 * self.evaluate_objective(
                        converter(res_string, num_nodes)
                    )
                    denominator += np.linalg.norm(amplitude) ** 2
            if numerator == 0 and denominator == 0:
                return math.inf
            else:
                return numerator / denominator
        else:
            for res_string, occurrence in eigenstate.items():
                if occurrence < 0:
                    raise ValueError(
                        "Negative probability or shot count encountered. "
                        "If results should be amplitudes, they aren't normalized"
                    )
                if self.is_feasible(res_string, num_nodes, encoding):
                    numerator += occurrence * self.evaluate_objective(
                        converter(res_string, num_nodes)
                    )
                    denominator += occurrence
            if numerator == 0 and denominator == 0:
                return math.inf
            else:
                return numerator / denominator

    @classmethod
    def convert_between_frames(cls, path: Sequence) -> list:
        """Converts a path from the time-frame to the city-frame or vice versa.

        Note: this only works if the given path is valid, i.e. all integers i in
        range(len(path)) occur only once.

        Parameter:
        path            list of integers

        Returns:
        converted_path  list of integers with switched interpretation
        """
        path = np.asarray(path, dtype=int)
        path = path - path.min()

        integer_sequence = np.zeros(len(path), dtype=int)
        for index, integer in enumerate(path):
            integer_sequence[integer] = index + 1

        return integer_sequence.tolist()

    def formulate_qubo(
        self, penalty_factor: float = 100, scaling_factor: float = 1
    ) -> tuple[float, np.ndarray]:
        """returns offset and qubo tensor for the TSP instance

        Parameters:
        penalty_factor      The factor penalizing the TSP constraints
        scaling_factor      A generic scaling factor for the cost function

        Returns:
        offset              The constant contribution to the objective function
        qubo_tensor         A qubo tensor encapsulating the interaction between the binary
                            variables
        """
        adjacency_matrix = np.array(self.adjacency_matrix, dtype=float)
        adjacency_matrix = (adjacency_matrix + adjacency_matrix.T) / 2
        n = adjacency_matrix.shape[0]

        i = np.arange(n)[:, None, None, None]
        j = np.arange(n)[None, :, None, None]
        k = np.arange(n)[None, None, :, None]
        el = np.arange(n)[None, None, None, :]

        position_diff = (j - el) % n
        position_adjacent = (position_diff == 1) | (position_diff == n - 1)
        upper_triangular = (n * i + j) < (n * k + el)
        mask = position_adjacent & upper_triangular

        qubo_tensor = np.where(mask, adjacency_matrix[i, k], 0.0)

        diag = np.arange(n)
        rows, cols = np.meshgrid(diag, diag, indexing="ij")
        qubo_tensor[rows, cols, rows, cols] -= 2 * penalty_factor

        j_idx, k_idx = np.triu_indices(n, k=1)
        for city in range(n):
            qubo_tensor[j_idx, city, k_idx, city] += 2 * penalty_factor
            qubo_tensor[city, j_idx, city, k_idx] += 2 * penalty_factor

        offset = 2 * n * penalty_factor

        return scaling_factor * offset, scaling_factor * qubo_tensor

    def formulate_qubo_condensed(
        self, penalty_factor: float = 100, scaling_factor: float = 1
    ) -> tuple[float, np.ndarray]:
        """returns offset and qubo tensor while removing the cyclic permutation
        freedom by fixing the first city to be visited first.

        Args and return values as for formulate_qubo method.
        """
        offset, wasteful_qubo = self.formulate_qubo(penalty_factor, scaling_factor)
        condensed_qubo = wasteful_qubo[1:, 1:, 1:, 1:].copy()
        num_nodes = np.shape(wasteful_qubo)[0]
        for city_index in range(1, num_nodes):
            condensed_qubo[city_index - 1, 0, city_index - 1, 0] += wasteful_qubo[
                0, 0, city_index, 1
            ]
            condensed_qubo[city_index - 1, -1, city_index - 1, -1] += wasteful_qubo[
                0, 0, city_index, -1
            ]
        offset = 2 * (num_nodes - 1) * penalty_factor * scaling_factor
        return offset, condensed_qubo

    def decode_result(
        self,
        sol_bitstring: str,
        mode: str = "QUBO_condensed",
    ) -> tuple[list, float]:
        """the inverse function for formulate_problem(). Takes a solution bitstring of the given
        formulation and converts it to an integer solution vector and the path length.
        """
        self._check_mode_support(mode)
        if mode.lower() not in ("qubo", "qubo_condensed"):
            raise ValueError(f"Unsupported mode {mode!r}.")

        tour = self._as_tour(sol_bitstring)
        return tour, self.evaluate_objective(tour)

    def formulate_problem(
        self, mode: str = "QUBO_condensed", penalty_factor: float = 100, scaling_factor: float = 1
    ) -> tuple[float, np.ndarray]:
        """returns an offset and QUBO tensor for the TSP instance"""

        self._check_mode_support(mode)
        mode = mode.lower()

        if mode == "qubo":
            return self.formulate_qubo(penalty_factor, scaling_factor)
        if mode == "qubo_condensed":
            return self.formulate_qubo_condensed(penalty_factor, scaling_factor)

        raise ValueError(f"Unsupported mode {mode!r}.")

    def to_dict(self) -> dict:
        """Converts the TSP instance to a dictionary.

        Returns:
            dict: A dictionary with problem_class, distance_matrix, and optionally
                coordinate_list.
        """
        out_dict = dict(problem_class="TSP")
        out_dict["distance_matrix"] = self.adjacency_matrix.tolist()
        if self.positions is not None:
            out_dict["coordinate_list"] = self.positions.tolist()

        return out_dict

    @classmethod
    def from_dict(cls, problem_dict: Mapping) -> "TSP":
        """constructs a TSP instance from a problem dictionary

        Either the key "distance_matrix" or "coordinate_list" must be present.
        """
        distance_matrix = problem_dict.get("distance_matrix")
        coordinate_list = problem_dict.get("coordinate_list")

        if distance_matrix is None and coordinate_list is None:
            raise InvalidProblemDictError(
                "Problem dictionary must contain 'distance_matrix' or 'coordinate_list'."
            )

        if distance_matrix is not None:
            distance_matrix = np.array(distance_matrix)
            if distance_matrix.ndim != 2 or distance_matrix.shape[0] != distance_matrix.shape[1]:
                raise InvalidProblemDictError(
                    f"Invalid distance matrix shape {distance_matrix.shape!r}."
                )
            try:
                ins = cls(distance_matrix)
            except TSPBadState as e:
                raise InvalidProblemDictError(str(e)) from e
            if coordinate_list is not None:
                ins.positions = coordinate_list
            return ins

        try:
            return cls.from_coordinate_list(np.array(coordinate_list))
        except TSPBadState as e:
            raise InvalidProblemDictError(str(e)) from e

    @classmethod
    def is_feasible(
        cls, solution_string: str, num_nodes: int | None = None, encoding: str = "one-hot"
    ) -> bool:
        """Determines whether a solution bitstring represents a feasible solution.

        For convenience, returns False if solution_string = None.
        """
        if solution_string is None:
            return False

        if encoding == "one-hot":
            try:
                integer_sequence = one_hot_to_integer(solution_string)
            except (ValueError, IndexError, InvalidBitstring):
                return False
        elif encoding == "binary":
            try:
                integer_sequence = binary_to_integer(solution_string, num_nodes)
            except InvalidBitstring:
                return False
        elif encoding == "edge":
            try:
                integer_sequence = edge_to_integer(solution_string)
            except (ValueError, IndexError, InvalidBitstring):
                return False
        else:
            raise ValueError("Unknown encoding.")

        integer_sequence = list(integer_sequence)
        if not integer_sequence:
            return False

        offset = min(integer_sequence)
        normalized = [i - offset for i in integer_sequence]
        return set(normalized) == set(range(len(normalized)))

    @classmethod
    def internal_loops(cls, x: Sequence) -> list:
        """Returns the loops contained in a (candidate) TSP solution given in edge encoding"""

        num_nodes = len(np.nonzero(x)[0])
        all_edges = np.array(list(combinations(range(num_nodes), 2)))
        occ_edges = all_edges[np.nonzero(x)[0], :]

        loops = []
        k = 1
        row_idx = 0
        integer_list = np.zeros(num_nodes, dtype=int)

        while len(occ_edges) > 0:
            row = occ_edges[row_idx]
            integer_list[k] = row[row != integer_list[k - 1]]
            occ_edges = np.delete(occ_edges, row_idx, axis=0)
            if len(occ_edges) == 0:
                loops.append(integer_list[:k])
                break
            row_idx = np.where(occ_edges == integer_list[k])[0]
            if len(row_idx) == 0:
                loops.append(integer_list[:k])
                integer_list = np.zeros(num_nodes, dtype=int)
                integer_list[0] = occ_edges[0, 0]
                k = 1
                row_idx = 0
                continue
            k += 1

        return loops

    def get_ising_offset(self, penalty_factor: float = 100, scaling_factor: float = 1) -> float:
        """calculates the offset between QUBO and Ising formulation of the problem.

        This is the analytic expression for the normal qubo formulation of TSP.
        """
        total_adjacency = np.sum(self.adjacency_matrix)
        num_nodes = len(self.adjacency_matrix)
        return (
            scaling_factor
            * (
                num_nodes * total_adjacency
                + 2 * penalty_factor * (num_nodes**3 - 3 * num_nodes**2)
            )
            / 4
        )

direct_encoding_modes class-attribute instance-attribute

direct_encoding_modes = ('QUBO_condensed', 'QUBO')

distance_tol class-attribute instance-attribute

distance_tol = 1e-05

optimal_path_length instance-attribute

optimal_path_length = None

graph instance-attribute

graph = from_numpy_array(self.adjacency_matrix)

adjacency_matrix property writable

adjacency_matrix

num_nodes property

num_nodes

positions property writable

positions

__init__

__init__(adjacency_matrix, positions=None)

Constructs a TSP instance from an adjacency matrix, possibly also setting a positions array.

Parameters:

Name	Type	Description	Default
adjacency_matrix	ndarray	The adjacency matrix containing the distance between the cities of the TSP.	required
positions	Optional[ndarray]	An optional list of positions with two coordinates for each city. It is not necessary since the TSP is fully defined via the adjacency matrix. Positions will be used for plotting. Defaults to None.	None

Source code in src/quast_decisiontree/problems/classes/tsp.py

def __init__(self, adjacency_matrix: np.ndarray, positions: np.ndarray | None = None) -> None:
    """Constructs a TSP instance from an adjacency matrix, possibly also setting a positions
    array.

    Args:
        adjacency_matrix (np.ndarray): The adjacency matrix containing the distance between
            the cities of the TSP.
        positions (Optional[np.ndarray], optional): An optional list of positions with two
            coordinates for each city. It is not necessary since the TSP is fully defined
            via the adjacency matrix. Positions will be used for plotting. Defaults to None.
    """
    self.optimal_path_length = None
    self.adjacency_matrix = adjacency_matrix
    self.positions = positions
    self.graph = from_numpy_array(self.adjacency_matrix)
    self._check_state()

tsp_value staticmethod

tsp_value(city_list, adjacency_matrix)

Calculate the total path length for a given city ordering.

Parameters:

Name	Type	Description	Default
city_list	Sequence[int]	Sequence of city indices representing the tour order.	required
adjacency_matrix	ndarray	Distance matrix between cities.	required

Returns:

Type	Description
float	Total round-trip path length.

Source code in src/quast_decisiontree/problems/classes/tsp.py

@staticmethod
def tsp_value(city_list: Sequence[int], adjacency_matrix: np.ndarray) -> float:
    """Calculate the total path length for a given city ordering.

    Args:
        city_list: Sequence of city indices representing the tour order.
        adjacency_matrix: Distance matrix between cities.

    Returns:
        Total round-trip path length.
    """
    n = len(city_list)
    total = 0.0
    for i in range(n):
        total += adjacency_matrix[city_list[i]][city_list[(i + 1) % n]]
    return total

__eq__

__eq__(other)

checks equality of the underlying graphs (understood as isomorphism respecting weights)

Source code in src/quast_decisiontree/problems/classes/tsp.py

def __eq__(self, other: object) -> bool:
    """checks equality of the underlying graphs (understood as isomorphism respecting
    weights)
    """
    if not isinstance(other, TSP):
        return NotImplemented
    return is_isomorphic(
        self.graph,
        other.graph,
        edge_match=lambda x, y: _edge_match(x, y, tol=self.distance_tol),
    )

from_coordinate_list classmethod

from_coordinate_list(coordinate_list)

Constructs a TSP instance from a coordinate list by inferring the adjacency matrix from the Euclidean distance of the coordinates.

Parameters:

Name	Type	Description	Default
coordinate_list	ndarray	A list of points with two coordinates each.	required

Returns:

Name	Type	Description
TSP	TSP	The TSP instance defined by the coordinate list.

Source code in src/quast_decisiontree/problems/classes/tsp.py

@classmethod
def from_coordinate_list(cls, coordinate_list: np.ndarray) -> "TSP":
    """Constructs a TSP instance from a coordinate list by inferring the adjacency matrix
    from the Euclidean distance of the coordinates.

    Args:
        coordinate_list (np.ndarray): A list of points with two coordinates each.

    Returns:
        TSP: The TSP instance defined by the coordinate list.
    """
    coordinates = np.asarray(coordinate_list, dtype=float)
    difference = coordinates[:, None, :] - coordinates[None, :, :]
    adjacency_matrix = np.linalg.norm(difference, axis=-1)

    ins = cls(adjacency_matrix)
    ins.positions = coordinates

    return ins

create_random_instance classmethod

create_random_instance(size, seed=None)

create a random TSP instance of the specified size.

Parameters: size how many nodes (cities) to create seed an optional seed for the random number generator

Returns: the newly created TSP instance

Source code in src/quast_decisiontree/problems/classes/tsp.py

@classmethod
def create_random_instance(cls, size: int, seed: Any = None):
    """create a random TSP instance of the specified size.

    Parameters:
    size       how many nodes (cities) to create
    seed            an optional seed for the random number generator

    Returns:
    the newly created TSP instance
    """
    rng = rd.default_rng(seed)
    coordinate_list = rng.uniform(0, 100, (size, 2))

    ins = cls.from_coordinate_list(coordinate_list)
    return ins

display

display()

draws the problem graph of the TSP instance

Source code in src/quast_decisiontree/problems/classes/tsp.py

def display(self) -> go.Figure:
    """draws the problem graph of the TSP instance"""

    if self.positions is None:
        raise AttributeError(
            "The display function cannot be called when the attribute positions is not"
            " specified"
        )

    edge_x = []
    edge_y = []

    for edge in self.graph.edges:
        x1, y1 = self.positions[edge[0]]
        x2, y2 = self.positions[edge[1]]
        edge_x.append(x1)
        edge_x.append(x2)
        edge_x.append(None)
        edge_y.append(y1)
        edge_y.append(y2)
        edge_y.append(None)

    edge_trace = go.Scatter(
        x=edge_x, y=edge_y, line=dict(width=1, color="#000"), hoverinfo="skip", mode="lines"
    )

    node_x = []
    node_y = []
    node_text = []

    for i, pos in enumerate(self.positions):
        x, y = pos
        node_x.append(x)
        node_y.append(y)
        node_text.append(str(i))

    node_trace = go.Scatter(
        x=node_x,
        y=node_y,
        mode="markers+text",
        hoverinfo="x+y",
        marker=dict(color="LightSkyBlue", size=30, line_width=2),
        text=node_text,
    )

    fig = go.Figure(
        data=[edge_trace, node_trace],
        layout=go.Layout(
            plot_bgcolor="#FFF",
            showlegend=False,
            hovermode="closest",
            margin=dict(b=20, l=5, r=5, t=40),
            xaxis=dict(showgrid=False, zeroline=False, showticklabels=False),
            yaxis=dict(showgrid=False, zeroline=False, showticklabels=False),
        ),
    )

    return fig

display_solution

display_solution(result)

Draws the TSP instance with a given solution.

Parameters:

Name	Type	Description	Default
result	Sequence	The sequence of cities to draw alongside the TSP instance.	required

Source code in src/quast_decisiontree/problems/classes/tsp.py

def display_solution(self, result: Sequence) -> None:
    """Draws the TSP instance with a given solution.

    Args:
        result (Sequence): The sequence of cities to draw alongside the TSP instance.
    """
    edges = [(result[i], result[(i + 1) % len(result)]) for i in range(len(result))]

    if self.positions is not None:
        draw_networkx(
            self.graph,
            pos={key: value for key, value in enumerate(self.positions)},
            edgelist=edges,
        )
    else:
        draw_networkx(self.graph, edgelist=edges)

evaluate_objective

evaluate_objective(result)

returns the path length of a proposed path.

Accepts a one-hot bitstring or an integer city sequence in condensed or full form (see :meth:_as_tour for the accepted conventions).

Source code in src/quast_decisiontree/problems/classes/tsp.py

def evaluate_objective(self, result: Sequence | str) -> float:
    """returns the path length of a proposed path.

    Accepts a one-hot bitstring or an integer city sequence in condensed or full form
    (see :meth:`_as_tour` for the accepted conventions).
    """
    tour = self._as_tour(result)
    return self.tsp_value(tour, self.adjacency_matrix)

tsp_length

tsp_length(eigenstate, num_nodes, encoding='one-hot')

Calculates the expectation value of the path length for a given state. Infeasible basis states are ignored, for feasible ones the TSP length is calculated and added to a weighted average.

Parameters:

Name	Type	Description	Default
eigenstate	Mapping	A dictionary of key : val pairs where the keys correspond to the bitstrings (computational basis states) and the vals may be amplitudes, probabilities or shot counts associated with the bitstrings.	required
num_nodes	int	The number of nodes in the associated TSP.	required
encoding	str	The encoding the bitstrings contained in the eigenstate come from. Defaults to "one-hot".	'one-hot'

Raises:

Type	Description
ValueError	If an unknown encoding is given.

Returns:

Name	Type	Description
float	float	The expectation value of the path length in the projection of the eigenstate on the feasible subspace.

Source code in src/quast_decisiontree/problems/classes/tsp.py

def tsp_length(self, eigenstate: Mapping, num_nodes: int, encoding: str = "one-hot") -> float:
    """Calculates the expectation value of the path length for a given state. Infeasible
    basis states are ignored, for feasible ones the TSP length is calculated and added
    to a weighted average.

    Args:
        eigenstate (Mapping): A dictionary of key : val pairs where the keys correspond to
            the bitstrings (computational basis states) and the vals may be amplitudes,
            probabilities or shot counts associated with the bitstrings.
        num_nodes (int): The number of nodes in the associated TSP.
        encoding (str, optional): The encoding the bitstrings contained in the eigenstate
            come from. Defaults to "one-hot".

    Raises:
        ValueError: If an unknown encoding is given.

    Returns:
        float: The expectation value of the path length in the projection of the eigenstate
            on the feasible subspace.
    """

    numerator = 0
    denominator = 0

    if encoding == "one-hot":
        converter = one_hot_to_integer
    elif encoding == "binary":
        converter = binary_to_integer
    elif encoding == "edge":
        converter = edge_to_integer
    else:
        raise ValueError("Unknown encoding.")

    if math.isclose(sum(np.linalg.norm(x) ** 2 for x in eigenstate.values()), 1):
        for res_string, amplitude in eigenstate.items():
            if self.is_feasible(res_string, num_nodes, encoding):
                numerator += np.linalg.norm(amplitude) ** 2 * self.evaluate_objective(
                    converter(res_string, num_nodes)
                )
                denominator += np.linalg.norm(amplitude) ** 2
        if numerator == 0 and denominator == 0:
            return math.inf
        else:
            return numerator / denominator
    else:
        for res_string, occurrence in eigenstate.items():
            if occurrence < 0:
                raise ValueError(
                    "Negative probability or shot count encountered. "
                    "If results should be amplitudes, they aren't normalized"
                )
            if self.is_feasible(res_string, num_nodes, encoding):
                numerator += occurrence * self.evaluate_objective(
                    converter(res_string, num_nodes)
                )
                denominator += occurrence
        if numerator == 0 and denominator == 0:
            return math.inf
        else:
            return numerator / denominator

convert_between_frames classmethod

convert_between_frames(path)

Converts a path from the time-frame to the city-frame or vice versa.

Note: this only works if the given path is valid, i.e. all integers i in range(len(path)) occur only once.

Parameter: path list of integers

Returns: converted_path list of integers with switched interpretation

Source code in src/quast_decisiontree/problems/classes/tsp.py

@classmethod
def convert_between_frames(cls, path: Sequence) -> list:
    """Converts a path from the time-frame to the city-frame or vice versa.

    Note: this only works if the given path is valid, i.e. all integers i in
    range(len(path)) occur only once.

    Parameter:
    path            list of integers

    Returns:
    converted_path  list of integers with switched interpretation
    """
    path = np.asarray(path, dtype=int)
    path = path - path.min()

    integer_sequence = np.zeros(len(path), dtype=int)
    for index, integer in enumerate(path):
        integer_sequence[integer] = index + 1

    return integer_sequence.tolist()

formulate_qubo

formulate_qubo(penalty_factor=100, scaling_factor=1)

returns offset and qubo tensor for the TSP instance

Parameters: penalty_factor The factor penalizing the TSP constraints scaling_factor A generic scaling factor for the cost function

Returns: offset The constant contribution to the objective function qubo_tensor A qubo tensor encapsulating the interaction between the binary variables

Source code in src/quast_decisiontree/problems/classes/tsp.py

def formulate_qubo(
    self, penalty_factor: float = 100, scaling_factor: float = 1
) -> tuple[float, np.ndarray]:
    """returns offset and qubo tensor for the TSP instance

    Parameters:
    penalty_factor      The factor penalizing the TSP constraints
    scaling_factor      A generic scaling factor for the cost function

    Returns:
    offset              The constant contribution to the objective function
    qubo_tensor         A qubo tensor encapsulating the interaction between the binary
                        variables
    """
    adjacency_matrix = np.array(self.adjacency_matrix, dtype=float)
    adjacency_matrix = (adjacency_matrix + adjacency_matrix.T) / 2
    n = adjacency_matrix.shape[0]

    i = np.arange(n)[:, None, None, None]
    j = np.arange(n)[None, :, None, None]
    k = np.arange(n)[None, None, :, None]
    el = np.arange(n)[None, None, None, :]

    position_diff = (j - el) % n
    position_adjacent = (position_diff == 1) | (position_diff == n - 1)
    upper_triangular = (n * i + j) < (n * k + el)
    mask = position_adjacent & upper_triangular

    qubo_tensor = np.where(mask, adjacency_matrix[i, k], 0.0)

    diag = np.arange(n)
    rows, cols = np.meshgrid(diag, diag, indexing="ij")
    qubo_tensor[rows, cols, rows, cols] -= 2 * penalty_factor

    j_idx, k_idx = np.triu_indices(n, k=1)
    for city in range(n):
        qubo_tensor[j_idx, city, k_idx, city] += 2 * penalty_factor
        qubo_tensor[city, j_idx, city, k_idx] += 2 * penalty_factor

    offset = 2 * n * penalty_factor

    return scaling_factor * offset, scaling_factor * qubo_tensor

formulate_qubo_condensed

formulate_qubo_condensed(
    penalty_factor=100, scaling_factor=1
)

returns offset and qubo tensor while removing the cyclic permutation freedom by fixing the first city to be visited first.

Args and return values as for formulate_qubo method.

Source code in src/quast_decisiontree/problems/classes/tsp.py

def formulate_qubo_condensed(
    self, penalty_factor: float = 100, scaling_factor: float = 1
) -> tuple[float, np.ndarray]:
    """returns offset and qubo tensor while removing the cyclic permutation
    freedom by fixing the first city to be visited first.

    Args and return values as for formulate_qubo method.
    """
    offset, wasteful_qubo = self.formulate_qubo(penalty_factor, scaling_factor)
    condensed_qubo = wasteful_qubo[1:, 1:, 1:, 1:].copy()
    num_nodes = np.shape(wasteful_qubo)[0]
    for city_index in range(1, num_nodes):
        condensed_qubo[city_index - 1, 0, city_index - 1, 0] += wasteful_qubo[
            0, 0, city_index, 1
        ]
        condensed_qubo[city_index - 1, -1, city_index - 1, -1] += wasteful_qubo[
            0, 0, city_index, -1
        ]
    offset = 2 * (num_nodes - 1) * penalty_factor * scaling_factor
    return offset, condensed_qubo

decode_result

decode_result(sol_bitstring, mode='QUBO_condensed')

the inverse function for formulate_problem(). Takes a solution bitstring of the given formulation and converts it to an integer solution vector and the path length.

Source code in src/quast_decisiontree/problems/classes/tsp.py

def decode_result(
    self,
    sol_bitstring: str,
    mode: str = "QUBO_condensed",
) -> tuple[list, float]:
    """the inverse function for formulate_problem(). Takes a solution bitstring of the given
    formulation and converts it to an integer solution vector and the path length.
    """
    self._check_mode_support(mode)
    if mode.lower() not in ("qubo", "qubo_condensed"):
        raise ValueError(f"Unsupported mode {mode!r}.")

    tour = self._as_tour(sol_bitstring)
    return tour, self.evaluate_objective(tour)

formulate_problem

formulate_problem(
    mode="QUBO_condensed",
    penalty_factor=100,
    scaling_factor=1,
)

returns an offset and QUBO tensor for the TSP instance

Source code in src/quast_decisiontree/problems/classes/tsp.py

def formulate_problem(
    self, mode: str = "QUBO_condensed", penalty_factor: float = 100, scaling_factor: float = 1
) -> tuple[float, np.ndarray]:
    """returns an offset and QUBO tensor for the TSP instance"""

    self._check_mode_support(mode)
    mode = mode.lower()

    if mode == "qubo":
        return self.formulate_qubo(penalty_factor, scaling_factor)
    if mode == "qubo_condensed":
        return self.formulate_qubo_condensed(penalty_factor, scaling_factor)

    raise ValueError(f"Unsupported mode {mode!r}.")

to_dict

to_dict()

Converts the TSP instance to a dictionary.

Returns:

Name	Type	Description
dict	dict	A dictionary with problem_class, distance_matrix, and optionally coordinate_list.

Source code in src/quast_decisiontree/problems/classes/tsp.py

def to_dict(self) -> dict:
    """Converts the TSP instance to a dictionary.

    Returns:
        dict: A dictionary with problem_class, distance_matrix, and optionally
            coordinate_list.
    """
    out_dict = dict(problem_class="TSP")
    out_dict["distance_matrix"] = self.adjacency_matrix.tolist()
    if self.positions is not None:
        out_dict["coordinate_list"] = self.positions.tolist()

    return out_dict

from_dict classmethod

from_dict(problem_dict)

constructs a TSP instance from a problem dictionary

Either the key "distance_matrix" or "coordinate_list" must be present.

Source code in src/quast_decisiontree/problems/classes/tsp.py

@classmethod
def from_dict(cls, problem_dict: Mapping) -> "TSP":
    """constructs a TSP instance from a problem dictionary

    Either the key "distance_matrix" or "coordinate_list" must be present.
    """
    distance_matrix = problem_dict.get("distance_matrix")
    coordinate_list = problem_dict.get("coordinate_list")

    if distance_matrix is None and coordinate_list is None:
        raise InvalidProblemDictError(
            "Problem dictionary must contain 'distance_matrix' or 'coordinate_list'."
        )

    if distance_matrix is not None:
        distance_matrix = np.array(distance_matrix)
        if distance_matrix.ndim != 2 or distance_matrix.shape[0] != distance_matrix.shape[1]:
            raise InvalidProblemDictError(
                f"Invalid distance matrix shape {distance_matrix.shape!r}."
            )
        try:
            ins = cls(distance_matrix)
        except TSPBadState as e:
            raise InvalidProblemDictError(str(e)) from e
        if coordinate_list is not None:
            ins.positions = coordinate_list
        return ins

    try:
        return cls.from_coordinate_list(np.array(coordinate_list))
    except TSPBadState as e:
        raise InvalidProblemDictError(str(e)) from e

is_feasible classmethod

is_feasible(
    solution_string, num_nodes=None, encoding="one-hot"
)

Determines whether a solution bitstring represents a feasible solution.

For convenience, returns False if solution_string = None.

Source code in src/quast_decisiontree/problems/classes/tsp.py

@classmethod
def is_feasible(
    cls, solution_string: str, num_nodes: int | None = None, encoding: str = "one-hot"
) -> bool:
    """Determines whether a solution bitstring represents a feasible solution.

    For convenience, returns False if solution_string = None.
    """
    if solution_string is None:
        return False

    if encoding == "one-hot":
        try:
            integer_sequence = one_hot_to_integer(solution_string)
        except (ValueError, IndexError, InvalidBitstring):
            return False
    elif encoding == "binary":
        try:
            integer_sequence = binary_to_integer(solution_string, num_nodes)
        except InvalidBitstring:
            return False
    elif encoding == "edge":
        try:
            integer_sequence = edge_to_integer(solution_string)
        except (ValueError, IndexError, InvalidBitstring):
            return False
    else:
        raise ValueError("Unknown encoding.")

    integer_sequence = list(integer_sequence)
    if not integer_sequence:
        return False

    offset = min(integer_sequence)
    normalized = [i - offset for i in integer_sequence]
    return set(normalized) == set(range(len(normalized)))

internal_loops classmethod

internal_loops(x)

Returns the loops contained in a (candidate) TSP solution given in edge encoding

Source code in src/quast_decisiontree/problems/classes/tsp.py

@classmethod
def internal_loops(cls, x: Sequence) -> list:
    """Returns the loops contained in a (candidate) TSP solution given in edge encoding"""

    num_nodes = len(np.nonzero(x)[0])
    all_edges = np.array(list(combinations(range(num_nodes), 2)))
    occ_edges = all_edges[np.nonzero(x)[0], :]

    loops = []
    k = 1
    row_idx = 0
    integer_list = np.zeros(num_nodes, dtype=int)

    while len(occ_edges) > 0:
        row = occ_edges[row_idx]
        integer_list[k] = row[row != integer_list[k - 1]]
        occ_edges = np.delete(occ_edges, row_idx, axis=0)
        if len(occ_edges) == 0:
            loops.append(integer_list[:k])
            break
        row_idx = np.where(occ_edges == integer_list[k])[0]
        if len(row_idx) == 0:
            loops.append(integer_list[:k])
            integer_list = np.zeros(num_nodes, dtype=int)
            integer_list[0] = occ_edges[0, 0]
            k = 1
            row_idx = 0
            continue
        k += 1

    return loops

get_ising_offset

get_ising_offset(penalty_factor=100, scaling_factor=1)

calculates the offset between QUBO and Ising formulation of the problem.

This is the analytic expression for the normal qubo formulation of TSP.

Source code in src/quast_decisiontree/problems/classes/tsp.py

def get_ising_offset(self, penalty_factor: float = 100, scaling_factor: float = 1) -> float:
    """calculates the offset between QUBO and Ising formulation of the problem.

    This is the analytic expression for the normal qubo formulation of TSP.
    """
    total_adjacency = np.sum(self.adjacency_matrix)
    num_nodes = len(self.adjacency_matrix)
    return (
        scaling_factor
        * (
            num_nodes * total_adjacency
            + 2 * penalty_factor * (num_nodes**3 - 3 * num_nodes**2)
        )
        / 4
    )