# RAMP: Reaction-Aware Motion Planning of Multi-Legged Robots for Locomotion
in Microgravity

Warley F. R. Ribeiro1, Kentaro Uno1, Masazumi Imai1, Koki Murase1 and Kazuya
Yoshida1 This work was not supported by any organization1Warley F. R. Ribeiro,
Kentaro Uno, Masazumi Imai, Koki Murase and Kazuya Yoshida are with the
Department of Aerospace Engineering, Graduate School of Engineering, Tohoku
University, Sendai, Japan. {warley, imai.masazumi.p2,
murase.koki.q2}[at]dc.tohoku.ac.jp and {unoken, yoshida.astro}[at]tohoku.ac.jp

###### Abstract

Robotic mobility in microgravity is necessary to expand human utilization and
exploration of outer space. Bio-inspired multi-legged robots are a possible
solution for safe and precise locomotion. However, a dynamic motion of a robot
in microgravity can lead to failures due to gripper detachment caused by
excessive motion reactions. We propose a novel Reaction-Aware Motion Planning
(RAMP) to improve locomotion safety in microgravity, decreasing the risk of
losing contact with the terrain surface by reducing the robot’s momentum
change. RAMP minimizes the swing momentum with a Low-Reaction Swing Trajectory
(LRST) while distributing this momentum to the whole body, ensuring zero
velocity for the supporting grippers and minimizing motion reactions. We
verify the proposed approach with dynamic simulations indicating the
capability of RAMP to generate a safe motion without detachment of the
supporting grippers, resulting in the robot reaching its specified location.
We further validate RAMP in experiments with an air-floating system,
demonstrating a significant reduction in reaction forces and improved mobility
in microgravity.

##  I Introduction

Mobile robots are part of the essential technologies necessary to increase
human presence and exploration of extraterrestrial bodies. Space has, however,
several extreme circumstances different from the regular mobile robots on
Earth, such as drastic temperature variation, high radiation, and low gravity.
A more specific condition is microgravity, present in two crucial situations:
orbiting space stations and exploration of asteroids.

By exploring Small Solar System Bodies (SSSBs), such as asteroids and comets,
we can obtain significant scientific data to understand the process of
formation and evolution of our Solar System. Those celestial objects also
showed substantial amounts of water and organic matter, which could give us
insights into the process of life formation. SSSBs can also threaten our
planet if they cross Earth’s orbit, and impact mitigation missions need to be
a concrete alternative for our safety. And lastly, SSSBs are also rich in
valuable materials for use on Earth and further space development [1]. All
those scenarios need to rely on mobile robots capable of performing
interactive activities on the surface of the SSSB, such as in-situ
investigation, terrain mapping, and resource sampling and analysis.

Another use of mobile robots is assisting astronauts in orbiting stations,
such as the International Space Station (ISS). Robots could achieve simple
tasks, such as helping with setting up experiments and recording the results
of the tests to optimize the valuable crew time [2]. High-capable robots can
also perform more complex activities, such as maintaining the components
necessary for habitable conditions, transporting and organizing the cargo
inside the station and assembling new structures and modules.

Figure 1: Scheme for the Reaction-Aware Motion Planning (RAMP). The trajectory
of the swing limb $\boldsymbol{x}_{e}$ is optimized to minimize motion
reactions using a Low-Reaction Swing Trajectory (LRST) algorithm. The motion
of the base $\boldsymbol{x}_{b}$ is generated to balance the swing reactions
with the Whole-Body Swing Momentum Distribution (MD) technique. A zero
constraint velocity is applied for all supporting limbs to mitigate motion
reactions at the contact points. RAMP can mitigate motion reactions to reduce
the risk of detachment from the surface.

###  I-A Related Works

The hopping mobility technique is the only robotic surface locomotion method
employed in actual missions to SSSBs up to this date. A robot uses an internal
actuator, such as a flywheel, to create a striking reaction force from the
surface to generate a jump with a ballistic trajectory. Due to the
microgravity condition, small and lightweight actuators can produce high hops
within the escape velocity of the celestial body [3]. The twin rovers MINERVA-
II-1 [4] and the lander MASCOT [5] successfully demonstrated this method in
the Japanese Hayabusa2 mission to the asteroid Ryugu [6]. A fourth robot,
MINERVA-II-2, was also aboard the Hayabusa2 mission, equipped with four
different actuators to produce hops [7]. However, the second MINERVA-II rover
did not accomplish its mobility objective due to a malfunction in the main
computer. The development of hopping robots continues to evolve, with more
robots incorporating new mobility concepts and strategies [8, 9].

When we consider past and current robots for assistance in a space station,
the primary locomotion method is that of the free-flying robots. By utilizing
propellant gases or fans to generate the reaction forces, small robots can
move inside the tridimensional space of the ISS thanks to the negligible
gravity. Astrobee [10], CIMON [11] and Int-Ball [12] are some successful
examples of free-flying robots with autonomous mobility and task-performing
capabilities operating in the ISS.

Bio-inspired legged robots were proposed as a viable general solution for
mobility in microgravity [13]. Considering uncontrollable hops inside a space
station can be dangerous for the astronauts, and propellant-based free-flying
motion is limited for SSSBs without an atmosphere, those robots do not work in
any microgravity locomotion context. Also, fast and accurate mobility is
challenging with robots based on free-flying and hopping techniques. Multi-
legged robots with surface attachment capabilities can achieve precise and
stable locomotion in space stations and small celestial bodies, although there
is additional complexity regarding fabrication and control.

Robonaut 2 was a humanoid robot with climbing legs that flew to the ISS as a
robot to assist astronauts but did not demonstrate mobility in microgravity
[14]. As an example of multi-legged robots for microgravity locomotion, we
have the LEMUR series developed by NASA JPL, where the latest version can have
microspine grippers for rough surfaces of SSSBs, or gecko-inspired adhesive
grippers for smooth surfaces in the ISS [15]. A similar robot with microspines
was also proposed, including stochastic grasping modelling [16]. Another
example is ReachBot, which includes extendable limbs to improve the
reachability of a compact designed robot [17].

###  I-B Objectives and Contributions

Besides the complexity of a system with multiple links and developing a
reliable grasping mechanism, the safe mobility of multi-legged robots in
microgravity is still an issue. Although the lack of gravity facilitates
static equilibrium, quasi-static legged mobility is not a suitable solution
for space robots, as the time constraints are significantly severe. It is
required to perform dynamic evaluations to incorporate the inertial
acceleration and forces, as the robot’s motion can cause the detachment of the
grippers. Chacin et al. proposed a compliant motion control considering the
friction forces conditions for equilibrium [18]. Yuguchi et al. proposed a
different approach, using a reactionless method for a bipedal robot to reduce
the ground reaction forces [19]. However, precise friction characteristics are
not available before exploration missions, complicating the estimation of
forces acting on the robot. Moreover, reactionless control was demonstrated
only for a single contact point and could cause considerable changes in the
robot’s pose, inducing kinematic singularities.

In this paper, we propose a novel Reaction-Aware Motion Planning (RAMP) to
improve the safety of multi-legged locomotion in microgravity by decreasing
the motion reactions. We also aim to provide low-reaction alternatives to the
reactionless control to prevent singular configurations. First, we introduce
our approach to mitigate the change in the robot’s momentum in two different
aspects. The first is to minimize the reactions generated from the swing leg
motion. And the second is the whole-body momentum distribution to partially or
fully compensate for the swing momentum, constraining the velocity of contact
points to zero. Later, we verify the increased performance of the proposed
motion planning in dynamic simulations of realistic scenarios of multi-legged
mobility in microgravity. And finally, we validate the approach in an emulated
microgravity facility.

##  II Reaction-Aware Motion Planning

In a microgravity condition, we can neglect the effects of the gravity force
when considering the elements that can cause the detachment of the grippers.
Assuming that only the ground reaction forces are acting on the robot, we can
limit the detachment analysis to the maximum pulling force the gripper can
endure before undesired separation from the surface. And the only source of
the pulling forces under this microgravity condition is the robot’s self-
movement, either due to the swing of the legs to the next grasping position or
the main body to reach the desired location. In this study, we define _Motion
Reaction_ as the pulling or pushing action of the supporting grippers induced
by the robot’s movement.

