m[#1]

# Spatial Assisted Human-Drone Collaborative Navigation and Interaction
through Immersive Mixed Reality

Luca Morando and Giuseppe Loianno The authors are with the New York
University, Tandon School of Engineering, Brooklyn, NY 11201, USA. email:
{luca.morando<EMAIL_ADDRESS>work was supported by the NSF CPS Grant
CNS-2121391, the NSF CAREER Award 2145277, the DARPA YFA Grant D22AP00156-00,
Qualcomm Research, Nokia, and NYU Wireless.

###### Abstract

Aerial robots have the potential to play a crucial role in assisting humans
with complex and dangerous tasks. Nevertheless, the future industry demands
innovative solutions to streamline the interaction process between humans and
drones to enable seamless collaboration and efficient co-working. In this
paper, we present a novel tele-immersive framework that promotes cognitive and
physical collaboration between humans and robots through Mixed Reality (MR).
This framework incorporates a novel bi-directional spatial awareness and a
multi-modal virtual-physical interaction approaches. The former seamlessly
integrates the physical and virtual worlds, offering bidirectional egocentric
and exocentric environmental representations. The latter, leveraging the
proposed spatial representation, further enhances the collaboration combining
a robot planning algorithm for obstacle avoidance with a variable admittance
control. This allows users to issue commands based on virtual forces while
maintaining compatibility with the environment map. We validate the proposed
approach by performing several collaborative planning and exploration tasks
involving a drone and an user equipped with a MR headset.

## Supplementary Material

Github: https://github.com/arplaboratory/mri˙ros˙public.git Video:
https://youtu.be/oWdqLuXe1Zw

##  I Introduction

The increased use of robotics systems in our everyday life raises the need to
design novel solutions that can facilitate the interaction between humans and
robots decreasing the physical and cognitive human workload especially when
executing complex or dangerous tasks. Human-Robot Collaboration and
Interaction (HRI) is a well-studied problem and being integrated in multiple
areas, such as Industry 4.0, construction, and emergency responses [1, 2].
Traditionally, human-centered robotic systems focus on collaboration between
humans and grounded robots. In aerial robotics, due to the additional
complexity introduced by the navigation requirements in 3D unstructured
environments, current state of the art approaches still limit the
collaboration to be mostly unidirectional (i.e., with commands or information
sent from the human to the robot) using classic teleoperation [3] with
representation of the output signals projected on flat surface monitors.
Aerial robotics will become essential in supporting humans in multiple tasks
including, but not limited to warehouses’ operations, vertical farming
inspection, space missions, and search and rescue. In these settings, robots’
role can no longer be limited to simply assist humans. Instead, they will need
to co-exist, collaborate, and co-work with them, sharing the same workspace
and participating in shared tasks’ execution. Therefore, these requirements
call for the novel approaches that enable a high level of mutual and seamless
human-robot interaction. We envision aerial robots to be elevated to the role
of humans’ teammates allowing a bi-directional flux of information to
concurrently augment both human and robot situational awareness.

In this work, we propose a novel tele-immersive framework for human-drone
cognitive and physical collaboration through Mixed Reality (MR) as illustrated
in Fig. LABEL:fig:fig1 to facilitate seamless interaction between humans and
robots. The emergence of innovative spatial computing techniques, such as
Virtual Reality (VR), Augmented Reality (AR), and Mixed Reality (MR), presents
a technological opportunity in robotics. These techniques facilitate enhanced
collaboration between humans and robots through multi-modal information
sharing within the human-robot team based on vision [4], gestures [5], natural
languages [6], and gaze [7, 8].

The contributions of this paper can be summarized as follows. First, we
present a novel bi-directional spatial awareness concept that enables co-
sharing of spatial information between a human and an aerial robot at
different levels of abstractions, making them co-aware of the surrounding
environment. Second, we propose a novel virtual-physical interaction that
further facilitates bi-directional collaboration. This approach uniquely
combines a Variable Admittance Control (VAC) [9, 10], implemented on top of
the drone’s guidance loop, with a planning algorithm [11] responsible for
robot’s obstacle avoidance. The proposed VAC allows through gesture
recognition to input a user virtual external force as an external command to
the robot while ensuring compatibility with the environment map. The user
force is further combined with a repulsive force field acting on the robot
when located in proximity of obstacles increasing the overall system safety.

To the best of our knowledge, the proposed framework is the first to enable
continuous spatial virtual-physical navigation and interaction with an aerial
robot via MR. We open source the framework to the community.

##  II Related Works

Teleoperation. Drone teleoperation, while arguably one of the simplest
interaction methods, remains an open research area especially in constrained
robot environments [12, 13]. Current approaches still rely on simple
joysticks, keyboards, and pointers receiving the feedback from the robots on
monitors. Overall, state of the art autopilots can enable drones to follow
aggressive trajectories [14], but they do not offer an easy and intuitive
communication paradigm between human and robot. Recent solutions try to
overcome these drawbacks by employing hand gesture [15] or impedance control
with vibrotactile feedback [16] to control multiple drones or impedance
control for human-machine [17, 18] or machine-environment interaction [19].

Extended Reality. Innovative body-machine interfaces or wearable devices
(e.g., head mounted displays) utilize body motion to control robots and
incorporate haptic feedback, along with VR or AR [20, 21, 22, 23, 24]. These
interfaces aim to address the constraints of basic teleoperation control. For
instance, AR proves to be a valuable interaction modality that facilitates
collaborative information visualization in robotics tasks, therefore enabling
efficient communication of robots’ intentions to human co-workers [25, 26]. In
[27], the authors study how to maximize the user ability to tele-operate a
particular drone projecting its camera field of view camera in the user
headset. In [28, 29, 30], an exocentric view of the drone pose in a virtual
scene that is built a-priori, helps the user to interact with the drone via a
pick and place gesture, commanding a new goal position to the robot. A gesture
based trajectory completion in MR is proposed in [31], where an automated
technique auto-complete the user intended motion. Conversely, the authors in
[32] show a collaborative MR interface for immersive planning where the user
is able, using simple gestures to draw 3D trajectories for quadrotors in its
workspace. While in [32] the main focus is on drone teleoperation via simple
gestures, in [33] the authors employ MR to analyze the information sharing
problem, in order to simplify the human-robot communication during an
exploration task. Information about the robot’s future actions and predicted
trajectory is visualized using AR. Finally, a recent work [34] proposes a
fascinating instance of enhancing human situational awareness in MR during
human-robot interaction. The MR experience merges a downsized satellite map
with the robot’s local environment representation. Users have the ability to
choose a location on the map and stream the locally captured robot’s point
cloud into the MR headset, providing them with a first-person view of the
scene.

Novelties. Our tele-immersive approach goes beyond simple tele-operation and
extended reality for robotics. It facilitates a bi-directional flux of
information to enhance the spatial awareness of both the human and the drone.
This aspect is further enhanced incorporating a novel human-robot virtual-
physical interaction in 3D based on a VAC [35]. Compared to [35], our work
goes beyond the manipulator case, it eliminates the need for prior knowledge
about the environment’s structure and seamlessly and uniquely designs and
deploys the VAC within the online 3D navigation framework.

##  III Methodology

As visible in Fig. LABEL:fig:fig1.a, the robot body frame is denoted with $B$
and the robot world frame, located at the drone take-off position, is denoted
with $W$ and it is represented in the MR environment as a red and white
cylinder. Referring to Fig. LABEL:fig:fig1.a, we denote the spatial
transformation between two frames $i$ and $j$ as with $\mathbf{T}_{j}^{i}$
that represents the relative pose of the $j$ frame with respect to the $i$.

###  III-A Spatial Awareness

1

Figure 1: Key framework’s components. The pointcloud set $P^{W}$ is collected
and sent to the SAM where a collective octomap is generated (left). A VAC
coupled with a planner provides assistance to the user during the drone
navigation (right).

We introduce the concept of human-drone spatial awareness. This cutting-edge
technique facilitates an improved representation of the robot’s perceived
environment on the user side, utilizing advanced mapping and visualization
methods. Moreover, on the robot side, the user’s perception of the surrounding
environment is translated in a map representation that is interpretable by the
robot navigation layer. The raw data obtained from the drone and the wearable
device sensors are processed through the Spatial Awareness Module (SAM) as
presented in Fig. LABEL:fig:fig1. Specifically, the robot point cloud
${P_{c}}_{robot}$ (in our case a stereo disparity as detailed in Section IV)
and the wearable device point cloud ${P_{c}}_{dev}$ are forwarded as input to
the robot SAM which provides as output two types of geometric representations
namely an octomap and a mesh visualization shown in Fig. LABEL:fig:fig1.b and
Fig. LABEL:fig:fig1.c respectively.

The set of the robot and device pointclouds
${P_{c}}^{W}=\\{{P_{c}}_{dev},{P_{c}}_{robot}\\}^{W}$ can be defined and
expressed in $W$ using the following transformations: i) $\mathbf{T}_{B}^{W}$
for ${P_{c}}_{robot}$ and ii) $\mathbf{T}_{SR}^{W}$ for ${P_{c}}_{dev}$, where
the frame $SR$ is initialized with the same pose of the device frame $H$, when
the MR application is launched. Leveraging the Voxblox framework [36], it is
possible to obtain an octomap set from the pointcloud ${P_{c}}^{W}$, denoted
as $V^{W}=\\{V_{robot},V_{dev}\\}^{W}$, where $V_{robot}$ is the voxel
representation on the robot side while $V_{dev}$ is similarly obtained from
the user perspective. The set $V^{W}$ define a unified global voxel
representation visible in Fig. LABEL:fig:fig1.c where the merged maps of two
adjacent rooms (one explored by the robot and the other one explored by the
user) are depicted. An example of the representation instantiated by the SAM
is depicted in Fig.LABEL:fig:fig1.b and Fig.LABEL:fig:fig1.c, showcasing a
dual visualization of the environments perceived by the user and the robot,
respectively.

In Fig. LABEL:fig:fig1.b, we show the user Ego-centric visualization that
allows the MR layer to translate the drone perceived environment in a spatial
language easily understandable by the human. The Exo-centric visualization
instead, depicted in Fig. LABEL:fig:fig1.c, merges the user and the drone maps
captured from two different rooms or perspectives, therefore augmenting the
drone perception capabilities. A further level of abstraction and
functionality is provided by the instantiation of the MR environment in the
user field of view, where both the wearable device frame $H$ and the world
frame $W$, are referred to the $SR$ frame via $\mathbf{T}_{W}^{SR}$ and
$\mathbf{T}_{H}^{SR}$ respectively, as shown in Fig. LABEL:fig:fig1.a. The
wearable device can extract information in the perceived environment about
various surfaces (i.e., walls, ceilings, floors etc.) or artificially create
additional elements in form of holograms around the user (see Fig.
LABEL:fig:fig1.b). This produces a seamless mix of real and virtual worlds
allowing a customized and reconfigurable spatial interaction with increased
configuration complexities between the human and/or the robot even in an empty
physical space. An example is provided in Fig. LABEL:fig:fig1.b, where the
drone mesh $M_{robot}$, obtained by ${P_{c}}_{robot}$ and visible in white, is
merged with other spatially defined holograms in the MR framework.

Finally, the drone pose is always represented in the MR environment through an
orange cube $C$ referred to the frame $W$ (see Fig. LABEL:fig:fig1.a).

###  III-B Virtual-Physical Human-Drone Interactive Navigation

We propose a VAC to link the physical and virtual worlds producing two types
of virtual-physical interaction

  * •

Assisted Physical and Virtual Interaction (APVI): In this modality, a planning
algorithm (in our case a Rapidly exploring Random Tree in its optimal
formulation (RRT*) [11]) is built on top of the generated octomap and combined
with a VAC algorithm. Concurrently, the VAC provides a modulated damping
feedback based on a desired function. This can modify the drone compliance
depending the user intention to modify its path from the original one proposed
by the planner.

  * •

Free Physical and Virtual Interaction (FPVI): This is a specific case of the
APVI modality when the RRT* planner is disabled by the user. This modality
provides additional freedom when driving the robot. Since the robot is not
anymore constrained to follow a specific path, an additional safety layer is
required. We propose the generation of an obstacle force field that is
injected in the admittance controller as an additional input.

A schematic representation of both modalities with their integration with the
SAM is visualized in Fig. 1. The two modalities offer varying levels of
intuitiveness and environmental safety during the human-robot interaction,
allowing the user to depend more on the robot’s perception of its surroundings
for effective exploration, especially in unfamiliar or obstacle-filled
environments. In the following, we describe the key components of the APVI and
FPVI namely the RRT*, the VAC and how this is couple with the RRT* for the
APVI modality, and the obstacle force field.

