A Study of Ising Formulations for Minimizing Setup Cost in the Two-Dimensional Cutting Stock Problem

: We proposed the method that translates the two-dimensional CSP for minimizing the number of cuts to the Ising model. After that, we conducted computer experiments of the proposed model using the benchmark problem. From the above, the following results are obtained. (1) The proposed Ising model adequately represents the target problem. (2) Acceptance rates were as low as 0.2% to 9.8% and from 21.8% to 49.4%. (3) Error rates from optimal solution were as broad as 0% to 25.9%. For future work, we propose the following changes: (1) Improve the Hamiltonian for constraints. (2) Improve the proposed model to adjust more complex two-dimensional CSP and reduce the number of spins when it deals with large materials and components. (3) Conduct experiments using a quantum annealer.


Introduction
In several industries, such as paper, steel, and textile, items of specified sizes are cut out from the raw material placed on a flat surface. Problems aimed at minimizing the cost of such cutting work are referred to as two-dimensional cutting stock problems (CSPs) and two-dimensional bin packing problems (BPP), on which several studies have been carried out so far.
In particular, when positioning the items to be cut from raw material of fixed width W and variable height H, if the purpose of the problem is to minimize the height H, it is formulated as the strip packing problem (SPP), and various heuristic solutions have been proposed. E.G. Coffman Jr. et al. proposed the typical solution method, the first-fit decreasing height (FFDH) method and the next-fit decreasing height (NFDH) method [1]; C. M. Valenzuela et al. arranged FFDH and proposed the best-fit decreasing height (BFDH) method [2]; A. Lodi, et al., also arranged FFDH and proposed the floor-ceiling no rotation (FCNR) method [3]. B. Baker et al., in addition to the FFDH method, proposed solutions including the bottom-left method [4], in which items are introduced from the upper right of a material and then stacked from the bottom left, and heuristic approaches such as the bottom-left fill method [5], which adds an algorithm for reducing the gap between items. There are also solutions proposed using a metaheuristic method, such as the hill-climbing method [6], which is a typical local search method; simulated annealing [7]; the genetic algorithm (GA) [8][9][10]; and the native evolution algorithm [11], which eliminates GA crossover and evolves only by mutation. The SPP and BPP presume that all members have been arranged as a precondition; however, solutions [12] have been proposed using a knapsack problem approach when not all members have been arranged for selecting members that are given gain so that their sum is maximized, although many of these deal with materials as cost objects. These previous studies are aimed at maximizing the yield of raw materials as the cost target. However, it is not just the material cost that must be considered [13]. For example, at a manufacturing site where special sheet-like materials are cut, a setup change cost is incurred for each cutting, which significantly impacts the total cost. Hence, a technique that can facilitate a reduction of the setup change cost is desired. In this paper, we propose an allocation method aimed at minimizing the number of setup changes by applying conventional methods, such as the FFDH and binary tree methods, for manufacturing sites that cut special sheet-like material. This method has helped achieve a reduction in the number of setup changes while ensuring a high yield [14].
Because a two-dimensional CSP is considered N P-hard, combinatorial explosion occurs when the complexity of the problem and, in turn, the scale of the solutions increases, leading to an exponential increase in the calculation time. In contrast, the actual problem represented by the two-dimensional CSP has a short lead time and requires a high-speed solution.
Quantum computers are next-generation computers that use quantum mechanics and are implemented through two types of methods: the quantum annealing and quantum gate methods. The quantum annealer is considered to be suitable for solving combinatorial optimization problems. Therefore, various combinatorial optimization problems have been studied using quantum annealing. In a previous study on the BPP, Terada et al. (2018) performed Ising modeling using sequence pairs for maximizing the yield of the rectangular packing problem [15]. Furthermore, Kida et al. (2000) performed Ising modeling with a formulation based on the NFDH method for the glass processing industry; they compared the performances of the quantum annealer and the classical bifurcation machine [16]. However, Ising models aimed at reducing setup change cost (= number of cuts) in the two-dimensional CSP have not yet been created.
Therefore, in this study, we propose an Ising model for the two-dimensional CSP that minimizes the number of cuts, we perform simulation experiments on benchmark problems, and we compare the results with the optimal solution determined using a solver [17].