In the case of quasi-static walking, the motion reactions are virtually null,
and no surface separation should occur. But responsive and agile robots are
needed for future space exploration and human assistance, and dynamic
locomotion is necessary to achieve our space utilization goals. The proposed
Reaction-Aware Motion Planning (RAMP) addresses the safety locomotion issue by
generating a swing motion for the leg with the minimum motion reaction.
Additionally, the robot performs an additional movement using the remaining
degrees of freedom to drive its main body, distributing the generated motion
reactions. The motion plan also adds zero velocity constraint for the
supporting grippers, assuring that no excessive pulling motion happens. The
momentum distribution with zero velocity at the contact points further
decreases motion reactions and ground reaction forces that could cause
detachment. Fig. 1 shows an outline of the proposed motion planning.

###  II-A Momentum Variation Approach

A simple way to quantify the pulling action of the grippers due to motion
reactions is by computing the total internal force induced by the planned
robot’s movement. This force is the derivative of the resultant momentum of
the motion generated by the robot. Equation (1) describes the momentum
$\boldsymbol{\mathcal{L}}$ of an articulated robot with multiple rigid links
and $n$ legs under microgravity. $\boldsymbol{H}_{b}$ and
$\boldsymbol{H}_{bm,i}$ are the inertia matrices of the base link and the
base-manipulator coupling, respectively. The vector with the tridimensional
position of the base link is $\boldsymbol{x}_{b}$, and $\boldsymbol{\phi}_{i}$
is the vector with the angular position of all articulated joints.

$\displaystyle\boldsymbol{\mathcal{L}}=\boldsymbol{H}_{b}\dot{\boldsymbol{x}}_{b}+\sum_{i=1}^{n}\boldsymbol{H}_{bm,i}\dot{\boldsymbol{\phi}}_{i}$
(1)

Therefore, by minimizing the variation in the robot’s momentum, we can reduce
the motion reactions causing the pulling movement of the grippers, which is
the primary reason for surface separation under microgravity, increasing the
locomotion stability.

###  II-B Low-Reaction Swing Trajectory

The swing motion of a leg from its current position to the location of the
next grasping point needs to be in the feasible space of the robot and should
avoid collisions with the environment. We include a third condition to the
trajectory generation algorithm: minimizing momentum variation to lower the
motion reactions, preventing the gripper detachment in microgravity
conditions.

A low-reaction trajectory generation was proposed based on optimizing the
coefficients of a polynomial curve to find the trajectory with minimum
momentum variation, reaching the desired step height to avoid collisions [20].
Equation (2) shows the details of the optimization problem that constrains the
end-effector trajectory $\boldsymbol{x}_{e}(t)$ from the initial time $t_{0}$
to the final time $t_{f}$ by a 7th-order Bezier curve of coefficients
$\boldsymbol{a}_{i}$. The optimization is also subject to the feasible angular
position of joints $\boldsymbol{\phi}$ bounded by the minimum and maximum
values, $\boldsymbol{\phi}_{\min}$ and $\boldsymbol{\phi}_{\max}$
respectively.

$\displaystyle\min_{\boldsymbol{a}_{3},\boldsymbol{a}_{4}}$ $\displaystyle
J_{1}({\dot{\mathcal{L}}}(t))+J_{2}(x_{z}(t))$ (2) s.t.
$\displaystyle\boldsymbol{x}_{e}(t)=\sum_{i=0}^{7}B_{i}(t)\boldsymbol{a}_{i}\quad\quad
t_{0}\leq t\leq t_{f}$
$\displaystyle\boldsymbol{\phi}_{\min}\leq\boldsymbol{\phi}(t)\leq\boldsymbol{\phi}_{\max}$

The objective functions $J_{1}$ and $J_{2}$ incorporate the momentum variation
and the step height $h$ of the desired trajectory, including selectable
coefficients $k$ to adjust the weight of each term.

$\displaystyle J_{1}({\dot{\mathcal{L}}}(t))=$ $\displaystyle
k_{1}\max({\dot{\mathcal{L}}}(t))$ (3) $\displaystyle J_{2}(x_{z}(t))=$
$\displaystyle k_{2}\left|h-\max(x_{z}(t))\right|$ (4) $\displaystyle+$
$\displaystyle k_{3}\left|h-\operatorname{ave}\left(x_{z}(t)\right)\right|$

###  II-C Whole-Body Swing Momentum Distribution

The Low-Reaction Swing Trajectory (LRST) generation method successfully
reduces motion reactions, improving the safety of legged locomotion in
microgravity. However, it does not suffice to guarantee stability during
mobility. Therefore, we propose adding the Whole-Body Swing Momentum
Distribution (MD) technique to avoid pulling forces in the grippers, forming
the RAMP framework.

To distribute the momentum from the swing motion and ensure zero velocity
conditions for the supporting grippers, we divide the total momentum from (1)
into the momenta of the base $\boldsymbol{\mathcal{L}}_{\text{base}}$,
supporting legs $\boldsymbol{\mathcal{L}}_{\text{sup}}$, and swing legs
$\boldsymbol{\mathcal{L}}_{\text{swing}}$, as shown by (5), where
$n_{\text{sup}}$ and $n_{\text{swing}}$ denote the number of supporting legs
and the number of swing legs, respectively.

$\displaystyle\boldsymbol{\mathcal{L}}=$
$\displaystyle\boldsymbol{\mathcal{L}}_{\text{base}}+\boldsymbol{\mathcal{L}}_{\text{sup}}+\boldsymbol{\mathcal{L}}_{\text{swing}}$
(5) $\displaystyle=$
$\displaystyle\boldsymbol{H}_{b}\dot{\boldsymbol{x}}_{b}+\sum_{i=1}^{n_{\text{sup}}}\boldsymbol{H}_{bm,i}\dot{\boldsymbol{\phi}}_{i}+\sum_{j=1}^{n_{\text{swing}}}\boldsymbol{H}_{bm,j}\dot{\boldsymbol{\phi}}_{j}$

We combine (5) with the forward kinematics equation, as shown in (6), to merge
the momenta of the base and support legs into two terms as a function of the
velocity of the base and the supporting grippers $\dot{\boldsymbol{x}}_{e,i}$.
By forcing the zero velocity condition to the supporting grippers, we obtain
the relation in (7), where $\boldsymbol{J}_{b}$ indicates the Jacobian matrix
of the base, and $\boldsymbol{J}_{mi}$ represents the Jacobian of the
manipulator part.

$\dot{\boldsymbol{x}}_{e,i}=\boldsymbol{J}_{b}\dot{\boldsymbol{x}}_{b}+\boldsymbol{J}_{mi}\dot{\boldsymbol{\phi}}_{i}$
(6)

$\boldsymbol{\mathcal{L}}=\left(\boldsymbol{H}_{b}-\sum_{i=1}^{n_{\text{sup}}}\boldsymbol{H}_{bm,i}\boldsymbol{J}_{mi}^{+}\boldsymbol{J}_{b}\right)\dot{\boldsymbol{x}}_{b}+\sum_{j=1}^{n_{\text{swing}}}\boldsymbol{H}_{bm,j}\dot{\boldsymbol{\phi}}_{j}$
(7)

If we assume that there is no change in the momentum, we can compute the
velocity of the base that compensates for a fraction $\alpha$ of the momentum
generated from the swing motion, as shown in (8). In summary, we distribute
the swing momentum to the base link and supporting legs, avoiding that high
motion reactions cause the detachment of the gripper by securing zero velocity
for supporting end-effectors.

$\left(\boldsymbol{H}_{b}-\sum_{i=1}^{n_{\text{sup}}}\boldsymbol{H}_{bm,i}\boldsymbol{J}_{mi}^{+}\boldsymbol{J}_{b}\right)\dot{\boldsymbol{x}}_{b}=-\alpha\sum_{i=1}^{n_{\text{swing}}}\boldsymbol{H}_{bm,j}\dot{\boldsymbol{\phi}}_{j}$
(8)

If the momentum distribution factor $\alpha$ is equal to zero in (8), we have
a base velocity equal to zero, which is the same as not performing any
momentum distribution. By selecting an appropriate value for the momentum
distribution factor greater than 0 and smaller than 1, it is possible to
reduce a significant amount of motion reactions in the gripper, keeping the
base attitude within safe limits to avoid singularities. We name this
procedure Reaction-Aware Motion Planning with Partial Momentum Distribution
(RAMP-PMD). If $\alpha$ equals 1, the swing momentum is fully distributed to
the remaining parts of the robot, mitigating the motion reactions. We call
this technique RAMP with Full Momentum Distribution (RAMP-FMD).