####  III-B1 RRT* assisted Human-Robot Navigation

The RRT* is enabled in the APVI mode, expanding its tree based on a local
robot map defined within a circular horizon centered on the drone body frame
origin $\mathcal{O_{B}}$. This avoids the problem to know a-priori the
environment map and local minima planning issues can be taken care by the user
intervention with the proposed VAC. The output of the RRT* is a sequence of
setpoints $\mathbf{P}_{S_{i}}$, for $i=0,\cdots N$, which are computed only
within the drone horizon. The path is updated on-line, while the drone is
moving along the proposed direction. As visible in Fig. 1 (right), the robot
reference trajectory is composed by a sequence of reference positions
$\mathbf{p}_{r}$ that lies on the line each defined by a segment between two
consecutive setpoints $\mathbf{P}_{S_{i}}$ and $\mathbf{P}_{S_{i+1}}$
generated by the RRT*. The sequence of reference position $\mathbf{p}_{r}$
compose the trajectory $\mathbf{\Delta}_{r}$.

During the random nodes exploration phase embedded in the RRT* , the total
cost of the new node $c^{*}$ is obtained by adding the linking cost $c_{l}$ to
the parent cost $c_{p}$.

####  III-B2 Coupling Planning and VAC for assisted virtual-physical HRI

We combine the RRT* planning algorithm, presented in Section III-B1 with a
VAC. Prior to describing the VAC, we design an admittance controller on top of
the drone inner control position pipeline. This controller, in its classic
formulation visible in eq. (1), is driven by the 3D virtual input force
$\mathbf{F}_{v}$ generating the position perturbation
$\mathcal{X}=\mathbf{p}_{c}-\mathbf{p}_{r}$ between a commanded position
$\mathbf{p}_{c}$ and a reference position $\mathbf{p}_{r}$. The dynamical
relationship between the input force $\mathbf{F}_{v}$ and the position
perturbation $\mathbf{\mathcal{X}}$ mimics a mass spring damper system as

$\mathbf{M}\left(\ddot{\mathbf{p}}_{c}-\ddot{\mathbf{p}}_{r}\right)+\mathbf{D}\left(\dot{\mathbf{p}}_{c}-\dot{\mathbf{p}}_{r}\right)+\mathbf{K}\left(\mathbf{p}_{c}-\mathbf{p}_{r}\right)=\mathbf{F}_{v}.\leavevmode\nobreak\
$ (1)

The terms $\mathbf{M}$, $\mathbf{D}$ and $\mathbf{K}$ are $3\times 3$
(diagonal matrices in our case) representing the mass, the damping and the
stiffness acting on each axis of the world frame $W$. The right input term to
the equation is the virtual force
$\mathbf{F}_{v}=\textbf{F}_{usr}+\textbf{F}_{rep}$ resulting from the sum of
the user interaction force $\textbf{F}_{usr}$ with the obstacle repulsive
force $\textbf{F}_{rep}$, which will be defined later in Section III-B3, while
$\mathbf{F}_{usr}=\mathbf{K}_{p}(\mathbf{p}-\mathbf{p}_{u})-\mathbf{K}_{d}(\dot{\mathbf{p}}-\dot{\mathbf{p}}_{u})$
with $\mathbf{p}$ and $\mathbf{p}_{u}$ are the drone position and the user
interactive marker position respectively. These are estimated using a linear
Kalman Filter (KF) for motion estimation considering a constant acceleration
model which takes as input the acceleration $\ddot{\mathbf{p}}_{u}$. The terms
$\mathbf{K}_{p}$ and $\mathbf{K}_{d}$ represent respectively the proportional
and the derivative terms. Compared to the classic admittance control, the
damping value $\mathbf{D}$, in the VAC definition, varies according to the
direction of the input force $\mathbf{F}_{v}$. A variable damping enables to
change the feeling of the user during the interaction since the drone is
constrained to follow the direction of minimum resistance.

Starting from eq. (1) for the classic admittance controller and referring to
the RRT* segment in Fig. 1, let $\theta=[0,\pi]$ be the angle determined by
the human force direction $\mathbf{F}_{v}$ with respect to
$\mathbf{l}_{i}:\theta=\cos^{-1}(\mathbf{l}_{i}^{\top}\mathbf{F}_{v}/||\mathbf{F}_{v}||)$
where $\mathbf{l}_{i}$ is the vector describing the direction, expressed in
$W$, of the segment between two consecutive setpoints $\mathbf{P}_{S_{i-1}}$
and $\mathbf{P}_{S_{i}}$ traversed by the drone. Focusing on a generic segment
in Fig. 1, $\mathbf{p}_{u}$ is projected on the segment defining the position
$\mathbf{p}_{r}$, forwarded as input to the Variable Admittance Control (VAC).
The sequence of reference positions $\mathbf{p}_{r}$ defines the reference
trajectory $\mathbf{\Delta}_{r}$. Inspired by [35], we define the VAC with
three functions $D_{k}=f_{k}(\theta)$, with $D_{k}$ a diagonal element of
$\mathbf{D}$, as $linear$, $squared$ and $exponential$ with

  * •

$D_{k,min}=f_{k}(\theta=0)$,

  * •

$D_{k,max}=[f_{k}(\theta=\pi)\land f_{k}(\theta=-\pi)]$,

  * •

$f_{k}(\theta)$ is continuous in $[0,\pi]$ and $[0,-\pi]$.

The three damping functions $D_{k}=f_{k}(\theta)$ are integrated within the
APVI mode. The functions modify the damping behaviour (visualized as the blue
area around the drone position in Fig. 1) depending on the pointing direction
of the virtual force vector $\mathbf{F}_{v}$ that is forwarded as input to the
admittance controller as shown in eq. (1).

The set of commanded positions $\mathbf{p}_{c}$ compose the trajectory
$\mathbf{\Delta}_{c}$, visible as a blue dashed line in Fig. 1 (right),
compliant with the direction of the force $\mathbf{F}_{v}$.

####  III-B3 Obstacle Perception as Force Field

The FPVI modality employs the same spatial representation of the APVI case as
visible in Fig. 1 (left). Given multiple repulsive forces per voxel
$\mathbf{F}_{r_{i}}$ with $i\in\\{0,...,N_{v}\\}$ with $N_{v}$ is the total
number of voxels $v_{i}\in V^{W}$ within the robot horizon $h$, it is possible
to define the direction of the resulting force $\mathbf{F}_{rep}$ acting on
the robot frame $B$ and expressed in the frame $W$. The normalized vector
$\mathbf{F}_{rep}=\mathbf{F}_{r}/\|\mathbf{F}_{r}\|$ defines the direction of
the repulsive force acting on the robot, where
$\mathbf{F}_{r}=\sum_{i=0}^{Nv}\mathbf{F}_{r_{i}}$ is composed by the sum of
all the force vectors $\mathbf{F}_{r_{i}}$ generated by each of the voxel
$v_{i}$ inside the robot horizon $h$. Each term $\mathbf{F}_{r_{i}}$ is
therefore defined following the exponential decay behaviour in function of the
euclidean distance $d$ from the voxel $v_{i}$ as

$\|\mathbf{F}_{r_{i}}\|=\frac{F_{s}}{k}e^{(-\lambda d)}(1-e^{(h-d)}),$ (2)

where $k=1-e^{h}$, $\|\mathbf{F}_{r_{i}}\|=F_{s}$ when $d=0$ and $\lambda$ is
the exponential decay constant. Finally, at each iteration, the magnitude of
the vector $\mathbf{F}_{rep}$ is defined as
$\|\mathbf{F}_{rep}\|=\|\mathbf{F}_{r}\|$.

Figure 2: Ego and Exo-centric visualization of the proposed experiment
captured at different salient moments. An enhanced safety layer is provided by
the robot and the SAM, with the generation of a safe trajectory through the
obstacles.

##  IV Experimental Results

We validate the components of the proposed framework and show their
application in joint human-drone exploration tasks through multiple real-world
experiments. These are conducted at the Agile Robotics and Perception Lab
(ARPL) at the New York University in the flying arena with size of $60$
$\mathrm{m}^{2}$ as well as in an adjacent room to show the fully capabilities
of the tele-immersive framework in collaborative tasks, as proposed in Fig.
LABEL:fig:fig1.b and Fig. LABEL:fig:fig1.c.

The quadrotor platform used in the experiments is equipped with a
$\text{Qualcomm}^{\circledR}\text{Snapdragon}^{\text{TM}}$ board for onboard
computing where an embedded stereo camera is used for a disparity
representation, obtaining the point cloud ${{P_{c}}_{robot}}^{B}$. Our
algorithm provides the state of the robot (pose and perception output)
directly in the $W$ frame via the transformation $\mathbf{T}_{B}^{W}$. In
addition, a mutual co-localization mechanism between the MR device and the
aerial robots can be used[37]. Both the robot and the AHRMD use onboard self-
localization mechanism. The Mixed Reality framework is implemented in C$\\#$
and it is executed on the Microsoft® Hololens 2.0. Combining different
features, the proposed framework offers four types of deployment modalities

  * •

Full Simulation (FS) The complete interaction is based on ARPL custom
quadrotor simulator and visualized in ROS RViz.

  * •

Mixed Simulation (MS): It is the same as FS, but the user can interact with
the robot using the MR headest.

  * •

Real Interaction (RI): The user employs ROS Rviz to send commands to the real
flying drone. In this case, the mapping and navigation functionalities are
enabled.

  * •

Real Interaction in MR (RIMR): In this case the framework is fully enabled.

We provide an experimental overview of the FPVI and APVI interaction
modalities proposed in Section III-B and shown in Fig. 1. These modalities are
tested in a lab environment by four individuals, all of whom provided positive
feedback regarding the projected amount of visual information and the ease of
interaction with an aerial robot. They consistently achieved the assigned task
to reach a target in a densely populated obstacle environment, relying on both
the robot perception feedback represented by the computed planner path and the
obstacle repulsive force. A visualization of two of the possible tasks
assigned to the subjects are described in the following sections.

###  IV-A VAC with RRT* exploration and Mapping

The interaction behaviour during the APVI mode is represented in Fig. 2
showing a sequence of frames with an Egocentric and Exocentric views at
different time steps during the same experiment whose results are in Fig. 3.

During the experiments the height of the robot is kept constant at
$p_{r,z}=0.8\leavevmode\nobreak\ $\mathrm{m}$$. The black dashed line at time
instant $t_{APVI}=19.5$ $\mathrm{s}$, in Fig. 3, represents the switches from
(FPVI) to (APVI) mode which happens when the user places a final goal to the
robot to reach, denoted as $\mathbf{G}^{W}$ and visualized in Fig. 2 with a
red circle. Once activated, the RRT* obstacle aware planning algorithm,
coupled with the VAC initialized with $M_{k}=2.4\leavevmode\nobreak\
$\mathrm{K}\mathrm{g}$$, $K_{k}=20\leavevmode\nobreak\
$\mathrm{K}\mathrm{g}\cdot\mathrm{s}^{-2}$$ and $1\leavevmode\nobreak\
$\mathrm{K}\mathrm{g}\cdot\mathrm{s}^{-1}$\leq D_{k}\leq
70\leavevmode\nobreak\ $\mathrm{K}\mathrm{g}\cdot\mathrm{s}^{-1}$$, where the
subscript $k$ represents the diagonal element of $\mathbf{M}$, $\mathbf{K}$,
and $\mathbf{D}$ respectively, constraints the drone to follow the proposed
path. Nonetheless, the user can still pull away the robot from the initial
trajectory despite this safety condition, applying an increasing virtual Force
$\mathbf{F}_{v}$ as input to the system, in order to define a maneuver to
bypass virtual and real obstacles. This makes the robot movement compliant to
the user intentions but still constrained to follow the reference trajectory
$\mathbf{\Delta}_{r}$, relying on the planner proposed path. The maneuver and
the human-drone assisted navigation are visible in detail also in the sequence
of frames in Fig. 2.