Ising Model
When solving a combinatorial optimization problem using a quantum annealer, it is necessary to transform it into an equation referred to as the Ising model. The Ising model is a statistical mechanics model that represents the properties of a magnetic material, and the state of the magnetic material is evaluated by calculating the upward or downward state of the spin. Spin has two forces: a local field that acts only on that spin, and an interaction that acts on the connected spins. The state of the spin is updated by these forces so that the total energy is stabilized. When the entire system is in the most stable state, the Hamiltonian H, which represents the total energy of the Ising model, is the minimum (ground state). Equation (1) expresses this Hamiltonian. Here, it should be noted that s i and s j are the directions (±1) of the spins i and j, respectively, J ij represents the interaction, h i is the local field of spin i, E is the number of edges, and V is the number of spins:

Target Overview
We consider a two-dimensional CSP as the target problem. However, if applied as it is, the number of spins in the quantum annealer may be insufficient, as there are several constraints. Hence, some simplification is required. Figure 1 presents an example of the allocation and specific conditions. The main conditions of this model are as follows: 1.
There are infinite base materials with fixed widths and lengths.

2.
The items shall be allocated from the lower left of the base material.

3.
The allocated items must not protrude from the base material. 4.
The allocated items must not overlap.

5.
All the items are allocated from one of the base materials only once. 6.
The items must not be rotated. 7.
Allocation in the vertical direction of the base material is known as a step, and the length of each step is the length of the longest piece in that step. 8.
Items shall not be allocated vertically in each step. 9.
If there are two or more items of the same length in each step, their horizontal cutting is carried out all at once. 10. If the allocated piece in each base material touches the upper side of the base material, horizontal cutting is omitted. 11. If the total width of the items allocated in each step matches the width of the base material, the vertical cutting of the rightmost piece is omitted. Table 1 lists the variables used in this paper. Furthermore, there are five types of spins as follows:

Ising Modeling
Total length of longest items of each step of base material i is p |0 Otherwise Equation (3) represents the number of cuts in the horizontal direction and satisfies Constraints 9 and 10 listed in Section 3. For each piece length type at each step, the number of cuts is incremented, and when the sum of lengths of the longest items at each step matches the length of the base material, the number of cuts is omitted once.
Equation (4) expresses the number of cuts in the vertical direction and satisfies Constraint 11 listed in Section 3. The number of cuts is incremented for each piece at each step, and when the total width of the items allocated to each step matches the width of the base material, the number of cuts is omitted once.
H A is expressed as follows: Equation (5) expresses a constraint on the number of times items allocated, and satisfies Constraint 5 listed in Section 3. Its value is the lowest when all the items are allocated to one of the base materials only once.
H H and H W are expressed as follows: Equation (6) represents the constraint for preventing the items from protruding vertically from the base material, and its value is the lowest when there is only one spin with sht = 1 in each base material and the subscript p of the spin with sht = 1 matches the sum of the subscripts m of the sh at each step. Equation (7) expresses the constraint for preventing the items from protruding horizontally from the base material, and its value is the lowest when there is only one spin with sw = 1 at each step for each base material, and the subscript l of the spin with sw = 1 matches the total width of the items allocated to the target step.
Equations (6) and (7) satisfy Condition 3 listed in Section 3. H Htype is expressed as follows: Although Equations (3)- (7) can be used to represent a two-dimensional CSP that minimizes the number of cuts, the spin shc cannot always correctly represent the piece length type at each step. Therefore, a constraint is set on shc by Equation (8). This equation yields the minimum value when the values of all the variables s i,j,k with h k = n become 0 at each step for each base material, leading to shc = 0, and when one of these values becomes 1, leading to shc = 1.
H Hlongest is expressed as follows: n · shc i,j,n sh i,j,m|m=n ∀i, j (9) Similar to shc, the spin sh cannot always correctly represent the longest piece at each step only through Equations (3)- (8). Therefore, a constraint is set on sh by Equation (9). This equation yields the minimum value only when the sh with the subscript m that matches the longest of the items allocated at each step for each base material becomes 1.