With RAMP, we calculate the position for the swing end-effectors using LRST in
(2). Moreover, we compute the necessary base velocity in (8) to distribute the
swing reactions based on the momentum calculated from the moving legs’ speed.
Therefore, considering the zero velocity condition for the supporting
grippers, RAMP outputs the trajectory for all end-effectors and the base,
considerably mitigating the detachment risk for mobility in microgravity.

##  III Simulation Studies

Testing new technologies and strategies for locomotion in microgravity is
extremely difficult due to the limitations in sending equipment for tests in
space. Therefore, numerical simulations are essential in developing new
mobility approaches for robots to explore SSSBs and assist in space station
activities.

In this paper, we use ClimbLab as a simulation platform. ClimbLab is a dynamic
simulator developed in MATLAB for climbing-legged robots [21]. ClimbLab uses
SpaceDyn to compute the dynamics of multi-body robots and other parameters,
such as inertia and Jacobian matrices [22].

The main criterion for failure analysis in the simulations is the detachment
of the gripper due to excessive pulling forces. In this model, if a ground
reaction force acting on the robot pulls it towards the surface, exceeding the
maximum holding force of the gripper, a detachment occurs. We use a compliant
contact model to compute ground reaction forces, assuming the end-effector as
a single point and the initial contact point as the neutral position.

###  III-A Quadrupedal Robot on an Asteroid (Rough Surface)

For this simulation study, we assume a rough surface to emulate an asteroid
with a gravity of $10^{-6}$ G. The surface coefficients for stiffness and
damping are 4000 N/m and 1 Ns/m, respectively. As for the robot model, we
assume the inertial and length parameters of Hubrobo, a quadruped insect-type
robot with three controllable articulations per leg [23]. Details of the
parameters of the robot are in Table I. For this simulation, we assume a
considerably small holding force of 0.9 N, as the fragile and porous surface
of an SSSB requires weak grasping.

TABLE I: Inertial Parameters of Quadrupedal Robot Link | Size (mm) | Mass (g) |  | Inertia (kgm2) |
---|---|---|---|---|---
|  |  | $I_{xx}$ | $I_{yy}$ | $I_{zz}$
Base | 108 x 108 | 441.7 | 9.2e-4 | 1.4e-3 | 1.3e-3
1 | 28.5 x 17.5 | 205.6 | 4.2e-5 | 2.8e-5 | 2.8e-5
2 | 107 | 27.3 | 1.5e-5 | 2.4e-5 | 3.7e-5
3 | 143 | 220.2 | 2.2e-5 | 1.8e-5 | 1.2e-5

The robot walks with a periodic crawl gait, moving first the hind legs and
then the front legs, with a fixed stride of 8 cm and a step height of 4 cm.
For this study, we decouple the onward motion of the base and the swing leg so
the RAMP can distribute the swing momentum to the whole body, including the
base link. Therefore, between swing movements, we have a phase where the base
moves 2 cm from its previous neutral position, adjusting the position and
attitude to complete a periodic cycle after the motion of the four legs. We
also implemented release and grasping vertical movements of 1 cm each to avoid
ground collisions with the swing trajectory. The total swing period of each
leg is 1.5 s, the same as the period of each base movement, adding to 12 s for
one cycle period of this quadruped robot.

Figure 2: Snapshots of the simulation study for rough surface locomotion with
a quadruped robot. Blue lines represent the swing trajectory, and pink is the
base trajectory. The baseline (BL) case fails due to detachment, inducing
flotation on the robot, while the proposed method RAMP with Partial Momentum
Distribution (RAMP-PMD) successfully reaches the goal position. RAMP with Full
Momentum Distribution (RAMP-FMD) fails due to a singularity.

Fig. 2 shows snapshots of three different cases of simulations. The first is a
baseline (BL) case without any RAMP element, simply computing swing
trajectories from the current to the next grasping position through a via
point at the desired height. The second case performs the Reaction-Aware
Motion Planning with a Partial Momentum Distribution (RAMP-PMD) with a
distribution factor of 0.5. And the last case is also a RAMP, but with Full
Momentum Distribution (RAMP-FMD), i.e., the distribution factor is equal to 1.
After RAMP computes the desired trajectories for the robot, an inverse
kinematics algorithm computes the angular positions of each joint, which are
the inputs of a PD torque controller to execute the motion.

The results show that after a few steps, walking with the regular baseline
motion plan leads to ground reaction forces higher than the gripper holding
limit, causing detachment and flotation of the robot. This failure case is
highly undesired in real scenarios, as this could lead to a complete loss of
the robot during the exploration mission. By planning the motion with RAMP-
PMD, the robot successfully advances 40 cm without exceeding force limits.
However, in the last case with RAMP-FMD, the robot moves its base excessively
to distribute the entire swing momentum, falling into a singularity
configuration.

The simulation results show that RAMP can effectively increase the safety of
legged mobility in microgravity, transforming a failure case into a successful
walk until the desired location. However, the distribution factor has to be
tuned appropriately to avoid a swing momentum distribution beyond the robot’s
capabilities, considering parameters such as inertia, moving speed, and
manipulability.

##  IV Experimental Verification

Although it is difficult to emulate microgravity conditions accurately on
Earth, it is possible to create a planar environment by eliminating the
friction between a robot and an extremely flat surface. A common technique is
to levitate the robot with pressurized air, avoiding direct contact with the
flat surface. We use such system to perform two-dimensional verification of
the RAMP method to improve mobility in microgravity.

Figure 3: Dual-Arm Air-Floating Robot. The robot floats with pressurized air
from the air-bearings, eliminating friction to create a two-dimensional
microgravity environment. Each limb has 3 degrees of freedom controllable by
DC motors, and a force/torque sensor at its tip. One of the limbs is directly
attached to a fixed plate to emulate the grasping of a surface.

###  IV-A Planar Robot with Fixed Leg

TABLE II: Inertial Parameters of Dual-Arm Air-Floating Robot Link | Size (mm) | Mass (g) | Inertia (kgm2)
---|---|---|---
|  |  | $I_{zz}$
Base | 160 x 295 | 8524.0 | 1.2e-1
1 | 25 | 608.0 | 1.6e-3
2 | 17.5 | 626.2 | 7.1e-4
3 | 87.25 | 248.1 | 2.0e-4

For this experimental study, we use the Dual-Arm Air-Floating Robot shown in
Fig. 3, with size and inertial parameters detailed in Table II. The robot has
three Degrees of Freedom per limb, with a force and torque (F/T) sensor at its
tip. We attach one of the limbs to a fixed plate to emulate grasping without
detachment to measure ground reactions. Since the robot cannot walk
continuously, we only perform a single-step analysis. The step has a 15 cm
stride and 7 cm height, with a 2 s swing period. We compute the trajectories
of the base and the limbs, and the inverse kinematics offline, while the robot
applies a velocity control based on the desired joint positions.

Figure 4: Snapshots of all experiments cases using the Dual-Arm Air-Floating
Robot. The left limb is fixed directly to a plate, and the right leg moves
along the planned trajectory (blue line) during 2 s. The baseline (BL) case
moves the leg with a standard spline trajectory, while the Low-Reaction Swing
Trajectory (LRST) case moves the swing limb in a trajectory that minimizes the
momentum variation. The partial Momentum Distribution (PMD) case also moves
the base and support limb to distribute half of the swing momentum to the
whole body. The Full Momentum Distribution (FMD) case distributes 100% of the
swing momentum.

Fig. 4 shows snapshots of four distinct experiment cases. The first is the
Baseline (BL) case with a standard trajectory and no momentum distribution,
while the other three use RAMP with different momentum distribution factors.
The second case is the Low-Reaction Swing Trajectory (LRST), which has only
the optimized trajectory but no momentum distribution. The third and fourth
cases are Partial Momentum Distribution (PMD) and Full Momentum Distribution
(FMD), with momentum distribution factor of 0.5 and 1, respectively.

We can observe the differences between low-reaction and standard trajectories
by comparing the four cases. The LRST minimizes the motion reactions by
reducing the limb’s inertia when bringing it closer to the main body.
Moreover, we notice how the robot moves its whole body to distribute the swing
momentum in the last two cases, showing an attitude change of $11^{\circ}$ and
$21^{\circ}$ for PMD and FMD, respectively. As this robot has higher inertia
and manipulability when compared to the robot used in the simulation study, we
do not observe any problem related to singular configurations, even with the
robot moving at high speeds.