To show the system’s compliance to the human intentions, in Fig. 2.b, we
observe that the user’s action can pull the robot away from the RRT* planned
path and drive it along a new trajectory passing through a narrow gap measured
in only $0.68\leavevmode\nobreak\ $\mathrm{m}$$ (schematized as a green
arrow). The estimated obstacle Truncated Signed Distance Field (TSDF) [38,
39], represented by the green mesh $M_{robot}$, visualizes the detected
surface of the obstacle as perceived by the drone sensors during the bypass
maneuver. Once this maneuver terminates, in Fig. 2.c, the RRT* restarts its
expansion towards $\mathbf{G}^{W}$, only when the drone is attracted back to
the last generated setpoint. Moreover, in Figs. 2.a, 2.b, and 2.c the
interaction marker (white cube), the drone camouflage (orange cube) and the
generated force (colorful string) are easily identifiable by the user. The
generated trajectory depending on the Force $\mathbf{F}_{v}$ is shown in Fig.
3. The Root Mean Squared Error between $\mathbf{\Delta}_{r}$ and
$\mathbf{\Delta}_{c}$ during the FPVI and APVI modalities is respectively of
$0.0373$ $\mathrm{m}$ and $2.0142$ $\mathrm{m}$. The higher RMSE value
obtained during the APVI phase is due to the obstacle bypass maneuver, where
the commanded trajectory is compliant with the user intentions rather than
following the original desired path.

###  IV-B FPVI with Obstacle Force Field

We propose a second experiment which provides a deeper insight of the behavior
of the system when subject to the action of the obstacle generated force
field, defined in Section III-B3 when the system is set to the FPVI mode.
Along the overall experiment, the gains of the admittance controller, except
for $D_{k}=20\leavevmode\nobreak\ $\mathrm{K}\mathrm{g}\cdot\mathrm{s}^{-1}$$,
are not modified with respect the previous scenario shown in Section 2. The
exponential decay behaviour of the force $\mathbf{F}_{rep}$ is defined
respectively with $Fs=8\leavevmode\nobreak\ $\mathrm{N}$$, $\lambda=1$ and
horizon $h=1.5\leavevmode\nobreak\ $\mathrm{m}$$. Referring to Section 2 and
following the exponential decay physical behavior, the magnitude of the force
$\|\mathbf{F}_{rep}\|=0$ when $d=h$ and $\|\mathbf{F}_{rep}\|=F_{s}$ when
$d=0$. The results of this experiment are provided in Fig. 4, with a detailed
Ego (”User view”) and Exo-centric (”Rviz view”) representations of the
interaction captured at time $t=57\leavevmode\nobreak\ $\mathrm{s}$$. As
visible from the location of the white interaction cube, visualized as a
sphere in the map, the user tries to impose an unsafe desired trajectory
$\mathbf{\Delta}_{r}$ passing through an obstacle. As the robot perceives the
environment a repulsive force $\mathbf{F}_{rep}$, defined as in III-B3, is
generated and visualized in the Exo-centric representation in Fig. 4 as a
green arrow. The user pulling force and the obstacle repulsive force are
forwarded as input to the admittance controller which generates new safer
commanded trajectory $\mathbf{\Delta}_{c}$, visualized in Fig. 4 as a green
line. The complete evolution of the trajectories $\mathbf{\Delta}_{c}$,
$\mathbf{\Delta}_{r}$ and the robot actual sequence of positions $\mathbf{p}$,
are represented in the plot in Fig. 4 where the components of the trajectories
along the $x,\leavevmode\nobreak\ y$ axis are visualized. The ”User View” in
Fig. 4 provides an enhanced visualization of the human robot interaction at
the time of the representation in the ”Rviz view” as seen through the head
mounted device.

Figure 3: VAC behavior during the APVI interaction.

###  IV-C Results Discussion

Figure 4: Quantitative visualization of the obstacle safety layer acting on
the system visualized in User and Rviz view

The primary objective of this work is to validate the technical modules of
proposed framework. The presented results validate each component and the
overall framework’s functionality, showcasing potential for natural and
intuitive interaction beyond traditional tele-operation and basic MR
solutions. Tests in a cluttered environment show that safety measures
introduced through the planner algorithm (Section IV-A) and obstacle force
field (Section IV-B) enable secure, gesture-based interaction with the drone,
adaptable to user needs. Tests performed by multiple users further validate
the methodologies and the benefits of the proposed solution. This approach
serves as a flexible building block, open to incorporation with other
modalities for enhanced bi-directional data interpretation and interaction.

##  V Conclusion

In this paper, we presented a tele-immersive human-robot collaborative
framework. The key features include a novel bi-directional spatial awareness
approach that provides seamless human-robot spatial awareness through MR and a
multi-modal virtual-physical interaction that further enhances this
collaboration experience. We open-source the framework with the goal to
promote research in this area and incorporate feedback from the community.

Despite four users already validated the usability of the system in the lab,
in the future, we expect to conduct a complete user case study, leveraging the
NASA Task Load Index [40], subjective questions and quantitative data
collected on an inspection target reaching based task, in order to have a
complete performance evaluation on multiple users. This will help to further
validate the presented approaches and results. Specifically, we would like to
compare the users’ feelings during a 2D based interaction and full MR
experience. We believe the results will guarantee an easier manipulation and a
continuity in the interaction between the robot and non expert pilots when
both agents share the same environment. Finally, we aim to extend the proposed
framework by allowing the incorporation of heterogeneous teams consisting of
multiple users and robots.

## References

  * [1] A. Ajoudani, A. M. Zanchettin, S. Ivaldi, A. Albu-Schäffer, K. Kosuge, and O. Khatib, “Progress and prospects of the human–robot collaboration,” _Autonomous Robots_ , vol. 42, no. 5, pp. 957–975, Jun 2018.
  * [2] L. Morando, C. T. Recchiuto, and A. Sgorbissa, “Social drone sharing to increase the uav patrolling autonomy in emergency scenarios,” in _29th IEEE International Conference on Robot and Human Interactive Communication (RO-MAN)_ , 2020, pp. 539–546.
  * [3] S. Kim, “Working with robots: Human resource development considerations in human–robot interaction,” _Human Resource Development Review_ , vol. 21, no. 1, pp. 48–74, 2022.
  * [4] R.-J. Halme, M. Lanz, J. Kämäräinen, R. Pieters, J. Latokartano, and A. Hietanen, “Review of vision-based safety systems for human-robot collaboration,” _Procedia CIRP_ , vol. 72, pp. 111–116, 2018, 51st CIRP Conference on Manufacturing Systems.
  * [5] H. Liu and L. Wang, “Gesture recognition for human-robot collaboration: A review,” _International Journal of Industrial Ergonomics_ , vol. 68, pp. 355–367, 2018.
  * [6] S. Pate, W. Xu, Z. Yang, M. Love, S. Ganguri, and L. L. S. Wong, “Natural language for human-robot collaboration: Problems beyond language grounding,” _CoRR_ , vol. abs/2110.04441, 2021. [Online]. Available: https://arxiv.org/abs/2110.04441
  * [7] L. Yuan, C. Reardon, G. Warnell, and G. Loianno, “Human gaze-driven spatial tasking of an autonomous mav,” _IEEE Robotics and Automation Letters_ , vol. 4, no. 2, pp. 1343–1350, 2019.
  * [8] M. Pavliv, F. Schiano, C. Reardon, D. Floreano, and G. Loianno, “Tracking and relative localization of drone swarms with a vision-based headset,” _IEEE Robotics and Automation Letters_ , vol. 6, no. 2, pp. 1455–1462, 2021.
  * [9] A. Lecours, B. Mayer-St-Onge, and C. Gosselin, “Variable admittance control of a four-degree-of-freedom intelligent assist device,” in _IEEE International Conference on Robotics and Automation (ICRA)_ , 2012, pp. 3903–3908.
  * [10] F. Ferraguti, C. T. Landi, L. Sabattini, M. Bonfè, C. Fantuzzi, and C. Secchi, “A variable admittance control strategy for stable physical human–robot interaction,” _The International Journal of Robotics Research_ , vol. 38, no. 6, pp. 747–765, 2019.
  * [11] X. Wang, X. Luo, B. Han, Y. Chen, G. Liang, and K. Zheng, “Collision-free path planning method for robots based on an improved rapidly-exploring random tree algorithm,” _Applied Sciences_ , vol. 10, no. 4, 2020.
  * [12] W. A. Isop, C. Gebhardt, T. Nägeli, F. Fraundorfer, O. Hilliges, and D. Schmalstieg, “High-level teleoperation system for aerial exploration of indoor environments,” _Frontiers in Robotics and AI_ , vol. 6, 2019.
  * [13] B. Rojas, H. Abaunza, P. Castillo, I. Thouvenin, and R. Lozano, “Teleoperation of a drone based on a virtual environment,” in _Workshop on Research, Education and Development of Unmanned Aerial Systems (RED UAS)_ , 2019, pp. 160–165.
  * [14] A. Saviolo, G. Li, and G. Loianno, “Physics-inspired temporal learning of quadrotor dynamics for accurate model predictive trajectory tracking,” _IEEE Robotics and Automation Letters_ , vol. 7, no. 4, pp. 10 256–10 263, 2022.
  * [15] J. Nagi, A. Giusti, L. M. Gambardella, and G. A. Di Caro, “Human-swarm interaction using spatial gestures,” in _IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_ , 2014, pp. 3834–3841.
  * [16] E. Tsykunov, R. Agishev, R. Ibrahimov, L. Labazanova, A. Tleugazy, and D. Tsetserukou, “Swarmtouch: Guiding a swarm of micro-quadrotors with impedance control using a wearable tactile interface,” _IEEE Transactions on Haptics_ , vol. 12, no. 3, pp. 363–374, 2019.
  * [17] A. Tagliabue, M. Kamel, S. Verling, R. Siegwart, and J. Nieto, “Collaborative transportation using mavs via passive force control,” in _IEEE International Conference on Robotics and Automation (ICRA)_ , 2017, pp. 5766–5773.
  * [18] F. Augugliaro and R. D’Andrea, “Admittance control for physical human-quadrocopter interaction,” 07 2013, pp. 1805–1810.
  * [19] M. Car, A. Ivanovic, M. Orsag, and S. Bogdan, “Position-based adaptive impedance control for a uav,” in _International Conference on Unmanned Aircraft Systems (ICUAS)_ , 2018, pp. 957–963.
  * [20] M. Macchini, F. Schiano, and D. Floreano, “Personalized telerobotics by fast machine learning of body-machine interfaces,” _IEEE Robotics and Automation Letters_ , vol. 5, no. 1, pp. 179–186, 2019.
  * [21] J. Bruce, J. Perron, and R. Vaughan, “Ready—Aim—Fly! Hands-Free Face-Based HRI for 3D Trajectory Control of UAVs,” in _14th Conference on Computer and Robot Vision_ , 2017.
  * [22] J. Nagi, A. Giusti, G. A. Di Caro, and L. M. Gambardella, “Human control of uavs using face pose estimates and hand gestures,” in _ACM/IEEE International Conference on Human-Robot Interaction_ , 2014, pp. 252–253.
  * [23] M. Walker, H. Hedayati, J. Lee, and D. Szafir, “Communicating robot motion intent with augmented reality,” in _ACM/IEEE International Conference on Human-Robot Interaction_ , 2018, pp. 316–324.
  * [24] C. Reardon, K. Lee, and J. Fink, “Come see this! augmented reality to enable human-robot cooperative search,” in _IEEE International Symposium on Safety, Security, and Rescue Robotics (SSRR)_ , 2018.
  * [25] M. D. Coovert, T. Lee, I. Shindev, and Y. Sun, “Spatial augmented reality as a method for a mobile robot to communicate intended movement,” _Computers in Human Behavior_ , vol. 34, pp. 241–248, 2014.
  * [26] E. Rosen, D. Whitney, E. Phillips, G. Chien, J. Tompkin, G. Konidaris, and S. Tellex, “Communicating and controlling robot arm motion intent through mixed-reality head-mounted displays,” _The International Journal of Robotics Research_ , vol. 38, no. 12-13, pp. 1513–1526, 2019.
  * [27] P. Fletcher, A. Luther, C. Detweiler, and B. Duncan, “Predicting visual differentiability for unmanned aerial vehicle gestures,” _IEEE Robotics and Automation Letters_ , vol. 7, no. 4, pp. 11 799–11 806, 2022.
  * [28] O. Erat, W. A. Isop, D. Kalkofen, and D. Schmalstieg, “Drone-augmented human vision: Exocentric control for drones exploring hidden areas,” _IEEE Transactions on Visualization and Computer Graphics_ , vol. 24, no. 4, pp. 1437–1446, 2018.
  * [29] C. Liu and S. Shen, “An augmented reality interaction interface for autonomous drone,” in _IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_ , 2020, pp. 11 419–11 424.
  * [30] J. Betancourt, B. Wojtkowski, P. Castillo, and I. Thouvenin, “Exocentric control scheme for robot applications: An immersive virtual reality approach,” _IEEE Transactions on Visualization and Computer Graphics_ , pp. 1–1, 2022.
  * [31] M. K. Zein, M. Al Aawar, D. Asmar, and I. H. Elhajj, “Deep learning and mixed reality to autocomplete teleoperation,” in _IEEE International Conference on Robotics and Automation (ICRA)_ , 2021, pp. 4523–4529.
  * [32] A. Angelopoulos, A. Hale, H. Shaik, A. Paruchuri, K. Liu, R. Tuggle, and D. Szafir, “Drone brush: Mixed reality drone path planning,” in _ACM/IEEE International Conference on Human-Robot Interaction (HRI)_ , 2022, p. 678–682.
  * [33] C. Reardon, K. Lee, J. G. Rogers, and J. Fink, “Augmented reality for human-robot teaming in field environments,” in _Virtual, Augmented and Mixed Reality. Applications and Case Studies_ , J. Y. Chen and G. Fragomeni, Eds. Cham: Springer International Publishing, 2019, pp. 79–92.
  * [34] M. Walker, Z. Chen, M. Whitlock, D. Blair, D. A. Szafir, C. Heckman, and D. Szafir, “A mixed reality supervision and telepresence interface for outdoor field robotics,” in _IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_ , 2021, pp. 2345–2352.
  * [35] D. Bazzi, G. Priora, A. M. Zanchettin, and P. Rocco, “Rrt* and goal-driven variable admittance control for obstacle avoidance in manual guidance applications,” _IEEE Robotics and Automation Letters_ , vol. 7, no. 2, pp. 1920–1927, 2022.
  * [36] H. Oleynikova, Z. Taylor, M. Fehr, R. Siegwart, and J. Nieto, “Voxblox: Incremental 3d euclidean signed distance fields for on-board mav planning,” in _IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_ , 2017, pp. 1366–1373.
  * [37] Microsoft, “azure_spatial_anchors_ro,” https://github.com/microsoft/azure_spatial_anchors_ros, 2020.
  * [38] B. Curless and M. Levoy, “A volumetric method for building complex models from range images,” in _Proceedings of the 23rd Annual Conference on Computer Graphics and Interactive Techniques_ , ser. SIGGRAPH ’96. New York, NY, USA: Association for Computing Machinery, 1996, p. 303–312.
  * [39] D. Werner, A. Al-Hamadi, and P. Werner, “Truncated signed distance function: Experiments on voxel size,” in _Image Analysis and Recognition: 11th International Conference, ICIAR 2014_. Springer-Verlag, 2014, p. 357–364.
  * [40] S. G. Hart and L. E. Staveland, “Development of nasa-tlx (task load index): Results of empirical and theoretical research,” _Human mental workload_ , vol. 1, no. 3, pp. 139–183, 1988.


