Quantum annealing is a metaheuristic for combinatorial optimization in quantum mechanics research(PhysRevE.58.5355),(qa_science). This new approach solves unconstrained optimization problems formulated as a Hamiltonian by evolving the time-dependent Schrödinger equation from a quantum mechanical superposition of all possible states to its ground state in physical systems(AIP_Nishimori). A quantum annealing machine is a special-purpose quantum computer that solves combinatorial optimization problems. D-wave is the first commercial quantum annealing machine(dwave).
Combinatorial optimization contributes to a more effective and reasonable life by minimizing the total costs or maximizing the benefits under certain constraints. Although sophisticated solvers can efficiently solve small problems, many practical problems are computationally intractable. Such problems are formally characterized as NP-hard; that is, there are no polynomial-time algorithms to solve the problems. Many researchers in computer science and operations research believe that there are no polynomial-time algorithms for NP-hard problems and continue to develop efficient heuristic algorithms to solve larger problems(10.5555/574848).
Quantum annealing is a metaheuristic that does not restrict the target problems. Thus far, some combinatorial optimization problems, but not many, have been formulated as Ising or Quadratic Unconstrained Binary Optimization (QUBO) models and solved using annealing machines(10.3389/fphy.2014.00005), (osti_1756458), (8643733), (NSP). Theoretical and experimental studies have contributed to improving the performance of annealing processes (PhysRevApplied.14.014100), (PhysRevLett.123.120501), (PhysRevLett.126.070505), (Hauke_2020). Useful software tools have been developed to utilize quantum annealing machines (pyqubo). Digital implementations of the quantum annealing process have also emerged because this new optimization technique has the potential to outperform traditional meta-heuristic algorithms. A graphics processing unit (GPU) based quantum annealing simulator has shown high-performance optimization by utilizing spatial and temporal parallelism during the annealing process (9057502). Therefore, we have high hopes for this new innovative optimizer.
However, obstacles to quantum annealing usability include the hardness of the parameter tuning and formulation of the target problems as Ising or QUBO models. The constraints of the original problem and objective functions are represented as a penalty function, minimized as an unconstrained optimization problem. Search processes can reach infeasible solutions, and to avoid such situations, it is necessary to carefully tune the parameters in pre-performing annealing process. The difficulty of the Ising model formulation is also an essential obstacle to the usability of quantum annealing. The new annealer requires users to create quadratic penalty functions, Ising, or QUBO models, for target optimization problems.
This study focuses solely on the latter problem, which is the difficulty of formulating the Ising and QUBO models. To overcome this problem, we propose a method for generating Ising or QUBO models from Petri net models representing the target optimization problems. Our proposal is a Petri net modeling approach in which, once we model the optimization problem with Petri nets, we can systematically obtain the Ising or QUBO formulation.
A Petri net is a mathematical and graphical modeling language for various systems, such as computer networks, manufacturing systems, transportation networks, biological systems, agricultural production processes, and other applications(24143)
. We only need the domain knowledge of the target systems for modeling because the rules and components of Petri nets are simple and intuitive. Petri net modeling is quite efficient for representing integer linear programming problems because Petri nets have mathematical forms, denoted as linear expressions(refId0). The authors proposed an algorithm for generating mixed-integer linear programming (MILP) based on Petri net models for practical scheduling and optimal resource assignment problems(MorikazuNAKAMURA2019). Our algorithm realizes the automatic formulation of MILPs for given combinatorial optimization problems.
This study extends our MILP version to quantum annealing models, where quadratic functions of binary variables need to be formulated. As far as we know, this is the first study to apply Petri net modeling to the quantum annealing. In this paper, we introduce a new notation calledbinary quadratic nets, which denote Ising or QUBO models, respectively. Our method incrementally targets binary quadratic nets from problem-domain Petri net models that represent target optimization problems.
This paper is organized as follows: Section 2 summarizes the basic knowledge on quantum annealing and Petri nets. Section 3 introduces a new class of Petri nets, binary quadratic nets, for representing Ising or QUBO models. Section 4 proposes a method for formulating Ising or QUBO models from problem-domain Petri nets and presents some examples. Finally, in Section 5, we present some concluding remarks and areas of future tasks.
This section provides the basic definitions and notations of Ising models and Petri nets.
2.1 Ising and QUBO Models
Quantum annealing is a metaheuristic for solving combinatorial optimization problems, where we try to find a value or for each of the Ising variables such that the following Hamiltonian of the Ising models gives the lowest energy state:
where is the magnetic field coefficient at site , and is the interaction coefficient between and . Ising variables correspond to discrete variables in metaheuristics. The lowest energy state, called the ground state of , provides the optimal solution. Figure 2 depicts a two-dimensional array Ising model. A circle represents a spin and is connected to four neighbor spins. The arrows and represent the up spin (+1) and down spin (-1), respectively.
A Hamiltonian may include an objective function, where the lower energy on the terms of the function contributes to a smaller objective value in minimization problems. For maximization problems, we reverse the sign of the Ising model parameters for the objective function and then minimize it. Constraints in optimization problems should also be composed of sets of terms in the Hamiltonian, where the lowest energy on the terms leads to the satisfaction of the corresponding constraints. Therefore, the annealing process needs to obtain feasible solutions that satisfy all constraints with a good quality of the objective value by exploring the lowest energy state of the total Hamiltonian.
Note that we can easily replace the type of decision variables from to , and vice versa.
2.2 Petri Net Fundamentals
A Petri net is a directed bipartite graph with initial marking . The bipartite graph has two types of vertex sets: a set of places and a set of transitions . The arc function defines the connection from places to transitions, and vice versa. The returned values of indicate the number of removed tokens from the starting place when or the number of generated tokens to the ending place when upon the firing of the transition.
Instead of an arc function, a matrix representation is often used for the connection. Here, and show the incident matrix of size to indicate the connection from places to transitions and from transitions to places, respectively. It should be noted that and .
Tokens are located in places such that the token distribution on
represents the status of the modeled system. The vectorof size represents the number of tokens in each place at step . In addition, , called the initial marking, is the marking at step , which corresponds to the initial status of the system. Therefore, the Petri net represents the structure and initial states of the modeled system.
The transition is enabled at step only when for each in , where and denote the set of input and output places of , respectively. Similarly, and are the sets of input and output transitions of , respectively. An enabled transition can fire. Firing implies the occurrence of an event in the system. The firing of transition removes tokens from each and adds tokens to each . Firing changes the token distribution, which indicates the change in status of the system by the corresponding event. The following mathematical form shows the change in status caused by the firing of transitions at step :
where is a vector of size showing the number of times each fires at step , and is called the firing count vector at step .
For a quantitative analysis of the dynamical behavior of a system, time was introduced to Petri nets (aalst). The timing methods can be categorized into three types: firing duration (FD), holding duration (HD), and enabling duration (ED) (aalst). The FD assigns time to transitions to represent the firing duration. The HD is referred to as place time Petri nets, where tokens are unavailable for firing for a particular period after being located in the place. In the last method, i.e., the ED, a transition cannot be fired for a given period after being enabled. Although other types can be allowed with our method, this paper focuses on timed Petri nets with the FD.
A timed Petri net is a six-tuple , where is the set of time values, and is a function that returns the firing duration time of transition . In this study, we assume that and is a set of natural numbers. A timestamp is also attached to tokens to record the token generation time. In the timed Petri net, the transition is enabled at time when each input place of has more than or equal to tokens, and its time stamp is no more than . By firing at time , the marking is changed according to the same rule as the Petri net described above, except that we attach the time stamp to each output token.
A colored Petri net is an extension of ordinary Petri nets, where we can introduce values (colors) to tokens, guard functions to transitions for their firing conditions, and arc functions to output arcs of transitions to calculate the colors of the output tokens(cpn). Therefore, tokens become much more informative and allow us flexible and efficient modeling.
This extended Petri net is denoted by , where , , and represent a set of places, transitions, and arcs, respectively. shows a set of colors, and is a color function for the places. Only tokens with colors specified by are located in place . In addition, is a set of arc variables, where for . Moreover, is the type of variable . Here, denotes an arc function where returns tokens on for each output place based on the binding values to each of the arc variables. The term represents a guard function that returns a Boolean value to determine whether the corresponding transition is enabled. A marking at step is a mapping of multiple sets of to each place. Transition is enabled at marking with binding when for each input place of ,
where denotes the result of the arc function for arc when we bind for each variable in At marking , is generated by firing with binding :
Colored Petri nets are extremely powerful in the sense that complicated systems can be easily modeled.
3 Binary Quadratic Nets
This section introduces a new Petri net class called binary quadratic nets to represent the Ising and QUBO models. The Ising model is a mathematical model in statistical mechanics. The model is a graph in which the vertices correspond to spins, and the edges interact between spins. Petri nets are also suitable for expressing Ising models because their fundamental components, places, and transitions can naturally represent the states of the spins and their interactions. A token in a place can show the corresponding spin state when we introduce colors to tokens. For QUBO models, the color type of tokens becomes .
In our Petri net modeling-based Ising and QUBO model formulation, Petri net models representing target optimization problems were converted into the corresponding binary quadratic nets. In other words, our proposed method is a net transformation from problem-domain Petri nets into binary quadratic nets. Binary quadratic nets are entirely equivalent to Ising or QUBO models.
3.1 Formal Definition
(Binary Quadratic Net) A binary quadratic net is a colored Petri net denoted by , , , , , :
There is one token in each place, where a place corresponds to a spin variable in the Ising models and a binary variable in the QUBO models, respectively. For simplicity, we denote the value of a variable associated with place by , that is, in the Ising models or in the QUBO models. Note that in the QUBO model, indicates that the place includes a token with color .
Current quantum annealing platforms have various physical topologies, chimera graphs for D-waves, and three-dimensional lattices for CMOS annealers (DaisukeOKU20192018EDP7411). Thus, we need to embed the logical Ising model to a specific graph topology to run the annealing process. However, we do not treat the embedding process in this study because we can use existing embedding algorithms for each platform. Moreover, Fujitsu Digital Annealer allows fully connected topologies. For reference, Fig. 2 depicts the binary quadratic net model for the 2D array Ising model shown in Fig.1.
To combine binary quadratic nets and Ising models, we introduce a measure to represent the state of binary quadratic nets with a marking . We define this as an energy function of binary quadratic nets corresponding to the Hamiltonian in the Ising models.
(Energy Function of Binary Quadratic Net) For a binary quadratic net with a marking , the energy function is defined as follows:
The first summed terms show the energy derived from token existence at each place. The second summed term represents the energy from the interaction between tokens on both sides of each transition. The weight parameters and are defined for and , respectively. The energy function in Definition 2 is equivalent to the Ising model shown in (1), and the energy function is uniquely defined from the corresponding binary quadratic net. Although it is straightforward, we simply summarize this fact as a proposition for the readability of the remaining parts.
For an Ising model, we have an equivalent binary quadratic net in the sense that the energy function is equivalent to the Hamiltonian of the Ising model.
In our binary quadratic nets, transitions represent interactive relations between tokens on both sides. There may be numerous interaction types. We choose interaction types carefully depending on the target optimization problems. In the Appendix, we summarize the primitive interaction for QUBO and the Ising model, respectively. Table A.1 shows the interaction primitives , for two binary variables in , where , and express all possible combinations of the two binary variables. In addition, correspond to the well-known logical functions AND, XOR, OR, and NOR, respectively. Moreover, indicate inconsistency and a tautology, respectively. Table A.2 lists the interactions for the Ising model converted from those in Table A.1 by applying the following relation.
where and denote the marking in the Ising and QUBO models, respectively.
3.2 Binary Quadratic Net Examples
As we describe in Proposition 1, binary quadratic net models are equivalent to the Ising models. We show the binary quadratic net models for well-known graph partitioning and minimum vertex cover problems. These problems are straightforwardly modeled with binary quadratic nets and formulated as marking problems in such nets.
Minimum Vertex Cover: For a given undirected graph , with vertex set and edge set , a vertex cover set satisfies the condition such that every edge in is incident on a vertex in the cover set. The problem is to find the minimum vertex cover set. The problem is known as NP-hard(10.5555/574848).
We choose the QUBO model for the problem, where we generate the place set and transition set of the binary quadratic net under the one-to-one correspondence with the vertex set and the edge set, respectively. We express a vertex cover by a marking, where we place a token with a color of 1 to show that the corresponding vertex is in the vertex cover, and where a color of 0 is not in the cover. To satisfy the feasibility of the vertex cover, we need to avoid a token with a color of 0 in both input places of each transition because the situation does not cover the corresponding edge. We formulate the feasibility of the condition corresponding to in Table A.1 as follows:
To minimize the objective function, we attempt to reduce the number of color 1 tokens in . Therefore, the following penalty function is suitable because only color 1 tokens increase the penalty.
where counts the number of transitions, both of which have a 0 color token. In addition, denotes the number of places colored with 1 to minimize the vertex cover. Moreover, and show the parameters for a trade-off between the constraint and the objective function. The formulation is equivalent to the QUBO model in (10.3389/fphy.2014.00005).
Graph Partitioning: Graph partitioning is an optimization problem, which is explained as follows: Let us consider an undirected graph with vertex set and edge set , where is the number of vertices. The problem is partitioning into two subsets whose sizes are equal to , minimizing the number of edges between the subsets. The problem is known as NP-hard(10.5555/574848).
We reduce the problem to a marking problem in binary quadratic nets, where each vertex in the given graph corresponds to a place in the binary quadratic net, and the marking such that shows the partition. The objective function minimizes the number of edges connecting the two partitioned groups. The transitions have a one-to-one correspondence with the edges. To minimize the objective function, we want to reduce the different color tokens in pairs of places that share the output transition. We design the objective function with in Table A.2 as follows:
Here, outputs 1 only when both places have different color tokens, that is, .
Second, we consider the constraint of the graph-partitioning problem and the equality in the vertex size of both partitioned groups. Concerning logic, we can design the following energy function for the constraint because minimizing the function leads to a satisfaction of the constraint:
The following penalty function shows the total energy function for graph partitioning and is equivalent to the Ising model presented in (10.3389/fphy.2014.00005).
Figure 5 shows the binary quadratic net for the graph partitioning example, where the subnet composed of the places and the transitions connected by the solid arcs correspond to , and the other subnet with all places and transitions connected by the dotted arcs are added based on .
One of the feasible solutions is such that and .
4 Binary Quadratic Net Construction from Problem Domain Petri Nets
In the previous section, we modeled combinatorial optimization problems directly, a minimum vertex cover problem, and a graph partitioning problem with binary quadratic nets and formulated them as marking problems. However, they are exceptional cases. In general, we can meet the conceptual gap between combinatorial optimization problems and target binary quadratic nets. Note that the problem is more serious in the direct formulation of the Ising or QUBO models.
Our approach attempts to minimize the gap by converting problem-domain Petri net models, represented by timed Petri nets and colored Petri nets, into the corresponding binary quadratic nets. This section proposes a method for constructing target binary quadratic nets from problem-domain Petri net models.
4.1 Incremental Construction based on Superposition Principle
A binary quadratic net can be composed incrementally by combining binary quadratic subnets, each corresponding to a constraint or an objective function in the original optimization problem. This superposition principle of the net structure and weight values simplifies the binary quadratic net construction. The composition is straightforward, where we add the weight values on the places and transitions if the same places and transitions in the subnets are to be combined; otherwise, add new places and transitions with their weight values. In the previous section, we observed this process in the graph partitioning example (Example 2). Even though we did not consider the weight values, we combined the two subnets. The following definition formally represents the binary quadratic net composition.
(Binary Quadratic Net Composition) For two given binary quadratic nets = (, , , , , and such that the model types of both nets, Ising or QUBO, are the same, the new binary quadratic net , , , , , is composed based on the following superposition principle:
The following proposition is a straightforward but essential property for validating our method.
(Superposition Property) Let BQN be composed from and . The following properties hold.
where is a scalar, and is the energy function of BQN, such that we replace the weight function with .
In Definition 3, the energy function is defined by the net structure and weight function . Based on the composition rules in (21), (22), and (23), all structural properties of and are transformed into BQN. Based on the rule for the weight function (25), we can confirm that and can be divided into + and + for the common place and transition , respectively. Therefore, the additivity in (26) holds.
The homogeneity property with a degree of 1 is given by Definition 2 if we replace the weight function with .
Owing to the additivity and homogeneity properties in Proposition 2, we have the following corollary:
Binary quadratic nets can be composed by the incremental application of Definition 3, in which we can scale the weight function in subnets with a constant factor.
In our approach, we construct a target binary quadratic net based on the incremental compositions of binary quadratic subnets. Each subnet is converted from a property of the Petri net model representing the optimization problems. We call the Petri net problem-domain Petri net. The properties of the problem-domain Petri net models are expressed with marking or firing sequences. To focus on markings or firing sequences, we employ marking-based or firing-based constructions. In the following subsections, we assume that binary quadratic nets are QUBO models unless otherwise stated, but can be converted into Ising models by using the conversion rule (58).
4.2 Marking-based Construction
Let us denote a problem-domain Petri net for a target optimization problem by with and . Marking represents a set of tokens in place at time step . In addition, is a vector of size , showing a token distribution on the Petri net at time step . The marking trajectory of a Petri net model, denotes the status changes triggered by the firing of transitions.
In the marking-based construction of binary quadratic nets, we represent a marking trajectory of the problem domain Petri net by the place set of the target binary quadratic net. If each place in the problem-domain Petri net has at most one token at any time and the maximum step is , we initially prepare the following place set of the target binary quadratic net.
If is a natural number, that is, more than one token or a color token with a natural number value is possible in each place in , we need to prepare more places because each place in the binary quadratic nets has one token with color in or .
where is the possible maximum value of . At the same time, we require a one-hot constraint because only one place between for for each and should be marked with , and the others should be marked with :
As the energy function representation, we need the following:
Note that we can obtain the corresponding binary quadratic net from an energy function.
The boundedness is an essential characteristic to ensure avoiding an overflow in the system behavior. In the Petri net theory, we can express this characteristic as follows.
where is the upper bound for . We can convert the constraint into a binary quadratic net and represent it as the energy function under the one-hot constraint (31):
Function (33) is sufficient for the equality constraint , but not for the upper bound. Therefore, we improve the function by introducing ancilla places, , . This technique is well known in the Ising model formulation (10.3389/fphy.2014.00005), (pyqubo).
The upper bound constraint appears in numerous optimization problems. The knapsack constraint is a well-known example of this constraint for the knapsack place. We can also express the boundedness of specific places by removing the other places from the function.
The boundedness shown in (32) denotes inequality constraints. The invariance leads to equality constraints based on markings. Note that the invariance is different from the structural invariance of the net theory. The behavioral invariance requires that the total weighted sum of the tokens be equal among the firing sequences. The sum may be the cost required for the resource to operate the system. The following constraint shows that the weighted sum becomes for all .
We can convert the constraint into a binary quadratic net and represent it as the energy function under the one-hot constraint (31):
4.3 Firing-based Construction
A firing count vector represents the firing counts of each transition in a Petri net model at step , where denotes the firing counts of transition at time step . Elements of are usually one of ; that is, transitions can fire only once at each time step. We call this single firing restriction single-server semantics in Petri net theory. However, any natural number values (allowing more than one) are also possible with infinite server semantics. For the single-server semantics, we generate a place set in the target binary quadratic net, where corresponds to :
Similarly, we can extend the single-firing semantics to the times firing semantics. Note that an infinite number of firing times is possible mathematically but impossible practically. Thus, we restrict the maximum number reasonably to .
As the energy function representation, we need the same one-hot constraint (31):
4.3.1 Resource Conflict
If transitions and have a common input place and conflict with a single token, cannot be allowed. Therefore, we need the following energy function:
where is a set of conflict transitions. In addition, can be extracted from the problem domain Petri net.
In timed Petri net models, we need to consider the firing duration. Let with and be a timed Petri net, where is a function that returns the firing duration. We then extend the resource conflict in the stepwise firing into the timed firing.
where is the timed conflict set such that conflict transitions cannot fire until the firing duration is complete. We can obtain from the given timed Petri net.
4.3.2 Firing Count
In some applications, we must specify the number of firing occurrences for each transition.
where is the specified number of firings of .
We can convert the constraint into a binary quadratic net and represent it as the energy function under the one-hot constraint (31):
Note that in the binary quadratic net corresponds to in the problem domain Petri net.
Assuming that each transition should fire exactly once during , , …, , the function (42) can be represented as follows:
In practical cases, this constraint is commonly used.
4.3.3 Precedence Relation
Let us assume again that each transition should fire exactly once during , , …, . We consider the precedence relation between the firing of transitions. If precedes structurally to , that is, and , the following penalty function is required.
where is the set of precedence relations. Here, can be extracted from the problem domain Petri net. Similar to the firing conflict, we consider timed Petri nets by introducing the firing duration.
4.4 Application Example 1 (Marking-based Construction): Traveling Salesman Problems
We consider the traveling salesman problem as an example of a marking-based construction. Let with and be a timed Petri net, where is a function that returns the firing duration. Each place corresponds to a city, and the transition denotes the path from to . The firing of indicates the movement from to . The initial marking includes only one token at .
Because the salesman visits each city only once from the problem definition,
In addition, because the salesman should be at one place for each step,
The objective function is to minimize the total time to visit all cities and return to the starting place.
By converting the constraints and objective function in the problem-domain Petri net into the binary quadratic nets, we have
The total binary quadratic net is as follows:
where , , and denote the scale factors of the corresponding subnets. The formulation is equivalent to the QUBO model in (10.3389/fphy.2014.00005).
4.5 Application Example 2 (Firing-based Construction): Job-shop Scheduling Problems
We implemented our method by incorporating CPNTools (cpn), SNAKES(10.1007/978-3-319-19488-2_13), and PyQUBO (pyqubo). As an example, we modeled a simple job-shop scheduling problem with three jobs, four tasks per job, and three shared resources. Figure 6 shows the model drawn by the GUI software of CPNTools. Three sequential systems,
represent the jobs, and places represent the resources. Our software reads the model, represents the problem-domain Petri net model with SNAKES, and then extracts the necessary information to convert it into the target binary quadratic net.
The following are the energy functions of the binary quadratic subnets corresponding to the precedence relation constraint between tasks in each job, the resource conflict constraint for each shared resource, and the firing count constraint for each task.
The precedence set in (54) is easily constructed by extracting the precedence relation in each sequential system (job). We can extract the timed conflict set in (55) from the connection between the transitions and resource places and the firing duration . MaxTime is the only parameter we need to set before the optimization process. This value indicates the delivery time deadline. The total binary quadratic net is as follows:
Note that this formulation is used to obtain a feasible solution; however, we can extend this to the optimization by incorporating it with a binary search to find the minimum MaxTime. This formulation is equivalent to the QUBO model in (venturelli2016quantum). Figure 7 shows the schedule obtained through the annealing process using a Fujitsu Digital Annealer.
This paper proposes a Petri net modeling approach to the Ising model formulation for quantum annealing. Although our method requires users to model their optimization problems with Petri nets, this process can be carried out in a relatively straightforward manner if we know the target problem and the simple Petri net modeling rules. Therefore, we can drastically relax the difficulty of the Ising model formulation. We implemented our method with Python incorporated using well-known Petri net tools, CPNTools and SNAKES. We can automatically generate the Ising models for optimization problems, such as scheduling pr