Figure 5: Time-history graphs of Ground Reaction Force and Moment of the
experimental study. Baseline (BL) is in red, Low-Reaction Swing Trajectory
(LRST) is in yellow, Partial Momentum Distribution (PDM) is in blue, and Full
Momentum Distribution (FMD) is in purple. We verify that by increasing the
degree of the Reaction-Aware Motion Planning (RAMP) with low-reaction swing
trajectory and momentum distribution, the reaction forces and moments acting
on the support limb are mitigated.

Figure 6: Summary of maximum and average ground reactions in relation to the
Baseline case. Both forces and moments are reduced by employing RAMP, with a
progressive decrease of maximum and average values as we increase the momentum
distribution factor. With FMD, it is possible to decrease up to 60% of the
maximum ground reactions and up to 75% of average ground reactions.

Fig. 5 and Fig. 6 show the data obtained from the F/T sensor of the fixed
limb. Compared to the baseline case shown in red, RAMP-FMD can reduce the
maximum ground reactions to 40% of the force and torque observed in the
Baseline case. RAMP also considerably decreases average force and moment,
ranging from 60% to 25% of the baseline average, depending on the momentum
distribution factor.

RAMP demonstrated its efficacy in decreasing ground reaction forces and
moments with the experimental results of a planar robot in an emulated
microgravity. From these results, we can expect that RAMP can effectively
increase the locomotion safety of multi-legged robots in a real scenario by
reducing the detachment risk of the gripper due to high pulling forces.

##  V Conclusions

We propose a novel Reaction-Aware Motion Planning (RAMP) for the mobility of
multi-legged robots in microgravity, computing the optimal trajectories for
the swing legs and base that mitigate motion reactions. The main contribution
of our approach is improving the safety of dynamic multi-legged locomotion in
microgravity environments by decreasing the risk of surface separation. RAMP
also provides low-reaction alternatives to reactionless control strategies
with an adjustable momentum distribution factor to avoid singular
configurations. We demonstrated the efficacy of our approach in preventing
surface detachment through numerical simulations and experiments, indicating
that real robotic mission scenarios can employ RAMP for safe mobility in
microgravity.

## References

  * [1] V. Badescu, _Asteroids: Prospective Energy and Material Resources_ , V. Badescu, Ed. Heidelberg: Springer Berlin, 2013.
  * [2] J. F. Russell, D. M. Klaus, and T. J. Mosher, “Applying analysis of international space station crew-time utilization to mission design,” _Journal of spacecraft and rockets_ , vol. 43, no. 1, pp. 130–136, 2006.
  * [3] T. Yoshimitsu, I. Nakatani, and T. Kubota, “New mobility system for small planetary body exploration,” in _Proceedings 1999 IEEE International Conference on Robotics and Automation (Cat. No. 99CH36288C)_ , vol. 2. IEEE, 1999, pp. 1404–1409.
  * [4] T. Yoshimitsu and T. Kubota, “Asteroid surface exploration by Minerva-II small rovers,” in _18th Annual Meeting of the Asia Oceania Geosciences Society: Proceedings of the 18th Annual Meeting of the Asia Oceania Geosciences Society (AOGS 2021)_. World Scientific, 2022, pp. 180–182.
  * [5] T.-M. Ho, V. Baturkin, C. Grimm, J. T. Grundmann, C. Hobbie, E. Ksenik, C. Lange, K. Sasaki, M. Schlotterer, M. Talapina, _et al._ , “MASCOT—the mobile asteroid surface scout onboard the Hayabusa2 mission,” _Space Science Reviews_ , vol. 208, no. 1, pp. 339–374, 2017.
  * [6] Y. Tsuda, T. Saiki, F. Terui, S. Nakazawa, M. Yoshikawa, S.-i. Watanabe, and H. P. Team, “Hayabusa2 mission status: Landing, roving and cratering on asteroid Ryugu,” _Acta Astronautica_ , vol. 171, pp. 42–54, 2020.
  * [7] K. Nagaoka, K. Yoshida, M. Kurisu, K. Osuka, K. Tadakuma, Y. Tsumaki, T. Mineta, S. Kimura, T. Narumi, T. Kubota, _et al._ , “Development of MINERVA-II2, a micro-robot for asteroid surface exploration with innovative mobility,” in _The 11th Low-Cost Planetary Missions Conference_ , 2015.
  * [8] B. J. Hockman, A. Frick, R. G. Reid, I. A. Nesnas, and M. Pavone, “Design, control, and experimentation of internally-actuated rovers for the exploration of low-gravity planetary bodies,” _Journal of Field Robotics_ , vol. 34, no. 1, pp. 5–24, 2017.
  * [9] K. Kobashi, A. Bando, K. Nagaoka, and K. Yoshida, “Tumbling and hopping locomotion control for a minor body exploration robot,” in _2020 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_. IEEE, 2020, pp. 1871–1878.
  * [10] M. G. Bualat, T. Smith, E. E. Smith, T. Fong, and D. Wheeler, “Astrobee: A new tool for ISS operations,” in _2018 SpaceOps Conference_ , 2018, p. 2517\.
  * [11] V. Schröder, R. Regele, J. Sommer, T. Eisenberg, and C. Karrasch, “GNC system design for the crew interactive mobile companion (CIMON),” in _Proceedings of IAF Human Spaceflight Symposium at 69th International Astronautical Congress_ , 2018.
  * [12] S. Mitani, M. Goto, R. Konomura, Y. Shoji, K. Hagiwara, S. Shigeto, and N. Tanishima, “Int-ball: Crew-supportive autonomous mobile camera robot on ISS/JEM,” in _2019 IEEE Aerospace Conference_. IEEE, 2019, pp. 1–15.
  * [13] K. Yoshida, T. Maruki, and H. Yano, “A novel strategy for asteroid exploration with a surface robot,” in _34th COSPAR Scientific Assembly_ , vol. 34, 2002, p. 1966.
  * [14] M. Diftler, T. Ahlstrom, R. Ambrose, N. Radford, C. Joyce, N. De La Pena, A. Parsons, and A. Noblitt, “Robonaut 2 — initial activities on-board the ISS,” in _2012 IEEE Aerospace Conference_. IEEE, 2012, pp. 1–12.
  * [15] A. Parness, N. Abcouwer, C. Fuller, N. Wiltsie, J. Nash, and B. Kennedy, “LEMUR 3: A limbed climbing robot for extreme terrain mobility in space,” in _2017 IEEE international conference on robotics and automation (ICRA)_. IEEE, 2017, pp. 5467–5473.
  * [16] Q. Wen, J. He, and F. Gao, “Stochastic models and control of anchoring mechanisms for grasping in microgravity,” _Applied Sciences_ , vol. 12, no. 6, p. 3196, 2022.
  * [17] S. Schneider, A. Bylard, T. G. Chen, P. Wang, M. Cutkosky, M. Lapôtre, and M. Pavone, “Reachbot: A small robot for large mobile manipulation tasks,” in _2022 IEEE Aerospace Conference (AERO)_. IEEE, 2022, pp. 1–12.
  * [18] M. Chacin, A. Mora, and K. Yoshida, “Motion control of multi-limbed robots for asteroid exploration missions,” in _2009 IEEE International Conference on Robotics and Automation_. IEEE, 2009, pp. 3037–3042.
  * [19] Y. Yuguchi, K. Nagaoka, and K. Yoshida, “Verification of gait control based on reaction null-space for ground-gripping robot in microgravity,” in _2016 IEEE International Conference on Robotics and Automation (ICRA)_. IEEE, 2016, pp. 2822–2827.
  * [20] W. F. R. Ribeiro, K. Uno, and K. Yoshida, “Low-reaction trajectory generation for a legged robot in microgravity,” in _2022 IEEE/SICE International Symposium on System Integration (SII)_. IEEE, 2022, pp. 505–510.
  * [21] K. Uno, W. F. R. Ribeiro, Y. Koizumi, K. Haji, K. Kurihara, W. Jones, and K. Yoshida, “ClimbLab: MATLAB simulation platform for legged climbing robotics,” in _Climbing and Walking Robots Conference_. Springer, 2021, pp. 229–241.
  * [22] K. Yoshida, “The SpaceDyn: a MATLAB toolbox for space and mobile robots,” in _Proceedings 1999 IEEE/RSJ International Conference on Intelligent Robots and Systems. Human and Environment Friendly Robots with High Intelligence and Emotional Quotients (Cat. No. 99CH36289)_ , vol. 3. IEEE, 1999, pp. 1633–1638.
  * [23] K. Uno, N. Takada, T. Okawara, K. Haji, A. Candalot, W. F. R. Ribeiro, K. Nagaoka, and K. Yoshida, “Hubrobo: a lightweight multi-limbed climbing robot for exploration in challenging terrain,” in _2020 IEEE-RAS 20th International Conference on Humanoid Robots (Humanoids)_. IEEE, 2021, pp. 209–215.