# On integer linear combinations of terms of rational cycles for the
generalized 3x+1 problem

Yagub N. Aliyev
School of IT and Engineering
ADA University
Ahmadbey Aghaoglu str. 61
Baku 1008, Azerbaijan
E-mail<EMAIL_ADDRESS>
###### Abstract

In the paper, some special linear combinations of the terms of rational cycles
of generalized Collatz sequences are studied. It is proved that if the
coefficients of the linear combinations satisfy some conditions then these
linear combinations are integers. The discussed results are demonstrated on
some examples. In some particular cases the obtained results can be used to
explain some patterns of digits in $p$-adic representations of the terms of
the rational cycles.

††2020 _Mathematics Subject Classification_ : Primary 11B75; Secondary 11A63,
11D88.†† _Key words and phrases_ : $3x+1$ problem, Collatz conjecture,
rational cycles, integer linear combinations.

##  1 Introduction

Collatz conjecture or $3x+1$ problem claims that for any positive integer
$x_{0}$, the recursive sequence defined for $n\geq 0$ by
$x_{n+1}=S(x_{n})=\frac{3x_{n}+1}{2}$ if $x_{n}$ is odd and
$x_{n+1}=T(x_{n})=\frac{x_{n}}{2}$ if $x_{n}$ is even, there is a positive
integer $N$ such that $x_{N}=1$ (see [7]). It is known that this holds true
for almost all $x_{0}$ in the sense of some density. See [15], [4], [5], [14],
and the references therein for the works in this direction. Another approach
is to find all cycles of this recursive sequence. It is conjectured that there
are only finitely many such cycles. The only known cycle is the one generated
by $x_{0}=1$. There are more cycles if $x_{0}$ is allowed to be zero or a
negative integer but it is conjectured that their number is also finite (The
Finite Cycles Conjecture, see [7], [8]). The only known non-positive cycles
are the ones generated by $x_{0}=0$, $x_{0}=-1$, $x_{0}=-5$, and $x_{0}=-17$.
If a sequence of functions consisted of several $S$ and $T$ is given, then one
can speak about rational cycles (see [2]). There is a rational number $x_{0}$
such that if the operations $S$ and $T$ are applied in the given order, then
the final result is again $x_{0}$. Rational cycles generated by such $x_{0}$
have some interesting properties (see [8]). Many generalizations of Collatz
conjecture were considered by replacing $S$ and $T$ operations by more general
$S_{k}(x)=\frac{p_{i}x+k}{q}$. One can find such generalizations in [11],
[12], [10]. See [3], [9], where similar generalizations are used to prove some
results related to undecidability properties. See [6], [13], where these
generalizations are considered in the context of $2$-adic numbers and $q$
based numeral systems.

In this paper we focus on properties of the terms of rational cycles $x_{i}$,
which show that these rational numbers are "integer like". A linear
combination with integer coefficients of two integers is again an integer. In
the paper it is proved that under some conditions over the coefficients of the
linear combination, the rational numbers $x_{i}$ and $x_{i+b}$ also form
integer linear combinations. Two worked out examples demonstarting the results
on particular cases are given. Some applications of these results explaining
peculiar patterns of digits in $p$-adic representations of these rational
numbers $x_{i}$ are also discussed.

##  2 Notations and lemmas

Consider composition $P=B_{0}\circ B_{1}\circ\ldots\circ B_{n-1}$ of functions
$B_{i}(x)=\frac{p_{i}x+k_{i}}{q}$, where $n>1$, $k_{i}$ are integers,
$p_{i},q$ are non-zero integers such that $(p_{i},q)=1$ for
$i=0,1,\ldots,n-1$. When it is necessary to extend the index $i$, beyond the
interval $[0,n-1]$, we suppose that $B_{i}=B_{j}$ if $i\equiv j\ (\bmod\ n)$.
Consider equation $B_{0}\circ B_{1}\circ\ldots\circ B_{n-1}(x)=x$, which can
also be written as

$\frac{p_{0}\frac{p_{1}\frac{\frac{p_{n-2}\frac{p_{n-1}x+k_{n-1}}{q}+k_{n-2}}{\begin{matrix}q\\\
\vdots\end{matrix}}}{q}+k_{1}}{q}+k_{0}}{q}=x.$

Note that its solution $x_{0}$ is a rational number (cf. formula 1.2 and 1.3,
[8]):

$x_{0}=\frac{p_{0}p_{1}\ldots p_{n-2}k_{n-1}+p_{0}p_{1}\ldots
p_{n-3}k_{n-2}q+\ldots+p_{0}k_{1}q^{n-2}+k_{0}q^{n-1}}{q^{n}-p_{0}p_{1}\ldots
p_{n-1}}.$

Similarly, consider equations $B_{i}\circ B_{i+1}\circ\ldots\circ
B_{i+n-1}(x)=x$ for $i=0,1,\ldots,n-1$. Note that their solutions $x_{i}$ are
also rational numbers:

$x_{i}=\frac{p_{i}p_{i+1}\ldots p_{i+n-2}k_{i-1}+p_{i}p_{i+1}\ldots
p_{i+n-3}k_{i-2}q+\ldots+p_{i}k_{i+1}q^{n-2}+k_{i}q^{n-1}}{q^{n}-p_{0}p_{1}\ldots
p_{n-1}},$

where all the indices are taken modulo $n$. In the following, it is assumed
that $x_{i}=x_{j}$ if $i\equiv j(\bmod n)$. Consider also numbers
$U_{i}=\frac{q^{i}}{q^{n}-p_{0}p_{1}\ldots p_{n-1}},$ for $i=0,1,\ldots,n$.
Note that $U_{n}=p_{0}p_{1}\ldots p_{n-1}U_{0}+1$. Note also that

$x_{i}=p_{i}p_{i+1}\ldots p_{i+n-2}k_{i-1}U_{0}+p_{i}p_{i+1}\ldots
p_{i+n-3}k_{i-2}U_{1}+\ldots+p_{i}k_{i+1}U_{n-2}+k_{i}U_{n-1}.$

Note that there are infinitely many pairs of non-zero integers $\alpha,\beta$
and integers $b$, such that $0<b<n$ and $\alpha U_{0}+\beta U_{b}$ is an
integer or equivalently, $q^{n}-p_{0}p_{1}\ldots p_{n-1}|\alpha+\beta q^{b}$.
Indeed, one can take, for example, $\alpha=k$,
$\beta=-kq^{\phi(|q^{n}-p_{0}p_{1}\ldots p_{n-1}|)-1}$, and $b=1$, where
$k=1,2,\ldots$ and $\phi$ is Euler’s totient function. If $\alpha$ and $\beta$
are fixed then one can ask if such $b$ exists. The answer to this question
depends on the choice of $\alpha$ and $\beta$ . For example, if $q=2$, $n=4$,
$p_{0}=3$, $p_{1}=p_{2}=p_{3}=1$, then $q^{n}-p_{0}p_{1}p_{2}p_{3}=13$. If
$\alpha=1,\beta=1$ then such $b$ ($0<b<4$) does not exist. But if
$\alpha=9,\beta=1$ then $b=2$ satisfies the condition.

######  Lemma 2.1.

If $\alpha U_{0}+\beta U_{b}$ is an integer, then $p_{0}p_{1}\ldots
p_{n-1}\beta U_{0}+\alpha U_{n-b}$ is also an integer.

###### Proof.

Our claim is equivalent to prove that if $q^{n}-p_{0}p_{1}\ldots
p_{n-1}|\alpha+\beta q^{b}$ then

$q^{n}-p_{0}p_{1}\ldots p_{n-1}|p_{0}p_{1}\ldots p_{n-1}\beta+\alpha q^{n-b}.$

Indeed, since $(q,q^{n}-p_{0}p_{1}\ldots p_{n-1})=1$ and

$q^{b}(p_{0}p_{1}\ldots p_{n-1}\beta+\alpha q^{n-b})=p_{0}p_{1}\ldots
p_{n-1}\beta q^{b}+\alpha q^{n},$

it is sufficient to show that

$q^{n}-p_{0}p_{1}\ldots p_{n-1}|p_{0}p_{1}\ldots p_{n-1}\beta q^{b}+\alpha
q^{n},$

which follows from

$p_{0}p_{1}\ldots p_{n-1}\beta q^{b}+\alpha q^{n}=p_{0}p_{1}\ldots
p_{n-1}(\alpha+\beta q^{b})+\alpha(q^{n}-p_{0}p_{1}\ldots p_{n-1}).$

∎

######  Lemma 2.2.

If $\alpha U_{0}+\beta U_{b}$ is an integer, then for $i=0,1,2,\ldots$ the
numbers $\alpha U_{i}+\beta U_{i+b}$ and $p_{0}p_{1}\ldots p_{n-1}\beta
U_{i}+\alpha U_{n+i-b}$ are also integers.

###### Proof.

By mutiplying $\alpha U_{0}+\beta U_{b}$ and $p_{0}p_{1}\ldots p_{n-1}\beta
U_{0}+\alpha U_{n-b}$, which are both integers, by integer $q^{i}$, we obtain
that both of the numbers $\alpha U_{i}+\beta U_{i+b}$ and $p_{0}p_{1}\ldots
p_{n-1}\beta U_{i}+\alpha U_{n+i-b}$ are integers. ∎

##  3 Main results

######  Theorem 3.1.

If $\alpha U_{0}+\beta U_{b}$ is an integer, then for any $i$, satisfying
$0\leq i<i+b<n$, the number $\alpha x_{i}+\beta p_{i}p_{i+1}\ldots
p_{i+b-1}x_{i+b}$ is also an integer.

###### Proof.

Note that we can write