• Experimental conditions
To confirm the validity of the proposed Ising model, we performed a computer experiment. Because the target problem is a benchmark problem [18], it is considered as an instance that satisfies the conditions specified in Table 2. However, for the sake of simplicity, the length of the base material Bin_H is considered as a variable, and it is assumed that a cut is always made in line with the longest piece among the allocated items. Thus, the full Ising model of this experiment is represented by Equation (10): For the experiment, Equation (10) is converted into a QUBO matrix, which is known to have a one-to-one correspondence with the Ising model. The QUBO matrix is applied to QBsolv and neal. QBsolv can optimize calculations through tabu search, and neal can optimize calculation using simulated annealing [19]. Table 3 shows each parameter setting that was set by the preliminary experiment. Under these experimental conditions, the number of logical spins are s i,j,k 10 × 10 × 20 = 2000, sw i,j,l : 11 × 10 × 20 = 2200, sh i,j,m : 11 × 10 × 20 = 2200, shc i,j,n : 10 × 10 × 20 = 2000 and sht i,p : 11 × 20 = 220. Since sw i,j,l and sh i,j,m need a spin to express height 0, these needs 11 spins for height.

• Experimental results
We now consider the acceptance rate and accuracy of the experimental results. Furthermore, we calculate the optimal solution of the target benchmark problem using Google OR-tools [20] and compare the two optimizer's results.
Tables 4 and 5 indicate that the relative error between the best solution using the proposed Ising model and the optimal solution using Google OR-tools were 0-25.9%, respectively, obtained using both QBsolv and neal. The calculation times are from 0.0 to 3.0 s in all instances. The optimal solutions were obtained in Instances 9 and 3 in QBsolv and neal, respectively. The acceptance rates were 0.2-9.8% and 21.8-49.4%, respectively. The difference in the acceptance rate is the algorithmic difference in each optimizer. As neal used SA, the acceptance rate of the average became higher because the acceptance rate of the initial search was high. About QBsolv, several solutions that violated the constraints were calculated. These optimizers are high-performance, but there is little room for parameter adjustment. For the improvement of the solution, the parameters can be adjusted or the spin's design of the Ising model can be modified. Although the parameters H Htype and H Hlongest of the Ising model are constraints, their minimum values change depending on the problem and other parameters, and they are effectively objective functions that impose constraints on other parameter settings. If a model could be constructed such that the values of these Hamiltonians could reach a minimum of 0, the range of acceptable values for the other parameters would become wider, and flexible response would become possible. Further, the reason for the difference in the acceptance rate depending on the instance, is that if the number of items with widths exceeding 50% of the width of the base material (number of high-width items) is high, the ratio of the solution space to the search space becomes small. The probability of obtaining an acceptable solution within the execution time becomes low. Figures 3 and 4 present the distribution maps of the acceptance rates and the numbers of high-width items at each instance obtained using QBsolv and neal, respectively. The broken lines in the figures represent the linear approximation curves. These figures indicate that the acceptance rate of the solution decreases proportionally with the number of high-width items.

Conclusions
In this paper, we propose an Ising model with the aim to minimizing the number of cuts in a two-dimensional CSP. This is first trying a formulation of a Ising model. In this study, five kinds of spins are prepared, and these change their state with interaction. The proposed Ising model worked appropriately, and as the number of cuts decreased, the energy consumed also reduced. However, in the experiment, the acceptance rates were 0.2-9.8% for QBsolv and 21.8-49.4% for neal. The relative error of the optimal solution was 0-25.9%, which indicates a wide range. Although the relative errors are depend on an instance condition, the parameter setting and the Ising modeling have to improve for versatility. When minimum value of parameters regarding Hamiltonian H Htype and H Hlongest become 0, the other parameter's values are able to change more flex. Furthermore, in the proposed method, the constraints for the allocation of items in the vertical and horizontal directions were realized by preparing spins based on the size of the base material. However, as the size of the base material increases, the number of spins to be secured also increases, making it difficult to apply to the quantum annealer. The following can be listed as future challenges: • Improving the constraints used in the Ising model; • Adjusting the parameters further; • Considering the omission of vertical cutting; • Devising an Ising model that can be applied to situations that use large-sized base materials; • Application to the quantum annealer.