# Solving break minimization problems
in mirrored double round-robin tournament
with QUBO solver ††thanks:  This work was supported by JSPS KAKENHI Grant
Number JP20K04973.

Koichi Fujii Department of Industrial Engineering and Economics, Tokyo
Institute of Technology  Tomomi Matsui 22footnotemark: 2

###### Abstract

The break minimization problem is a fundamental problem in sports scheduling.
Recently, its quadratic unconstrained binary optimization (QUBO) formulation
has been proposed, which has gained much interest with the rapidly growing
field of quantum computing. In this paper, we demonstrate that the state-of-
the-art QUBO solver outperforms the general mixed integer quadratic
programming (MIQP) solver on break minimization problems in a mirrored double
round-robin tournament. Moreover, we demonstrate that it still outperforms or
is competitive even if we add practical constraints, such as consecutive
constraints, to the break minimization problem.

Keywords: sports scheduling, round robin tournament, QUBO

##  1 Introduction

The construction of a suitable schedule for sports competitions is a crucial
issue in sports scheduling. The round-robin tournament is a competition
adopted in many professional sports, such as soccer and basketball. This paper
deals with a home-away assignment problem, which assigns home or away to each
match of a double round-robin tournament schedule to minimize the number of
breaks.

In recent years, Urdaneta, Yuan, and Sinqueria [1] formulated the break
minimization problem in a double round-robin tournament as a mixed integer
quadratic programming (MIQP) with no constraint or a quadratic unconstrained
binary optimization problem (QUBO). They demonstrated the advantage of QUBO
formulation against other formulations through numerical experiments.
Kuramata, Katsuki, and Nakata [2] successfully applied quantum annealing to
the QUBO formulation, obtaining better solutions quickly.

QUBO has recently gained much attention owing to its relationship with quantum
computing. One of the best QUBO solvers was developed by Rehfeldt, Koch and
Shinano [3]. They showed that their solver generally outperforms quantum
techniques, including quantum annealing, although it is designed to prove
optimality, contrary to heuristic approaches of quantum techniques.

In this paper, we demonstrate that the QUBO solver outperforms general MIQP
solver for break minimization problems in a mirrored double round-robin
tournament schedule. We also extend the break minimization problem by adding
consecutive constraints, which prohibit more than two consecutive home or away
matches. Formulating these extended problems as a QUBO model through
penalization, we demonstrate that the QUBO solver can also outperform the MIQP
solver.

The remainder of the paper is organized as follows. Section 2 describes the
break minimization problem in a double round-robin tournament schedule.
Section 3 gives a mathematical programming formulation of the break
minimization problem. We also show how to insert consecutive constraints into
our formulations of the break minimization problem. Section 4 provides
computational results comparing the general MIQP and QUBO solvers for the
break minimization problem. It is shown that the QUBO solver outperforms or is
competitive with the general MIQP solver, not only for the break minimization
problem but also for problems with consecutive constraints.

##  2 Break Minimization

In this section, we introduce round-robin tournament schedules, which
frequently arise in many sports scheduling problems.

A single round-robin tournament (SRRT) is a simple sports competition, where
each team plays against every other team once. A double round-robin tournament
(DRRT) is an extended version of an SRRT, where each team plays against all
others twice. In this paper, we consider a DRRT schedule with the following
properties:

  * •

Each team plays one match in each slot, the day when a match is held.

  * •

Each team has a home stadium, and each match is held at the home stadium of
one of the two playing teams. When a team plays a match at the home of the
opponent, we say that the team plays away.

  * •

Each team plays each other team exactly once at its home stadium and once
away.

A mirrored DRRT (MDRRT) is a DRRT with the same combination of matches in the
first and second halves.

###  2.1 Preliminaries

We use the following notations and symbols throughout this paper:

  * •

$2n$: number of teams, where $n\geq 2$,

  * •

$T=\\{1,2,\ldots,2n\\}$: set of teams,

  * •

$S=\\{1,2,\ldots,2(2n-1)\\}$: set of slots.

A schedule of a DRRT is described as a pair of a timetable and a home–away
assignment defined below. A timetable $\tau$ of a DRRT is a matrix whose rows
and columns are indexed by $T$ and $S$, respectively. An element $\tau(t,s)$
shows the opponent that plays against team $t$ at slot $s$. A timetable $\tau$
(of a DRRT schedule) should meet the following conditions: (i) a row of $\tau$
indexed by team $t\in T$ is obtained by concatenating two permutations of
$T\setminus\\{t\\}$, and (ii) $\tau(\tau(t,s),s)=t\;(\forall(t,s)\in T\times
S).$ In the case of an MDRRT, the corresponding timetable satisfies that each
row of $\tau$ indexed by $t$ is a repetition of a mutual permutation; i.e.,
$\tau(t,s)=\tau(t,s+2n-1)$ holds for each $(t,s)\in
T\times\\{1,\ldots,2n-1\\}$. Table 2 shows a timetable for an MDRRT schedule.

Table 1: Timetable Slot | $\mathbf{1}$ | $\mathbf{2}$ | $\mathbf{3}$ | $\mathbf{4}$ | $\mathbf{5}$ | $\mathbf{6}$
---|---|---|---|---|---|---
team 1 | 2 | 3 | 4 | 2 | 3 | 4
team 2 | 1 | 4 | 3 | 1 | 4 | 3
team 3 | 4 | 1 | 2 | 4 | 1 | 2
team 4 | 3 | 2 | 1 | 3 | 2 | 1

Table 2: Howe-away assignment Slot | $\mathbf{1}$ | $\mathbf{2}$ | $\mathbf{3}$ | $\mathbf{4}$ | $\mathbf{5}$ | $\mathbf{6}$
---|---|---|---|---|---|---
team 1 | 1 | 0 | 1 | 0 | 1 | 0
team 2 | 0 | 0 | 1 | 1 | 1 | 0
team 3 | 1 | 1 | 0 | 0 | 0 | 1
team 4 | 0 | 1 | 0 | 1 | 0 | 1

Given a timetable $\tau$ of a DRRT, we introduce a set defined by

$\Xi(\tau)=\\{(t_{1},t_{2},s_{1},s_{2})\in T^{2}\times S^{2}\mid
t_{1}<t_{2},s_{1}<s_{2},t_{2}=\tau(t_{1},s_{1})=\tau(t_{1},s_{2})\\}.$

We denote $\Xi(\tau)$ by $\Xi$ if there is no ambiguity. Each element
$(t_{1},t_{2},s_{1},s_{2})\in\Xi$ satisfies the following properties:

  * •

$t_{1},t_{2}$: a pair of teams with $t_{1}<t_{2}$.

  * •

$s_{1}$: the slot of the first match between $t_{1}$ and $t_{2}$.

  * •

$s_{2}$: the slot of the second match between $t_{1}$ and $t_{2}$.

Clearly, the cardinality of the set $\Xi$ is equal to $n(2n-1)$. For an MDRRT
schedule, $s_{2}=s_{1}+2n-1$ holds for each $(t_{1},t_{2},s_{1},s_{2})\in\Xi$.
The set $\Xi$ corresponding to the timetable in Table 2 becomes

$\Xi=\\{(1,2,1,4),(1,3,2,5),(1,4,3,6),(2,3,3,6),(2,4,2,5),(3,4,1,4)\\}.$

A team is said to be at home in slot $s$ if the team plays a match at its home
stadium in $s$; otherwise away in $s$. A home-away assignment (HA-assignment)
is a matrix $Y$ whose rows are indexed by $T$, and columns are indexed by $S$.
An element $y(t,s)$ of $Y$ is equal to 1 if team $t$ plays a match in $s$ at
its home stadium; otherwise it is 0. Tables 2 and 2 together show an example
of an MDRRT schedule. For a given timetable $\tau$, an HA-assignment $Y$ is
consistent with $\tau$ when the following holds:

$\forall(t_{1},t_{2},s_{1},s_{2})\in\Xi(\tau),\left(\begin{array}[]{cc}y(t_{1},s_{1})&y(t_{1},s_{2})\\\
y(t_{2},s_{1})&y(t_{2},s_{2})\end{array}\right)=\left(\begin{array}[]{cc}1&0\\\
0&1\end{array}\right)\mbox{ or }\left(\begin{array}[]{cc}0&1\\\
1&0\end{array}\right).$ (1)

Given an HA-assignment $Y,$ we say that team $t$ has a break at slot $s\in
S\setminus\\{1\\}$ if $y(t,s-1)=y(t,s).$ The number of breaks in a home–away
assignment is defined as the total number of breaks belonging to all teams.
For example, in Table 2, team $2$ has a break at slot $2$, as there are
consecutive away matches $(y(2,1)=y(2,2)=0)$. Similarly, team $2$ has a break
in slot $4$ because of the consecutive home matches $(y(2,3)=y(2,4)=1)$. Given
a timetable of a round-robin tournament, we need to find a consistent HA-
assignment to complete a schedule. In practical sports timetabling such as
[4], a break is considered undesirable.

###  2.2 Previous Studies

Most methods for constructing a tournament schedule use one of the following
two decomposition approaches [5]: (i) first-break-then-schedule and (ii)
first-schedule-then-break. In the first-break-then-schedule approach, a matrix
of HA-assignment is generated at the first stage. The second stage finds a
timetable consistent with the the fixed HA-assignment, if it exists. In the
first-schedule-then-break approach, in contrast, a timetable is constructed at
the first stage. The second stage determines the home team and the away team
of each match. In this approach, the home advantage is further determined in
the second stage. In this paper, we address a problem of the second stage in
the first-schedule-then-break approach:

Break Minimization Problem: Given a timetable $\tau$, the break minimization
problem finds an HA-assignment consistent to $\tau$ that minimizes the number
of breaks.

The break minimization problem has been widely studied [6, 1, 7]. In
particular, De Werra [8] proved that the number of breaks in an SRRT schedule
is more than $2n-2$, and the number of breaks in an MDRRT schedule is more
than $6n-6$. Post and Woeginger [9] proved that every timetable of an SRRT has
an HA-assignment such that the number of breaks is at most $n(n-1)$ if $n=2k$
and $(n-1)^{2}$ if $n=2k+1$. Though the break minimization problems in SRRTs
or MDRRTs are known to be hard to solve, they have not been proved as NP-hard
to the best of our knowledge.

Trick [6] proposed an integer linear programming model for the break
minimization problem in an SRRT. The model includes highly symmetry, which
makes the problem difficult to solve using an integer programming approach.
Régin [10] proposed a constraint programming model for the break minimization
problem in an SRRT. Elf, Jünger and Rinaldi [11] showed that solving the break
minimization problem in an SRRT can be transformed into a maximum cut problem
(MAX CUT). Miyashiro and Matsui [7] formulated the break minimization in an
SRRT as a special case of MAX RES CUT discussed by Goemans and Williamson
[12]; thus, it is solvable by an algorithm based on positive semidefinite
programming. Suzuka, et al. [13] addressed the HA-assignment problem and break
minimization/maximization problem and gave a unified view to the three
problems. Their computational results showed that the break minimization
problem in an MDRRT schedule is more difficult than that in an SRRT schedule.
More precisely, if we solve the problems by an integer programming solver, the
average computational time for the break minimization problem in an MDRRT is
apparently longer than that in an SRRT. Urdaneta, Yuan and Sinqueria [1]
demonstrated a QUBO formulation of the break minimization problem in a DRRT.
They demonstrated that it is superior to the mixed integer linear programming
formulation [6] for a general mathematical optimization solver. The QUBO
formulation proposed in [1] can also be regarded as an extension of the MAX
CUT formulation proposed in [11]. Kuramata, Katsuki and Nakata [2]
successfully applied quantum annealing to a QUBO formulation. They used the
quantum annealer D-wave Advantage and demonstrated that it could obtain better
solutions quickly compared to a state-of-the-art general solver.

The schedules of most major sports leagues contain restrictions on the number
of consecutive home and away matches. The most common restriction used is an
upper bound of consecutive home or away matches [14]. In practice, the
consecutive constraints appear as tournament problems of the Danish soccer
league [15]. In this paper, we also address a break minimization problem with
consecutive constraints, which prohibits consecutive home or away matches
exceeding two and/or three.

##  3 Mathematical Formulation

In this section, we introduce the QUBO formulation of the break minimization
problem studied by Urdaneta et al. [1].

For each pair $(t,s)\in T\times S$, we introduce a binary variable
$y_{t,s}\in\\{0,1\\}$, which takes 1 if team $t$ plays at home in slot $s$ and
0 if it plays away. A 0-1 matrix $Y=(y_{t,s})$ is consistent with a given
timetable $\tau$ if and only if $Y$ satisfies the following equations:

$\displaystyle
y_{t_{1},s_{1}}+y_{t_{2},s_{1}}=1\;\;(\forall(t_{1},t_{2},s_{1},s_{2})\in\Xi(\tau)),$
$\displaystyle
y_{t_{1},s_{1}}+y_{t_{1},s_{2}}=1\;\;(\forall(t_{1},t_{2},s_{1},s_{2})\in\Xi(\tau)),$
(2) $\displaystyle-
y_{t_{1},s_{1}}+y_{t_{2},s_{2}}=0\;\;(\forall(t_{1},t_{2},s_{1},s_{2})\in\Xi(\tau)).$

The above constraints are essentially equivalent to the condition in (1).
These variables and constraints enable the break minimization problem to be
formulated as a constrained quadratic binary optimization problem defined as

$\displaystyle\text{CQ}:$ minimize $\displaystyle\quad\sum_{t\in T}\sum_{s\in
S\backslash\\{4n-2\\}}\Bigl{(}y_{t,s}y_{t,s+1}+\left(1-y_{t,s}\right)\left(1-y_{t,s+1}\right)\Bigr{)}$
subject to $\displaystyle\quad
y_{t_{1},s_{1}}+y_{t_{2},s_{1}}=1\quad(\forall(t_{1},t_{2},s_{1},s_{2})\in\Xi(\tau)),$
(3) $\displaystyle\quad
y_{t_{1},s_{1}}+y_{t_{1},s_{2}}=1\quad(\forall(t_{1},t_{2},s_{1},s_{2})\in\Xi(\tau)),$
(4) $\displaystyle\quad-
y_{t_{1},s_{1}}+y_{t_{2},s_{2}}=0\quad(\forall(t_{1},t_{2},s_{1},s_{2})\in\Xi(\tau)),$
(5) $\displaystyle\quad y_{t,s}\in\\{0,1\\}\quad(\forall(t,s)\in T\times S).$

The objective function of CQ consists of the total sums of $y_{t,s}y_{t,s+1}$
and $(1-y_{t,s})(1-y_{t,s+1})$, which represent the number of breaks
corresponding to the two consecutive home and away matches, respectively.

Urdaneta et al. [1] demonstrated that the constrained optimization problem CQ,
which represents a break minimization problem in a DRRT, can be expressed as
an unconstrained quadratic binary optimization problem. For each
$\xi=(t_{1},t_{2},s_{1},s_{2})\in\Xi$, we introduce a variable $z_{\xi}$ and
eliminate four variables
$y_{t_{1},s_{1}},y_{t_{2},s_{1}},y_{t_{1},s_{2}},y_{t_{2},s_{2}}$ by the
following equalities:

$\displaystyle y_{t_{1},s_{1}}=z_{\xi},\quad y_{t_{2},s_{1}}=1-z_{\xi},\quad
y_{t_{1},s_{2}}=1-z_{\xi},\quad y_{t_{2},s_{2}}=z_{\xi}.$ (6)