$\alpha x_{i}+\beta p_{i}p_{i+1}\ldots p_{i+b-1}x_{i+b}$
$=\alpha(p_{i}p_{i+1}\ldots p_{i+n-2}k_{i-1}U_{0}+p_{i}p_{i+1}\ldots
p_{i+n-3}k_{i-2}U_{1}+\ldots+p_{i}k_{i+1}U_{n-2}+k_{i}U_{n-1})$ $+\beta
p_{i}p_{i+1}\ldots p_{i+b-1}(p_{i+b}p_{i+b+1}\ldots p_{i+b+n-2}k_{i+b-1}U_{0}$
$+p_{i+b}p_{i+b+1}\ldots
p_{i+b+n-3}k_{i+b-2}U_{1}+\ldots+p_{i+b}k_{i+b+1}U_{n-2}+k_{i+b}U_{n-1})$
$=k_{0}M_{0}+k_{1}M_{1}+\ldots+k_{n-1}M_{n-1},$

where

$M_{j}=\begin{cases}p_{i}p_{i+1}\ldots p_{j-1}(\alpha U_{n+i-1-j}+\beta
p_{0}p_{1}\ldots p_{n-1}U_{i+b-1-j})&\text{if }i\leq j<i+b,\\\
p_{i}p_{i+1}\ldots p_{n+j-1}(\alpha U_{i-1-j}+\beta U_{i+b-1-j})&\text{if
}0\leq j<i,\\\ p_{i}p_{i+1}\ldots p_{j-1}(\alpha U_{n+i-1-j}+\beta
U_{n+i+b-1-j})&\text{if }i+b\leq j<n.\end{cases}$

By Lemma 2.1 and Lemma 2.2, all $M_{j}$, for $j=0,1,\ldots n-1$, are integers
and therefore the claim is true. ∎

######  Corollary 3.2.

If $\alpha U_{0}+\beta U_{b}$ is an integer, then for any $i$, the number
$\alpha x_{i}+\beta p_{i}p_{i+1}\ldots p_{i+b-1}x_{i+b}$ is also an integer.

###### Proof.

Since it was assumed that $x_{i}=x_{j}$ if $i\equiv j(\bmod n)$, without loss
of generality, we can suppose that $0\leq i<n$. The case $i+b<n$ was
considered in Theorem 2.1. So, we can suppose that $i+b\geq n$. Since $0<b<n$,
we have $0\leq i+b-n<i$, and therefore $x_{i+b}=x_{i+b-n}$. Consequently,

$\alpha x_{i}+\beta p_{i}p_{i+1}\ldots p_{i+b-1}x_{i+b}=\alpha x_{i}+\beta
p_{i}p_{i+1}\ldots p_{i+b-1}x_{i+b-n}.$

Multiply this number by $p_{i+b-n}p_{i+b-n+1}\ldots p_{i-1}$, which is
relatively prime to $q^{n}-p_{0}p_{1}\ldots p_{n-1}$, and therefore can not
change the property of being or not being an integer for the number $\alpha
x_{i}+\beta p_{i}p_{i+1}\ldots p_{i+b-1}x_{i+b}$, we obtain

$\beta p_{0}p_{1}\ldots p_{n-1}x_{i+b-n}+\alpha p_{i+b-n}p_{i+b-n+1}\ldots
p_{i-1}x_{i},$

which is an integer by Lemma 2.1 and Theorem 3.1. ∎

######  Remark 3.3.

The above results are trivially true for the cases when $b=0$ and $b=n$.
Indeed, for the case when $b=0$, if $(\alpha+\beta)U_{0}$ is an integer, then
$(q^{n}-p_{0}p_{1}\ldots p_{n-1})|(\alpha+\beta)$, and therefore
$(\alpha+\beta)x_{i}$ is also an integer for $i=0,1,\ldots,n-1$. For the case
when $b=n$, if $\alpha U_{0}+\beta U_{n}=\alpha U_{0}+\beta p_{0}p_{1}\ldots
p_{n-1}U_{0}+\beta$ is an integer, then $(\alpha+\beta p_{0}p_{1}\ldots
p_{n-1})U_{0}$ is also an integer ($\beta$ is an integer), and therefore
$(q^{n}-p_{0}p_{1}\ldots p_{n-1})|(\alpha+p_{0}p_{1}\ldots p_{n-1}\beta)$.
Consequently $(\alpha+\beta p_{0}p_{1}\ldots p_{n-1})x_{i}$ is an integer for
$i=0,1,\ldots,n-1$.

######  Corollary 3.4.

If one of the numbers $x_{j}$ $(j\in\\{0,1,\ldots,n-1\\})$ is a fraction with
denominator $d$ in its simplest form then all of $x_{i}$ for
$i=0,1,\ldots,n-1$, are like fractions with the same denominator $d$.

###### Proof.

It was mentioned earlier that one can always take $\alpha=1$,
$\beta=-q^{\phi(|q^{n}-p_{0}p_{1}\ldots p_{n-1}|)-1}$, and $b=1$. Then by the
main result, the number $x_{i}+\beta p_{i}x_{i+1}$ is an integer for
$i=0,1,\ldots,n-1$. Note that $d$ is a divisor of $q^{n}-p_{0}p_{1}\ldots
p_{n-1}$, and therefore $d$ is relatively prime to $\beta$ and all $p_{i}$ for
$i=0,1,\ldots,n-1$. Since $x_{j}$ is a fraction with denominator $d$, the
other numbers $x_{j-1}$, $x_{j-2}$, $x_{j-3},...$ are all like fractions with
the same denominator $d$ in their simplest forms. ∎

######  Corollary 3.5.

If one of the numbers $x_{j}$ $(j\in\\{0,1,\ldots,n-1\\})$ is an integer then
all of $x_{i}$ for $i=0,1,\ldots,n-1$, are integers.

##  4 Examples

Let us take $q=3$. Consider composition of functions $P=B_{0}\circ B_{1}\circ
B_{2}\circ B_{3}$, where $B_{0}(x)=\frac{-5x-2}{3}$,
$B_{1}(x)=\frac{2x+1}{3}$, $B_{2}(x)=\frac{7x+6}{3}$,
$B_{3}(x)=\frac{-x+3}{3}$. Here $n=4$, $p_{0}=-5$, $k_{0}=-2$, $p_{1}=2$,
$k_{1}=1$, $p_{2}=7$, $k_{2}=6$, $p_{3}=-1$, $k_{3}=3$, and
$q^{n}-p_{0}p_{1}p_{2}p_{3}=3^{4}-(-5)\cdot 2\cdot 7\cdot(-1)=11$.

The solution of equation $B_{0}\circ B_{1}\circ B_{2}\circ B_{3}(x)=x$ is the
number $x_{0}=-69/11$. Note that $x_{0}=x_{4}$. We can also find the other
numbers $x_{1}=x_{5}=37/11$, $x_{2}=50/11$, $x_{3}=12/11$, by solving the
equations $B_{1}\circ B_{2}\circ B_{3}\circ B_{0}(x)=x$, $B_{2}\circ
B_{3}\circ B_{0}\circ B_{1}(x)=x$, $B_{3}\circ B_{0}\circ B_{1}\circ
B_{2}(x)=x$, respectively. We also find the numbers $U_{i}=2^{i}/11$
($i=0,1,2,3,4$). Note that $4U_{0}+2U_{2}=2$ is an integer, which is
equivalent to say that $11|(4+2\cdot 3^{2})$. So, we can take $\alpha=4$,
$\beta=2$, and $b=2$. We observe that $4x_{i}+2p_{i}p_{i+1}x_{i+2}$ is an
integer for each of $i=0,1,2,4$. Indeed,

$4x_{0}+2p_{0}p_{1}x_{2}$ | $=4\cdot(-69/11)+2\cdot(-5)\cdot 2\cdot(50/11)$ | $=-116,$
---|---|---
$4x_{1}+2p_{1}p_{2}x_{3}$ | $=4\cdot(37/11)+2\cdot 2\cdot 7\cdot(12/11)$ | $=44,$
$4x_{2}+2p_{2}p_{3}x_{4}$ | $=4\cdot(50/11)+2\cdot 7\cdot(-1)\cdot(-69/11)$ | $=106,$
$4x_{3}+2p_{3}p_{4}x_{5}$ | $=4\cdot(12/11)+2\cdot(-1)\cdot(-5)\cdot(37/11)$ | $=38,$

are all integers.

Now, note that $-5U_{0}-13U_{1}=-4$ is also an integer, which is equivalent to
say that $11|-44=-5-13\cdot 3^{1}$. This means that we can take $\alpha=-5$,
$\beta=-13$, and $b=1$. So, $-5x_{i}-13p_{i}x_{i+1}$ should be an integer for
each of $i=0,1,2,3$. Indeed,

$-5x_{0}-13p_{0}x_{1}$ | $=-5\cdot(-69/11)-13\cdot(-5)\cdot(37/11)$ | $=250,$
---|---|---
$-5x_{1}-13p_{1}x_{2}$ | $=-5\cdot(37/11)-13\cdot 2\cdot(50/11)$ | $=-135,$
$-5x_{2}-13p_{2}x_{3}$ | $=-5\cdot(50/11)-13\cdot 7\cdot(12/11)$ | $=-122,$
$-5x_{3}-13p_{3}x_{4}$ | $=-5\cdot(12/11)-13\cdot(-1)\cdot(-69/11)$ | $=-87,$

are all integers. These observations are in perfect agreement with the main
results of the current paper.

##  5 Applications

If $p_{i}\in\\{1,p\\}$, where $p$ is a nonzero integer and $(p,q)=1$ then
there are two type of functions $S_{k}(x)=\frac{px+k}{q}$ and
$T_{k}(x)=\frac{x+k}{q}$. Let us denote by $m$ the number of $S$ functions in
$P$. Then $U_{i}=\frac{q^{i}}{q^{n}-p^{m}},$ for $i=0,1,\ldots,n$. Denote by
$\sigma(i,j)$ the number of $S$ functions in the fragment $B_{i}B_{i+1}\ldots
B_{j-1}$ of $P$. In particular, $\sigma(i,i)=0$, because it corresponds to the
empty fragment of $P$. Let $x_{i}$ be the solution of the equation $B_{i}\circ
B_{i+1}\circ B_{i+2}\circ\ldots\circ B_{i+n-1}(x)=x$, where all the indices
are taken modulo $n$. Take $\alpha=p^{l}$ for some non-negative integer $l$,
and $\beta=-1$. For this special case the main result of the current paper can
be written as $p^{l}x_{i}-p^{\sigma(i,i+b)}x_{i+b}\in\mathbb{Z}$. This can be
visualized by writing $x_{i}$ as $p$-adic numbers in a table and noting that
$p$-adic digits at the corresponding place values of $p^{l}x_{i}$ and
$p^{\sigma(i,i+b)}x_{i+b}$ are identical, except for finitely many digits at
lower place values. Let us demonstrate this on an example. Let $q=2$, $p=11$,
$P=B_{0}\circ B_{1}\circ B_{2}\circ B_{3}\circ B_{4}\circ B_{5}\circ B_{6}$,
where $B_{0}(x)=B_{1}(x)=B_{2}(x)=B_{3}(x)=B_{5}(x)=T_{0}(x)$,
$B_{4}(x)=S_{5}(x)$, and, $B_{6}(x)=S_{3}(x)$. Here $n=7$, $m=2$,
$q^{n}-p^{m}=2^{7}-11^{2}=7$ and $U_{i}=2^{i}/7$ ($i=0,1,2,\ldots$). Note that

$U_{0}-U_{3}=-1,\ \ 11U_{0}-U_{2}=1,\ \ 11^{2}U_{0}-U_{1}=17,\ \
11^{3}U_{0}-U_{0}=190,$

are all integers, which is equivalent to say that

$7|(1-2^{3}),\ \ 7|(11-2^{2}),\ \ 7|(11^{2}-2^{1}),\ \ 7|(11^{3}-2^{0}),$

respectively. By the main result of the current paper we can say that

$x_{i}-11^{\sigma(i,i+3)}x_{i+3},\ \ 11x_{i}-11^{\sigma(i,i+2)}x_{i+2},\ \
11^{2}x_{i}-11^{\sigma(i,i+1)}x_{i+1},\ \ 11^{3}x_{i}-x_{i},$

are also integers for $i=0,1,2,\ldots$. The functions $B_{i}$, the numbers
$x_{i}$, their 11-adic representations (the letter A means digit 10), and the
patterns formed by the digits can be seen in the following table:

$x_{0}=53/7$ | $=\ldots$ | 7 | 9 | 4 | 7 | 9 | 4 | 8 | 6 |  |  |  |  | $S_{3}=B_{6}$
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---
$x_{6}=302/7$ | $=\ldots$ | 9 | 4 | 7 | 9 | 4 | 7 | 9 | 8 | 7 |  |  |  | $T_{0}=B_{5}$
$x_{5}=151/7$ | $=\ldots$ | 4 | 7 | 9 | 4 | 7 | 9 | 4 | 9 | 9 |  |  |  | $S_{5}=B_{4}$
$x_{4}=848/7$ | $=\ldots$ | 7 | 9 | 4 | 7 | 9 | 4 | 7 | A | 4 | 8 |  |  | $T_{0}=B_{3}$
$x_{3}=424/7$ | $=\ldots$ | 9 | 4 | 7 | 9 | 4 | 7 | 9 | 5 | 2 | 4 |  |  | $T_{0}=B_{2}$
$x_{2}=212/7$ | $=\ldots$ | 4 | 7 | 9 | 4 | 7 | 9 | 4 | 8 | 1 | 2 |  |  | $T_{0}=B_{1}$
$x_{1}=106/7$ | $=\ldots$ | 7 | 9 | 4 | 7 | 9 | 4 | 7 | 9 | 6 | 1 |  |  | $T_{0}=B_{0}$
$x_{0}=53/7$ | $=\ldots$ | 9 | 4 | 7 | 9 | 4 | 7 | 9 | 4 | 8 | 6 |  |  |

More examples of such patterns and applications to the original $3x+1$ problem
are given in [2] and the references therein.

The Finite Cycles Conjecture mentioned at the beginning of this paper claims
that the only integer cycles for $3x+1$ problem are the ones generated by
$x_{0}=0$, $x_{0}=-1$, $x_{0}=1$, $x_{0}=-5$, and $x_{0}=-17$. These numbers
correspond to compositions $P_{1}=T$, $P_{2}=S$, $P_{3}=T\circ S$,
$P_{4}=T\circ S\circ S$, and $P_{5}=T\circ T\circ T\circ S\circ S\circ S\circ
T\circ S\circ S\circ S\circ S$. For these compositions $q=2$, $p_{i}=1$ or
$3$, and the numbers in the following table are calculated.

$x_{0}$ | $P$ | $n$ | $q^{n}-p_{0}p_{1}\ldots p_{n-1}$
---|---|---|---
$0$ | $P_{1}$ | 1 | $2^{1}-1=1$
$-1$ | $P_{2}$ | 1 | $2^{1}-3=-1$
$1$ | $P_{3}$ | 2 | $2^{2}-3=1$
$-5$ | $P_{4}$ | 3 | $2^{3}-3\cdot 3=-1$
$-17$ | $P_{5}$ | 11 | $2^{11}-3^{7}=-139$

These compositions with integer $x_{0}$ also show that the main results of the
current paper, namely Theorem 3.1 and Corollary 3.2 can not be written as
"iff" statements. Indeed, if the numbers $x_{i}$ $(i\in\\{0,1,\ldots,n-1\\})$
are integers then $\alpha x_{i}+\beta p_{i}p_{i+1}\ldots p_{i+b-1}x_{i+b}$ is
an integer for any choice of integers $\alpha,\beta,b$, which is not the case
for $\alpha U_{0}+\beta U_{b}$. Nevertheless, Lemma 2.1 can be written as an
"iff" statement.

The composition $P_{5}$ is different from the others in the sense that
$q^{n}-p_{0}p_{1}\ldots p_{n-1}\neq\pm 1$ but $x_{0}$ is still an integer. The
compositions $P_{i}^{k}$ ($i=1,2,\ldots,5$; $k=1,2,\ldots$), defined
recursively by $P_{i}^{1}=P_{i}$ and $P_{i}^{k+1}=P_{i}^{k}\circ P_{i}$
($i=1,2,\ldots,5$; $k=1,2,\ldots$) also have integer $x_{0}$. Determination of
all such compositions with integer $x_{0}$ or proving that all other
compositions correspond to non-integer rational $x_{0}$ will be helpful for
the solution of $3x+1$ problem. The results of the current paper shed some
light on the structure of rational cycles and therefore they might be useful
for this purpose.

##  6 Conclusion

In the paper some generalizations of Collatz conjecture or $3x+1$ problem are
studied. Some results are obtained proving that special linear combinations of
the terms of rational cycles are integers. Demonstrations of these results on
some concrete examples are given. These results are then used to explain some
patterns of digits in $p$-adic representations of the rational cycles.

### Acknowledgements

## References

  * [1]
  * [2] Y.N. Aliyev, The $3x+1$ Problem For Rational Numbers : Invariance of Periodic Sequences in $3x+1$ Problem, 2020 IEEE 14th International Conference on Application of Information and Communication Technologies (AICT), Tashkent, Uzbekistan, 2020, 1-4. https://doi.org/10.1109/AICT50176.2020.9368585
  * [3] J.H. Conway, On Unsettleable Arithmetical Problems, The American Mathematical Monthly 120(3), 192-198 (2013). https://doi.org/10.4169/amer.math.monthly.120.03.192
  * [4] C.J. Everett, Iteration of the number-theoretic function $f(2n)=n,f(2n+1)=3n+2$, Advances in Mathematics, 25(1) (1977) 42-45. https://doi.org/10.1016/0001-8708(77)90087-1
  * [5] I. Korec, A density estimate for the $3x+1$ problem, Math. Slovaca 44 (1994), no. 1, 85–89. https://eudml.org/doc/32414
  * [6] T. Laarhoven, B. de Weger, The Collatz conjecture and De Bruijn graphs, Indagationes Mathematicae, 24(4), 2013, 971-983. https://doi.org/10.1016/j.indag.2013.03.003
  * [7] J.C. Lagarias (Ed.), The Ultimate Challenge: The $3x+1$ Problem, 344 pp., AMS, 2010.
  * [8] J.C. Lagarias, The set of rational cycles for the $3x+1$ problem, Acta Arith. 56 (1990) 33–53. http://eudml.org/doc/206298
  * [9] E. Lehtonen, Two undecidable variants of Collatz’s problems, Theor. Comput. Sci., 407, 1–3 (2008) 596–600. https://doi.org/10.1016/j.tcs.2008.08.029
  * [10] D. Levy, Injectivity and surjectivity of Collatz functions, Discrete Mathematics, 285, Issues 1–3, 2004, 191-199. https://doi.org/10.1016/j.disc.2004.03.008
  * [11] K. Matthews, A. Watts, A generalization of Hasse’s generalization of the Syracuse algorithm, Acta Arithmetica 43.2 (1984) 167-175. http://eudml.org/doc/205897
  * [12] K. Matthews, A. Watts, A Markov approach to the generalized Syracuse algorithm, Acta Arithmetica 45.1 (1985): 29-42. http://eudml.org/doc/205953
  * [13] J.L. Rouet, M.R. Feix, A Generalization of the Collatz Problem. Building Cycles and a Stochastic Approach, Journal of Statistical Physics 107, 1283–1298 (2002). https://doi.org/10.1023/A:1015122111180
  * [14] T. Tao, Almost all orbits of the Collatz map attain almost bounded value, Forum of Mathematics, Pi, 10, E12, 2022. https://doi.org/10.1017/fmp.2022.8
  * [15] R. Terras, A stopping time problem on the positive integers, Acta Arith.30.3 (1976), 241. http://eudml.org/doc/205476


# SOLIS: Autonomous Solubility Screening using Deep Neural Networks
††thanks: This work was supported by the Leverhulme Trust through the
Leverhulme Research Centre for Functional Materials Design and the H2020 ERC
Synergy Grant Autonomous Discovery of Advanced Materials under grant agreement
no. 856405. Correspondence should be directed to
<EMAIL_ADDRESS>
Gabriella Pizzuto  Department of Chemistry, University of Liverpool, United
Kingdom  Jacopo de Berardinis  Department of Informatics, King’s College
London, United Kingdom  Louis Longley  Department of Chemistry, University of
Liverpool, United Kingdom  Hatem Fakhruldeen  Leverhulme Research Centre for
Functional Materials Design, University of Liverpool, United Kingdom  Andrew
I. Cooper  Department of Chemistry, University of Liverpool, United Kingdom
Leverhulme Research Centre for Functional Materials Design, University of
Liverpool, United Kingdom

###### Abstract

Accelerating material discovery has tremendous societal and industrial impact,
particularly for pharmaceuticals and clean energy production. Many
experimental instruments have some degree of automation, facilitating
continuous running and higher throughput. However, it is common that sample
preparation is still carried out manually. This can result in researchers
spending a significant amount of their time on repetitive tasks, which
introduces errors and can prohibit production of statistically relevant data.
Crystallisation experiments are common in many chemical fields, both for
purification and in polymorph screening experiments. The initial step often
involves a solubility screen of the molecule; that is, understanding whether
molecular compounds have dissolved in a particular solvent. This usually can
be time consuming and work intensive. Moreover, accurate knowledge of the
precise solubility limit of the molecule is often not required, and simply
measuring a threshold of solubility in each solvent would be sufficient. To
address this, we propose a novel cascaded deep model that is inspired by how a
human chemist would visually assess a sample to determine whether the solid
has completely dissolved in the solution. In this paper, we design, develop,
and evaluate the first fully autonomous solubility screening framework, which
leverages state-of-the-art methods for image segmentation and convolutional
neural networks for image classification. To realise that, we first create a
dataset comprising different molecules and solvents, which is collected in a
real-world chemistry laboratory. We then evaluated our method on the data
recorded through an eye-in-hand camera mounted on a seven degree-of-freedom
robotic manipulator, and show that our model can achieve 99.13% test accuracy
across various setups, while being simple and fast to train and, as a result,
easily transferable to a robotic platform.

###### Index Terms:

solubility screening, autonomous material discovery, laboratory automation

##  I Introduction

Laboratory robotics working collaboratively with researchers using modular and
reconfigurable lab equipment brings us closer to the vision of Industry 5.0
[1], which has the overarching goal of adopting human-centric approaches to
digital technologies. Self-driving laboratories that are capable of planning,
executing and evaluating the results of experiments have the potential to make
faster, more efficient progress in materials research. The automation specific
devices, such as sample changers, is generally confined to a specific task and
has a high cost of entry, making it inaccessible to many research groups.
Also, autonomous labs can only be realised through linking the abstract
concepts of chemical processes with the hardware responsible for its execution
[2], making use of algorithms to facilitate decision making for scientists.
These algorithms should be robust to real laboratory conditions and
importantly, generalisable to different experiments. This would allow
scientists to trust the predictions and in turn facilitate automation of
chemical discovery.

Currently, human scientists spend a significant amount of time on mundane
tasks because robotic platforms are only capable of performing a limited set
of operations in comparison to humans, and most of these require the scientist
to do all of the decision making, even for the most trivial scenarios [3].
Hence, a long-term goal is to reduce both human labour with respect to tedious
tasks and to have algorithms that carry out mundane resolutions, bringing us
closer to a digitalisation of chemical manufacturing which is a critical
technology path toward a sustainable society [4].

Figure 1: The robotic setup used for our model architecture to be trained on
and evaluated for solubility screening.

Material discovery laboratories rely on diverse experimental workflows that
are key to development of new, sustainable materials. An important process in
such laboratories is solubility screening. Solubility screens can be performed
in a variety of different ways, such as by titration, ultraviolet-spectroscopy
or high performance liquid chromatography (HPLC) methods [5], all of which
involve expensive, customised equipment.

A recent method has evaluated low-cost equipment (cameras) with a simple
computer vision method that measures turbidity [3], where the captured image
is analysed to measure the average brightness of the solute-solvent mixture.
The latter is then compared to a reference image of the solvent to conclude
whether the solution is fully dissolved or not. The brighter the image, the
more undissolved solute is present because the clear solvent appears black
when placed against a dark background. However, this approach is limited to
simple conditions within highly controlled environments that may not be easily
generalisable in complex real-world laboratory environments. Specifically,
initial researcher input is required every time the workflow runs, a black
background is needed, and the method does not tolerate different lighting
conditions.

Recently, vision methods based on deep neural networks and convolutional
neural networks (CNNs) have proven effective in laboratory and medical
environments [6]. Such methods have gained traction for visual recognition of
lab hardware [7] and have shown promising results in challenging real-world
settings. This is remarkable since laboratory glassware and many solvents are
in fact transparent and hence pose a significant challenge for classical
machine vision algorithms. Deep models pave the way for future laboratory
robotics in autonomous laboratories and will be key to overcome current
limitations in methods relying on more traditional approaches.

Inspired by how human chemists visually assess solutions, we propose a novel
optimal solubility screening method that does not require analytic
instruments, such as HPLC-based methods, and that by using deep models as
feature extractors can generalise to standard laboratory environments. Since
accurate knowledge of the precise solubility limit of a molecule is often not
required, it can be enough to measure a threshold of solubility in a solvent,
for example to facilitate solvent-antisolvent polymorph screening. As
illustrated in Fig. 1, our novel model architecture for solubility screening
is evaluated in a real-world laboratory with all the challenges associated
with this environment. The input images are recorded directly from the Franka
Emika Panda robot and these are used to train and evaluate our cascaded system
combining the state-of-the-art mask regional-CNN (R-CNN), which feeds the
features from the cropped region of interest (RoI) around the sample vial into
the classifier trained to determine whether the solute-solvent has fully
dissolved.

In summary, our main contributions in this paper include:

  1. 1.

the first cascaded deep model for solubility screening, relying on state-of-
the-art image segmentation and classifier;

  2. 2.

an end-to-end architecture without human intervention or pre-calibration;

  3. 3.

the evaluation of the proposed model on a novel dataset collected in a a
conventional, unmodified materials chemistry laboratory.

The rest of the paper is structured as follows: Section II gives a review of
related work on solubility screening, and Section III reports the novel
computational model for automated solubility screening and the dataset
collected. Section IV details the experimental evaluation, including the
results obtained. Finally, Section V draws conclusions and gives direction for
future work.

##  II Related Work

Solubility screening has been used in research and development laboratories
for a wide array of experimental studies and it remains a key process in many
workflows. In most cases, the process is a case of adding ‘excess solid’ or
‘excess solvent’. For the former, different analytical methods can be used
once the solution reaches equilibrium, such as titration techniques,
gravimetric methods, HPLC, or NMR [8]. Another method for the ‘excess solid’
process involves filtering, drying, and weighing the undissolved solid [8]. On
the other hand, for the ‘excess solvent’ process the common practice is to add
solvent to the solute until the latter fully dissolved in it [8]. In this
case, techniques such as Fourier transform infrared spectroscopy (FTIR) [9]
are used to measure the point at which the solute dissolves. While both
methods have their respective strengths and shortcomings, the ‘excess solvent’
process is best suited for integrating with robotic platforms as it requires
less specialist equipment, and does not require additional
filtration/separation steps. It also happens to be the one that in practice is
most widely used.

Given that the method of using a light source and a detector aims at
replicating the human eye, methods relying on measuring ‘turbidity’ or a
‘clear point’ were introduced. A recent work proposed a method to measure the
brightness of the solution of captured images [3]. While this method uses low-
cost and non-invasive equipment, the main challenge that it fails to address
is robustness to dynamic laboratory environments and being an end-to-end model
that does not require human input. Moreover, as the method relies on the user
selecting the RoI, an inexperienced scientist might not select the correct
area and in turn this would affect the results. In fact, while here the RoI is
is a manually determined slice of the vial that would vary depending on the
human chemist/experiment, our method uses the whole vial and hence is the same
every time.

Our proposed approach addresses this by adopting the mask R-CNN to focus on
the RoI, and hence remove the need of human input. Together with using deep
models for extracting features from the RoI, our method provides flexibility
in using the network on images without the need of adding a black background
or fixing the illumination. From the technical point of view, our approach is
novel because the proposed framework combines mask R-CNN with deep model
classifiers for decision-making within a laboratory environment, which is
still an open problem. Our architecture is able to deal with data coming from
images recorded either from the human chemist or the ‘robotic chemist’ due to
its ability in dealing with background and illumination variance.

Figure 2: Overview of the proposed approach for determining whether a molecule
dissolves in a given solvent. An input image recorded from the Franka Emika
Panda robot is fed into our cascaded system combining the mask R-CNN, which
feeds the features from the cropped RoI around the vial into the classifier
trained to determine whether the solute-solvent has fully dissolved.

##  III SOLIS Architecture

Our model here aims at addressing the challenge of determining which molecules
dissolve in a given solvent. Our primary goal and main interest is to obtain
an architecture that can be transferred to a robotic platform in a laboratory
environment that has been designed by and for human chemists. We present our
cascaded system combining a mask R-CNN, with a deep model classifier for
solubility screening. An overview of the proposed method is illustrated in
Fig. 2. Section III-A gives an in-depth overview of the cascaded model, while
further details about the dataset we collected to evaluate our models is given
in Section III-B.

###  III-A Cascaded Model

Our system relies on the state-of-the-art CNNs tailored for autonomous
solubility screening. The architecture can detect glass laboratory vials
captured directly from eye-in-hand camera mounted on the robotic manipulator,
without the need of a human scientist to highlight the RoI. Given the
coordinates of the bounding box predicted around the mask of the vial, we use
the features from this region to determine whether the molecule has dissolved.
As shown in Fig. 2, our architecture is mainly divided into three stages. In
the first part, we collect and preprocess the image. Following this, the data
is fed into an image segmentation model and a CNN-based image classifier. Our
system can be deployed on a robotic platform as a “robotic scientist” to help
lab-based human scientists to avoid repetitive tasks.

####  III-A1 Image Acquisition

For the image acquisition stage, we use a Franka Panda Emika robot equipped
with an Intel RealSense D435i camera on its end-effector; the latter attached
via a 3D printed mount. The robot perceived its workspace from above the IKA
Plate, as illustrated in Fig. 1. All images were recorded using a framerate of
30 frames per second. In contrast to [3], we opted to put the camera on the
end-effector of the robot rather than as its own station since this gives us
the possibility to use the camera also for scene understanding and visual
servoing. It also offers the future potential to mount this camera on a mobile
“robotic chemist” [10]. For these instances, the mask R-CNN from the first
stage could easily be adapted to realise such tasks.

####  III-A2 Image Segmentation

Following successful image acquisition, the next step is to extract meaningful
features for training the solubility classifier. As a result, the first stage
of the cascaded model is a mask R-CNN [11]. Given that the mask R-CNN is
currently the state-of-the-art method for instance segmentation with respect
to a large array of image segmentation tasks [11], this work adopts the pre-
trained model on the TransProteus dataset [7]. The role of this model is to
find the region and boundaries of the glass laboratory vial in the image that
is recorded directly from the end-effector of the Panda Emika Franka robotic
arm. As we are operating in a challenging environmental setup, it would be
ideal for the classifier to focus on the output segment that is contained
within the RoI mask. By taking this approach, there is still the important
feature of limiting the region of the image in which the output segment is
found while removing the need for human input in selecting the RoI as in
related works [3]. Specifically for our task, the mask R-CNN uses ResNet as a
backbone [12]. From the candidates proposed by the mask R-CNN, the instance
bounding box and class are predicted by the box head. Here, we use
$vial\\_bbox$, which are the coordinates of the bounding box of the glass vial
that are predicted by the mask R-CNN. For the next step, the underlying image
features obtained from the bounding box are used as input to the CNN. While
masks are also provided by the network for both the vessel (glass vial) and
the material (solution), our work only requires the image segmentation to
produce a RoI with useful features of the complete vial for the classifier;
hence, we only used the bounding box around the segmented mask.

####  III-A3 Solubility Classification

Here, we view the problem of solubility screening as a classification process,
since for many material discovery experiments, we only need a threshold of
solubility. Convolutional neural networks are the current state-of-the-art
methods for image classification. Therefore, given the challenging data that
is involved in our application, this framework was the obvious choice. Given
the small dataset, we opted to use transfer learning from a pre-trained
network on the Imagenet dataset. For a thorough evaluation, we looked into
both finetuning and feature extraction strategies for training. For
finetuning, the pretrained model is retrained and hence all of its parameters,
that is, the higher order features and the classifier are updated with respect
to the current task of classifying whether the dissolved state has been
reached. For feature extraction, the model is frozen such that only the output
layer was updated; hence, the feature maps from the pre-trained model are
retained. Another important criteria we considered was image augmentation,
where we adopted random resizing, cropping, and horizontal flipping for the
training set and resizing and center cropping for the validation set. All the
data was normalised prior to feeding it into the model architecture. We have
tested our cascaded model using different CNN architectures, and we give
further detail of these setups in the experimental evaluation (Section IV).

###  III-B Dataset

Figure 3: An overview of the recorded dataset for benzimidazole with acetone.
(a) Initial solid dispensing is carried out, where benzimidazole is weighed
into the glass vial. A magnetic stir bar is also placed in the vial. (b)
Initial solvent dispensing occurs, and the stirring plate is also switched on
to mix the solution. (c) Mixing of the solution continues and the mixture
becomes partially dissolved. (d) The solution becomes fully dissolved.

Given the novelty of this application, there exists no dataset that is
publicly available for automated solubility screening. Moreover, the images
recorded in a relevant study [3] were recorded using a black background and a
reference vial. While we too considered a setup with a black background, we
specifically focused on having non-sterile conditions to facilitate easy
deployment in a real laboratory environment. Hence, it was crucial for us to
build a specific dataset for this problem. We collected data directly from a
materials discovery laboratory where a human chemist manually carried out the
‘excess solvent’ procedure. After preparing the solute in 20 mL vials with a
magnetic stir bar inside, the chemist pipetted the solvent into the sample
using increments of initially 1 mL but then 0.5 mL when it was judged by eye
that sample turbidity was dropping. During the experiment, the sample was
continuously stirred with a magnetic stirrer111IKA Plate (RCT digital) was
used for our experiments. to prevent the powder dropping to the bottom of the
vial. Stirring also prevented sample clumping The full dataset comprises 2
solids: caffeine and benzimidazole and 3 solvents: water, ethanol and acetone;
hence six different setups were used. The choice of these chemicals pivoted
around the fact that these can be used safely on an open-bench and could be
integrated in a crystallisation workflow. The experiments were supervised by a
human chemist, who was also responsible for annotating the dataset. Hence, the
dataset collection process was terminated once the human chemist visually
decided that the solute had dissolved. As mentioned, the RGB images were
recorded using an eye-in-hand camera on a Panda Emika robotic manipulator.

##  IV Experimental Evaluation

We posed two questions to evaluate this method. First, how does SOLIS (our
proposed approach) perform on the autonomous solubility screening task? More
precisely, this first experiment aims to compare the performance of different
state-of-the-art artificial neural networks to implement the CNN
classifier222We generally refer to this part of the pipeline as the CNN
classifier, regardless of the specific architecture implementing it. While the
previous stage, the mask R-CNN, is also a CNN classifier, we usually refer to
this as the image segmentation part. in our pipeline (c.f. Fig. 2). The
criteria for this comparison are the average classification performance on the
solubility screening task and the complexity of the classifier (number of
parameters). These are both significant for the selection of the specific
classification model that will be deployed in SOLIS, as we aim at minimising
inference times and the risk of failure. Second, given that solubility
screening is usually integrated within a larger workflow, when does it fail?
From the first experiment, we identify the most accurate and parsimonious
model to deploy in SOLIS, which is then used in the second experiment. This
comparison aims to illustrate that, although SOLIS can produce
misclassifications, our understanding of when this happens will equip us with
the right approach to compensate for these predictions and, as a result avoid,
interrupting expensive and time-consuming experiments.

###  IV-A Experimental Setup

For both the experiments, we evaluated all models on the solubility screening
task, quantifying their ability to predict the state of the solution given a
raw image recorded directly from the robotic platform. This was done by
computing the cross-entropy loss and accuracy of each model on the test set.
Cross entropy (CE) is the dominant measure for training and evaluating
classification models. Intuitively, it measures the difference between two
probability distributions: the target distribution (the actual solubility
status), and the predicted one (the model’s predictions of the solubility
status). Accuracy is derived from the model’s predictions by considering the
class associated with the maximum activation for a given input.

The framework is powered by PyTorch, using a machine equipped with an AMD
Ryzen Threadripper 3970X 32 Core CPU and a single NVIDIA GeForce RTX 3090
Graphical Processing Unit (GPU).

TABLE I: Summary of the CNN architectures tested as classifiers for solubility screening. The number of trainable parameters makes distinction between the backbone of each architecture and the last fully-connected layer used for classification (in brackets). CNN model | Input size | No. of parameters
---|---|---
VGG | 224x224 | 128 M (8.2K)
ResNet18 | 224x224 | 11.2 M (1K)
InceptionV3 | 224x224 | 24.3 M (5.6K)
Densenet | 299x299 | 7 M (2K)

###  IV-B Solubility Screening Classification

In this section, we explore different networks for the solubility classifier.
The architectures that were used include VGG [13], ResNet [12], Inception [14]
and Densenet [15]. The input size expected by each architecture and the
overall complexity are reported in Table I. Although there are different
models that could have been applicable, our decision was also motivated by
having an architecture that has proved to work in different domains and hence
using the reference fields to address our challenging task was the most
sensible approach. Moreover, given that the method presented in [3] requires
human input, we could not directly compare this approach to ours since our
proposed does not require human interventions.