By substituting equations (6) into CQ, constraint (3) in CQ becomes redundant,
as $y_{t_{1},s_{1}}+y_{t_{2},s_{1}}=z_{\xi}+1-z_{\xi}=1.$ In a similar manner,
we can show the redundancy of constraints (4) and (5). By deleting theses
redundant constraints, we obtain an unconstrained optimization problem
$\\{\mbox{\boldmath$z$}^{\top}H\mbox{\boldmath$z$}\mid\mbox{\boldmath$z$}\in\\{0,1\\}^{\Xi}\\}$
for the break minimization problem.

Next, we introduce additional consecutive constraints. A sequence of
consecutive away matches of a team is called a road trip, and a sequence of
consecutive home matches is called a home stand. The consecutive constraints
for a positive integer $u\geq 2$ is defined as follows:

CC($u$):

the length of any home stand is at most $u$, and the length of any road trip
is at most $u$.

In this paper, we deal with cases that $u\in\\{2,3\\}$. For example, in Table
2, team 2 plays at away from slot 3 to slot 5, which indicates a violation of
the constraint CC(2). The constraint CC(3) has been employed in many studies
[14, 16, 17, 18]. The constraint CC(2) is discussed in [19, 20, 21].

It is well-known that every single round robin timetable has a consistent HA-
assignment satisfying the constraint CC($u$) for any $u\geq 2$. For example,
the proof of Remark 7 in [22] gives a procedure to construct an HA-assignment
(of an SRRT) in which each team has no breaks at slots $2,4,\ldots,2n-2$. This
procedure directly shows that every timetable of MDRRT has an HA-assignment
satisfying CC(3). The hardness of the problem of finding an HA-assignment
satisfying CC(2) remains unknown in case of an MDRRT schedule.

We obtain a formulation of the break minimization problem in MDRRT with the
constraint CC($u$) by adding the following linear inequalities to CQ:

$1\leq y_{t,s}+y_{t,s+1}+\cdots+y_{t,s+u}\leq u\;\;(\forall(t,s)\in
T\times\\{1,2,\ldots,2n-1\\}).$ (7)

The first inequality indicates that the number of home matches should be
positive, which is equivalent to having less than or equal to $u$ consecutive
away matches. We need not consider cases that $(t,s)\in
T\times\\{2n,2n+1,\ldots\\}$, as we discuss an MDRRT schedule. In the case of
CC(2), the constraint (7) is equivalently expressed in quadratic equality
constraints as

$\displaystyle(y_{t,s}+y_{t,s+1}+y_{t,s+2}-1)(y_{t,s}+y_{t,s+1}+y_{t,s+2}-2)=0\;\;\;$
(8) $\displaystyle(\forall(t,s)\in T\times\\{1,2,\ldots,2n-1\\}).$

When we address a problem with CC(3), we introduce artificial binary variables
$w_{t,s}\in\\{0,1\\}$ and employ the following quadratic equality constraints

$\displaystyle(y_{t,s}+y_{t,s+1}+y_{t,s+2}+y_{t,s+3}+w_{t,s}-2)$ (9)
$\displaystyle\;\;\;(y_{t,s}+y_{t,s+1}+y_{t,s+2}+y_{t,s+3}+w_{t,s}-3)=0$
$\displaystyle\;\;\;\;\;\;\;\mbox{\qquad}(\forall(t,s)\in
T\times\\{1,2,\ldots,2n-1\\}).$

We employ a penalty method to solve the break minimization problem with
constraints CC(2) and/or CC(3). We construct a QUBO formulation by adding a
penalty function (to the objective function) that consists of penalty
parameters multiplied by the left-hand side of (8) or (9). It is easy to see
that the left-hand sides of both (8) and (9) are nonnegative when the
corresponding variables are 0-1 valued. Thus, by introducing sufficiently
large penalty parameters, every optimal solution of the QUBO formulation
obtained by the above penalty method is also optimal to the break minimization
problem with constraints CC(2) and/or CC(3), when it is feasible. Fred Glover
et al. [23] discussed a general technique to incorporate equality constrains
in a QUBO formulation.

##  4 Computational Results

This section presents computational results on break minimization problems in
MDRRTs with consecutive constraints. We report computational results obtained
by a general MIQP solver, Gurobi Optimizer (version 9.5.2), and a QUBO solver,
QuBowl [3]. The computational experiments were performed on Intel (R) Xeon(R)
CPU E3-1270 with 32 GB RAM.

QuBowl is a state-of-the-art QUBO solver developed by Rehfeldt, Koch and
Shinano [3]. It converts QUBO problems to MaxCut problems and solve them by
conducting a branch-and-cut method. Enhancing the algorithm with presolving
and cutting-plane generation, it outperforms state-of-the-art commercial
solvers on general QUBO instances.

### Testset

We generate break minimization problem instances by the following method from
[2].

  1. Step 1:

Create an RRT timetable using the Kirkman method [24].

  2. Step 2:

Shuffle the order of the slots in the timetable created in Step 1.

  3. Step 3:

Make a copy of the timetable obtained in Step 2 and concatenate two identical
timetables.

### Computational comparison of MIQP and QUBO

In our experiments, we generated and solved 5 instances for each
$n\in\\{3,4,\ldots,15\\}$ using Gurobi and QuBowl.

In the first experiment, we solved (unconstrained) break minimization
problems. Table 3 provides the results of our first experiment for a time-
limit of one hour using a single thread. The 2nd and 3rd columns show the
numbers of instances completely solved (within one hour) by Gurobi and QuBowl,
respectively. The 5th and 6th columns show the average execution times of
Gurobi and QuBowl, where time is set as the time limit when the solvers failed
to solve.

Both Gurobi and QuBowl took longer to compute as the number $n$ increased.
Table 3 clearly shows that QuBowl outperforms Gurobi. Gurobi failed to solve 3
instances (within one hour) with 28 teams $(n=14)$ and 30 teams $(n=15)$,
while QuBowl successfully solved all instances. For 26 teams $(n=13)$, QuBowl
is faster than Gurobi by a factor of more than 40 on average.

In the second experiment, we solved the break minimization problems with
consecutive constraints CC(2) and CC(3) using Gurobi and QuBowl. When we use
Gurobi, we solve nonconvex quadratic integer problems with linear constraints
(3). If we use QuBowl, we penalize quadratic constraints (8) and/or (9) by
setting penalty parameters to 10 to translate them into QUBOs. Tables 4 and 5
list the results of the second experiment with a time-limit of one hour. It
can be observed that QuBowl requires more computational time than the first
experiment. We verified that the solutions obtained by QuBowl satisfy the
constraints CC(2) and/or CC(3). Table 4 shows that for 26 teams $(n=13)$,
QuBowl is faster than Gurobi by a factor of more than 15 on average. Though
the running time gaps between QuBowl and Gurobi get close compared with the
first experiment, QuBowl still outperforms Gurobi for those problems. In case
of CC(3), Table 5 shows that QuBowl is competitive with Gurobi.

Table 3: MIQP solver vs. QUBO solver for break minimization problems | Solved | run time [s]
---|---|---
$n$ | Grb | QuBowl | Grb | QuBowl
3 | 5 | 5 | 0.01 | 0.00
4 | 5 | 5 | 0.03 | 0.01
5 | 5 | 5 | 0.06 | 0.02
6 | 5 | 5 | 0.14 | 0.04
7 | 5 | 5 | 0.57 | 0.10
8 | 5 | 5 | 1.09 | 0.14
9 | 5 | 5 | 2.26 | 0.31
10 | 5 | 5 | 4.98 | 0.76
11 | 5 | 5 | 36.48 | 2.20
12 | 5 | 5 | 155.71 | 4.53
13 | 5 | 5 | 492.31 | 10.43
14 | 2 | 5 | 2734.31 | 752.48
15 | 2 | 5 | 2677.92 | 155.88
Table 4: MIQP solver vs. QUBO solver for break minimization problems with CC(2). | Solved | run time [s]
---|---|---
$n$ | Grb | QuBowl | Grb | QuBowl
3 | 5 | 5 | 0.02 | 0.00
4 | 5 | 5 | 0.05 | 0.03
5 | 5 | 5 | 0.08 | 0.05
6 | 5 | 5 | 0.18 | 0.11
7 | 5 | 5 | 0.59 | 0.26
8 | 5 | 5 | 0.93 | 0.28
9 | 5 | 5 | 7.51 | 0.79
10 | 5 | 5 | 20.18 | 4.89
11 | 5 | 5 | 128.43 | 9.33
12 | 5 | 5 | 288.53 | 23.57
13 | 5 | 5 | 687.90 | 61.73
14 | 2 | 5 | 2797.60 | 688.87
15 | 2 | 4 | 2839.15 | 1002.10
Table 5: MIQP solver vs. QUBO solver for break minimization problems with CC(3). | Solved | run time [s]
---|---|---
$n$ | Grb | QuBowl | Grb | QuBowl
3 | 5 | 5 | 0.02 | 0.26
4 | 5 | 5 | 0.04 | 0.39
5 | 5 | 5 | 0.08 | 0.56
6 | 5 | 5 | 0.15 | 0.49
7 | 5 | 5 | 0.27 | 1.35
8 | 5 | 5 | 0.73 | 1.80
9 | 5 | 5 | 3.01 | 5.54
10 | 5 | 5 | 7.87 | 12.32
11 | 5 | 5 | 54.09 | 50.21
12 | 5 | 5 | 165.97 | 107.89
13 | 5 | 5 | 543.06 | 338.84
14 | 3 | 2 | 2523.42 | 2323.49
15 | 3 | 2 | 2403.28 | 2866.82