As the RoI from the output of the mask R-CNN varied, we scaled the features to
the input size each architecture expects (c.f. Table I). Nonetheless, a
different transformation strategy was followed for the training sets in order
to add stochasticity during training. In particular, for each epoch, we select
a random region from the resulting segmentation – matching the input size of
the classifier. Although this approach may result in some samples where the
solubility status is hard to estimate/detect (e.g., the liquid appears only in
a relatively small area), we expect the regularising effect to improve the
generalisation capabilities of the resulting models.

All CNN models had been pre-trained on Imagenet [16], then considered for
transfer learning on the solubility screening task. To achieve this, the last
fully-connected layer of each architecture was discarded (as it was originally
designed for ImageNet’s classes) and replaced with a linear fully-connected
layer with 2 output units - one per solubility outcome (undissolved, and
dissolved). Each CNN classifier was then trained to minimise the CE loss
according to two different strategies: feature extraction, where the gradient
is back-propagated only up to the last (new) fully-connected layer - meaning
that the parameters of the rest of the network (the backbone) are left
unchanged, and fine-tuning, where the gradient is fully back-propagated and
all the parameters of the network are adjusted accordingly.

We used the Adam optimiser with a learning rate equal to 0.001, together with
an early stopping policy (with 10 epochs of patience, and a delta of
$10^{-2}$) to prevent over-fitting. Considering the limited size of the
dataset, we performed 5-fold cross validation. This allowed to evaluate each
architecture without exposing our results to the potential bias of a single
(fixed) split. The dataset was split into 5 partitions (folds) and each of
them was used as the test set of a distinct experiment (for the same
architecture). The predictions of models on each test fold are then
aggregated, so that test losses and accuracy can be computed. Although the
splits were produced apriori for all the experiments, we also set a random
seed (1992) for a number of Python modules (numpy, pandas, torch, math) to
ensure reproducibility.

The results are presented in Table II. All of the models give an overall good
performance that would be suitable to deploy in laboratory settings. However,
the best performing models are the ones that are finetuned. When looking at
the test loss, ResNet18 achieves the best performance in terms of average CE
loss and standard deviation. We believe that although transfer learning aided
our training methods as the dataset was relatively small (as the task was
quite unique and the images are very different to the ones found in Imagenet)
training only the last fully connected layer was not enough in all models.
Hence, the finetuning approach, while requiring higher training times in
comparison, achieved the best performance.

TABLE II: Overview of the performance of our proposed architecture with the CNN classifiers. For each model, we use finetuning for training the network and this is represented as FT, and we also used feature extraction during training which is illustrated through FE. Model | Method | Test Accuracy $\uparrow$ | Cross Entropy Loss $\downarrow$
---|---|---|---
VGG | FT | 0.9846 | 0.0694 ± 0.6645
VGG | FE | 0.9675 | 0.0872 ± 0.3853
ResNet18 | FT | 0.9904 | 0.0264 ± 0.2486
ResNet18 | FE | 0.9483 | 0.1274 ± 0.4345
InceptionV3 | FT | 0.9913 | 0.0295 ± 0.3989
InceptionV3 | FE | 0.9334 | 0.1759 ± 0.4137
Densenet | FT | 0.9911 | 0.0363 ± 0.5648
Densenet | FE | 0.9528 | 0.1188 ± 0.4177

###  IV-C Robustness Evaluation

For our application, solubility screening requires a high level of accuracy to
have a successful completion of the laboratory experiment. For example,
solubility screening is a precursor for widely used pharmaceutical
applications, such as crystallisation. Hence, should the solubility screen
produce incorrect results this could be costly. As a result, understanding the
cases where it could fail is important to ensure workflow success and increase
robustness of the overall system.

Since we opted for ResNet18 as the final model for deployment, we analysed in-
depth the misclassifications and these are illustrated in Fig. 4. Fig. 4
illustrates the five worst performances for each cross-validation fold. At a
closer look, we believe that a large proportion of the errors, for example
Fig. 4.b (1) and (V), result from reflections or spots on the glass vial
walls. In such cases, even a human chemist might misclassify the solution and
to mitigate these errors they would take another look at the sample,
potentially in different settings or under different light conditions.
Motivated by this, our future plan is that when the model predicts the
dissolved state, we create a buffer of images rather than use a single one to
base our decision on.

Another issue is that at the early stages of the experiment, the vial does not
have a large proportion of solution and if the majority of the crop takes into
account the top part of the transparent glass vial, it would be easy for the
network to confuse this with a clear dissolved liquid. Future work could also
explore the use of a third classification category indicating cases that are
borderline between dissolved and undissolved. This could be used as a basis
for subsequent experiments with the same solvent and solute molecule aiming at
more accurate determination of solubility limits.

To mitigate the challenge this brings (as illustrated in some cases in Fig.
4), in future work we plan to explore how different image crops can affect the
overall performance of our model.

Figure 4: An in-depth analysis of the misclassifications across the five
different folds for the ResNet18 model, which is finetuned. For each fold, we
report the top five worst predictions.

###  IV-D Discussion

All of our experiments were carried out using the recorded dataset (Section
III-B). While we emphasise that the dataset was carried out in real-world
laboratory conditions, in our preliminary experiments we explored the option
of including the same black background, as in [3], and also went further to
explore different lighting conditions. In our case, as the mask R-CNN in the
first step puts attention on the glass vial through the RoI selected within
the bounding box, this step made our model robust to different environmental
conditions and, as a result, we did not observe a better performance when
using a black background. The use of the mask R-CNN in the first stage also
gives the possibility of having the vial in different positions, not
necessarily within the centre of the image, and the proposed pipeline would
still work. Again, this is quite advantageous as it never is the case that
human or robotic scientists place glassware exactly in the same location.

Our long-standing goal is to have an end-to-end model that is not only
accurate, but is also computationally efficient for deployment on a real
robot. As such, we remained faithful towards models that have a smaller number
of parameters (such as opting for a ResNet18 against a ResNet152). From the
results obtained in Table II, it is evident that the difference between the
finetuned models ResNet18, InceptionV3 and Densenet is very small. In fact, it
is our understanding that when deployed on a robotic platform all of these
models would perform similarly. Therefore, knowing that InceptionV3 has in
reality significantly more parameters than ResNet18 or Densenet, our model
choice for deployment is either of the latter as this would result in faster
training and inference times.

In our results, we illustrated how this cascaded model can be used for
solubility screening; however, different processes in material chemistry also
rely on visual observation of solutions for multiphase materials systems.
Given the robustness and modularity of our method, we believe our model has
huge potential to change visual assessment of samples in pharmaceutical and
clean energy applications.

##  V Conclusion

In the fields of materials discovery and pharmaceutical chemistry, the
automatic screening of solutions is an important task and a staple procedure
in various workflows. We addressed this challenge of solubility screening in a
real-world laboratory environment by introducing and evaluating SOLIS, a novel
end-to-end cascaded model capable of determining whether a solute has
dissolved in a particular solvent. Specifically, we achieved this by combining
image segmentation using a mask R-CNN for extracting the RoI with CNNs as
image classifiers. We demonstrated how our method possesses key advantages
over existent ones, mainly in terms of not requiring human input and also by
evaluating it on a dataset collected in an unmodified laboratory, without any
inconvenient modifications to the laboratory. Our results indicate that the
proposed methodology can achieve excellent performance on this challenging
dataset, thus illustrating its potential to be used in autonomous laboratory
experiments.

As the model has already been evaluated on images captured directly from the
robotic manipulator, the evident next step is to deploy this model on the
robot in a closed-loop material discovery workflow. The laboratory instruments
and the robotic platform are already integrated in a robotic architecture for
chemistry laboratories [17] and hence this should be straightforward. Other
interesting research directions would be to look in more depth at using the
system with organic solvents: this would require the workflow to be carried
out in a fumehood. We believe our approach of training and evaluating the
proposed model in real-world laboratory conditions will ease this.

In the future, additional refinements to the method could facilitate a more
accurate estimation of the solubility, while still building on the same
classifier outlined here; for example, in a titration like method involving
several repeat measurements with the same solvent and solute combination.
These repeats would employ smaller increments of solvent addition, longer
equilibration times between liquid addition steps, and larger image buffers.
This may offer the potential for accurate measure of solubility limits, above
and beyond the solubility thresholds that have been investigated here.

## References

  * [1] K. A. Demir, G. Döven, and B. Sezen, “Industry 5.0 and human-robot co-working,” Procedia Computer Science, vol. 158, pp. 688–695, 2019.  3rd World Conference on Technology, Innovation and Entrepreunership.
  * [2] J. Bai, L. Cao, S. Mosbach, J. Akroyd, A. A. Lapkin, and M. Kraft, “From platform to knowledge graph: Evolution of laboratory automation,” JACS Au, vol. 0, no. 0, p. null, 0.
  * [3] P. Shiri, V. Lai, T. Zepel, D. Griffin, J. Reifman, S. Clark, S. Grunert, L. P. Yunker, S. Steiner, H. Situ, F. Yang, P. L. Prieto, and J. E. Hein, “Automated solubility screening platform using computer vision,” iScience, vol. 24, no. 3, p. 102176, 2021.
  * [4] O. Inderwildi, C. Zhang, X. Wang, and M. Kraft, “The impact of intelligent cyber-physical systems on the decarbonization of energy,” Energy Environ. Sci., vol. 13, pp. 744–771, 2020.
  * [5] A. Petereit and C. Saal, “What is the solubility of my compound? assessing solubility for pharmaceutical research and development compounds,” American Pharmaceutical Review, vol. 14, pp. 68–73, 07 2011.
  * [6] S. Eppel, H. Xu, M. Bismuth, and A. Aspuru-Guzik, “Computer vision for recognition of materials and vessels in chemistry lab settings and the vector-labpics data set,” ACS Central Science, vol. 6, no. 10, pp. 1743–1752, 2020.  PMID: 33145411.
  * [7] S. Eppel, H. Xu, Y. R. Wang, and A. Aspuru-Guzik, “Predicting 3d shapes, masks, and properties of materials inside transparent containers, using the transproteus cgi dataset,” Digital Discovery, pp. –, 2022.
  * [8] S. Black, L. Dang, C. Liu, and H. Wei, “On the measurement of solubility,” Organic Process Research & Development, vol. 17, no. 3, pp. 486–492, 2013\.
  * [9] A. d. P. M. da Silva and J. F. Cajaiba da Silva, “Determination of the adipic acid solubility curve in acetone by using atr-ftir and heat flow calorimetry,” Organic Process Research & Development, vol. 15, no. 4, pp. 893–897, 2011.
  * [10] B. Burger, P. M. Maffettone, V. V. Gusev, C. M. Aitchison, Y. Bai, X. Wang, X. Li, B. M. Alston, B. Li, R. Clowes, N. Rankin, B. Harris, R. S. Sprick, and A. I. Cooper, “A mobile robotic chemist,” Nature, vol. 583, pp. 237–241, Jul 2020.
  * [11] K. He, G. Gkioxari, P. Dollár, and R. B. Girshick, “Mask r-cnn,” 2017 IEEE International Conference on Computer Vision (ICCV), pp. 2980–2988, 2017.
  * [12] K. He, X. Zhang, S. Ren, and J. Sun, “Deep residual learning for image recognition,” 2016 IEEE Conference on Computer Vision and Pattern Recognition (CVPR), pp. 770–778, 2016.
  * [13] K. Simonyan and A. Zisserman, “Very deep convolutional networks for large-scale image recognition,” in International Conference on Learning Representations, 2015.
  * [14] C. Szegedy, V. Vanhoucke, S. Ioffe, J. Shlens, and Z. Wojna, “Rethinking the inception architecture for computer vision,” in 2016 IEEE Conference on Computer Vision and Pattern Recognition (CVPR), pp. 2818–2826, 2016.
  * [15] G. Huang, Z. Liu, L. Van Der Maaten, and K. Q. Weinberger, “Densely connected convolutional networks,” in 2017 IEEE Conference on Computer Vision and Pattern Recognition (CVPR), pp. 2261–2269, 2017.
  * [16] J. Deng, W. Dong, R. Socher, L.-J. Li, K. Li, and L. Fei-Fei, “Imagenet: A large-scale hierarchical image database,” in 2009 IEEE conference on computer vision and pattern recognition, pp. 248–255, Ieee, 2009.
  * [17] H. Fakhruldeen, G. Pizzuto, J. Glawucki, and A. I. Cooper, “Archemist: Autonomous robotic chemistry system architecture,” IEEE International Conference on Robotics and Automation, 2022.