##  5 Conclusion

This study demonstrated that the QUBO solver outperforms a general MIQP solver
for break minimization problems. This indicates an important role of the QUBO
solver in sport scheduling problems. Further, the QUBO solver can also
outperform in break minimization problems by transforming them into QUBOs.
Although penalizing the quadratic constraints makes the objective function
more dense, the enhancing techniques in the QUBO solver enables more efficient
solving than a general MIQP solver.

Our results indicate a further potential of the QUBO solver and QUBO
formulations. In recent years, QUBO formulations have been proposed for
industry applications, such as nurse scheduling [25], job shop scheduling
[26], and vehicle routing [27]. Although they assumed quantum techniques such
as quantum annealing, our results suggest that the state-of-the-art QUBO
solver could be competitive with the general MIQP solver for them.

## References

  * [1] H. L. Urdaneta, J. Yuan, A. S. Siqueira, Alternative integer linear and quadratic programming formulations for ha-assignment problems, Proceeding Series of the Brazilian Society of Computational and Applied Mathematics 6 (1) (2018).  doi:10.5540/03.2018.006.01.0311.
  * [2] M. Kuramata, R. Katsuki, K. Nakata, Solving large break minimization problems in a mirrored double round-robin tournament using quantum annealing, Plos one 17 (4) (2022) e0266846.  doi:10.1371/journal.pone.0266846.
  * [3] D. Rehfeldt, T. Koch, Y. Shinano, Faster exact solution of sparse maxcut and qubo problems, Mathematical Programming Computation (2023) 1–26. doi:10.1007/s12532-023-00236-6.
  * [4] G. L. Nemhauser, M. A. Trick, Scheduling a major college basketball conference, Operations research 46 (1) (1998) 1–8.  doi:10.1287/opre.46.1.1.
  * [5] A. Drexl, S. Knust, Sports league scheduling: Graph-and resource-based models, Omega 35 (5) (2007) 465–471.  doi:10.1016/j.omega.2005.08.002.
  * [6] M. A. Trick, A schedule-then-break approach to sports timetabling, in: Practice and Theory of Automated Timetabling III. PATAT 2000. Lecture Notes in Computer Science, Vol. 2079, Springer, 2001, pp. 242–253.  doi:10.1007/3-540-44629-x_15.
  * [7] R. Miyashiro, T. Matsui, Semidefinite programming based approaches to the break minimization problem, Computers & Operations Research 33 (7) (2006) 1975–1982.  doi:10.1016/j.cor.2004.09.030.
  * [8] D. De Werra, Scheduling in sports, Studies on graphs and discrete programming 11 (1981) 381–395.  doi:10.1057/jors.1982.188.
  * [9] G. Post, G. J. Woeginger, Sports tournaments, home–away assignments, and the break minimization problem, Discrete Optimization 3 (2) (2006) 165–173.  doi:10.1016/j.disopt.2005.08.009.
  * [10] J.-C. Régin, Minimization of the number of breaks in sports scheduling problems using constraint programming, in: DIMACS workshop on on Constraint programming and large scale discrete optimization, 2000, pp. 115–130.  doi:10.1090/dimacs/057/07.
  * [11] M. Elf, M. Jünger, G. Rinaldi, Minimizing breaks by maximizing cuts, Operations Research Letters 31 (5) (2003) 343–349.  doi:10.1016/S0167-6377(03)00025-7.
  * [12] M. X. Goemans, D. P. Williamson, Improved approximation algorithms for maximum cut and satisfiability problems using semidefinite programming, Journal of the ACM (JACM) 42 (6) (1995) 1115–1145.  doi:10.1145/227683.227684.
  * [13] A. Suzuka, R. Miyashiro, A. Yoshise, T. Matsui, The home-away assignment problems and break minimization/maximization problems in sports scheduling, Pacific Journal of Optimization 3 (2007) 113–133.
  * [14] K. Easton, G. Nemhauser, M. Trick, Solving the travelling tournament problem: A combined integer programming and constraint programming approach, Lecture notes in computer science (2003) 100–112. doi:10.1007/978-3-540-45157-0_6.
  * [15] R. V. Rasmussen, Scheduling a triple round robin tournament for the best danish soccer league, European Journal of Operational Research 185 (2) (2008) 795–810.  doi:10.1016/j.ejor.2006.12.050.
  * [16] D. Yamaguchi, S. Imahori, R. Miyashiro, T. Matsui, An improved approximation algorithm for the traveling tournament problem, Algorithmica 61 (4) (2011) 1077–1091.  doi:10.1007/978-3-642-10631-6_69.
  * [17] R. Miyashiro, T. Matsui, S. Imahori, An approximation algorithm for the traveling tournament problem, Annals of Operations Research 194 (2012) 317–324.  doi:10.1007/s10479-010-0742-x.
  * [18] C. Thielen, S. Westphal, Complexity of the traveling tournament problem, Theoretical Computer Science 412 (4-5) (2011) 345–351.  doi:10.1016/j.tcs.2010.10.001.
  * [19] C. Thielen, S. Westphal, Approximation algorithms for TTP(2), Mathematical Methods of Operations Research 76 (1) (2012) 1–20.  doi:10.1007/s00186-012-0387-4.
  * [20] S. Imahori, A $1+\mbox{O}(1/n)$ approximation algorithm for TTP(2), in: Proceedings of the International Symposium on Scheduling 2015, 2015, pp. 186–191.
  * [21] M. Xiao, S. Kou, An improved approximation algorithm for the traveling tournament problem with maximum trip length two, in: 41st International Symposium on Mathematical Foundations of Computer Science (MFCS 2016), Schloss Dagstuhl-Leibniz-Zentrum fuer Informatik, 2016.  doi:10.4230/LIPIcs.MFCS.2016.89.
  * [22] R. Miyashiro, T. Matsui, A polynomial-time algorithm to find an equitable home–away assignment, Operations Research Letters 33 (3) (2005) 235–241.  doi:10.1016/j.orl.2004.06.004.
  * [23] F. Glover, G. Kochenberger, R. Hennig, Y. Du, Quantum bridge analytics I: a tutorial on formulating and using qubo models, Annals of Operations Research 314 (1) (2022) 141–183.  doi:10.1007/s10288-019-00424-y.
  * [24] T. P. Kirkman, On a problem in combinations, Cambridge and Dublin Mathematical Journal 2 (1847) 191–204.  doi:0.1112/PLMS/S1-2.1.69.
  * [25] K. Ikeda, Y. Nakamura, T. S. Humble, Application of quantum annealing to nurse scheduling problem, Scientific reports 9 (1) (2019) 1–10.  doi:10.1038/s41598-019-49172-3.
  * [26] C. Carugno, M. Ferrari Dacrema, P. Cremonesi, Evaluating the job shop scheduling problem on a d-wave quantum annealer, Scientific Reports 12 (1) (2022) 6539.  doi:10.1038/s41598-022-10169-0.
  * [27] H. Irie, G. Wongpaisarnsin, M. Terabe, A. Miki, S. Taguchi, Quantum annealing of vehicle routing problem with time, state and capacity, in: Quantum Technology and Optimization Problems: First International Workshop, QTOP 2019, Munich, Germany, March 18, 2019, Proceedings 1, Springer, 2019, pp. 145–156.  doi:10.1007/978-3-030-14082-3.