# Scheduling and control of high throughput screening systems with uncertainties and disturbances

Adetola Oke, Laurent Hardouin, Xin Chen, Ying Shang

► To cite this version:

Adetola Oke, Laurent Hardouin, Xin Chen, Ying Shang. Scheduling and control of high throughput screening systems with uncertainties and disturbances. Production & Manufacturing Research: An Open Access Journal, 2022, 10 (1), pp.450-469. 10.1080/21693277.2022.2091679. hal-03794671

HAL Id: hal-03794671

https://univ-angers.hal.science/hal-03794671

Submitted on 13 Feb 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Scheduling and control of high throughput screening systems with uncertainties and disturbances

Adetola Oke, Laurent Hardouin, Xin Chen & Ying Shang

**To cite this article:** Adetola Oke, Laurent Hardouin, Xin Chen & Ying Shang (2022) Scheduling and control of high throughput screening systems with uncertainties and disturbances, Production & Manufacturing Research, 10:1, 450-469, DOI: 10.1080/21693277.2022.2091679

**To link to this article:** https://doi.org/10.1080/21693277.2022.2091679

<figure><img src="image_3.png" /></figure>

© 2022 The Author(s). Published by Informa UK Limited, trading as Taylor & Francis Group.

<figure><img src="image_4.png" /></figure>

Published online: 23 Jun 2022.

<figure><img src="image_5.png" /></figure>

Submit your article to this journal

<figure><img src="image_6.png" /></figure>

Article views: 316

<figure><img src="image_7.png" /></figure>

View related articles

<figure><img src="image_8.png" /></figure>

View Crossmark data

CrossMark
OPEN ACCESS

<figure><img src="image_9.png" /></figure>

Check for updates

# Scheduling and control of high throughput screening systems with uncertainties and disturbances

Adetola Okeª, Laurent Hardouinᵇ, Xin Chenᶜ and Ying Shangᵈ

ªDepartment of Electrical and Computer Engineering, Southern Illinois University Edwardsville, Edwardsville, Illinois, USA; ᵇLaboratoire Angevin de Recherche en Ingénierie des Systèmes (LARIS), Université d'Angers, France; ᶜDepartment of Industrial Engineering, Southern Illinois University Edwardsville, Edwardsville, Illinois, USA; ᵈDepartment of Electrical Engineering, Indiana Institute of Technology, Fort Wayne, Indiana, USA

## ABSTRACT

Based on the existing research work on the modeling and control of cyclically operated high-throughput screening systems, this paper presents the scheduling and control of high-throughput screening systems with uncertainties and disturbances. Different definitions of disturbance decoupling problems are considered in order to achieve the optimal scheduling of high-throughput screening systems. Online scheduling was achieved by computing optimal pre-filters and output feedback controllers, which can be computed using residue theory. The control strategies achieved by the controllers are optimal in the sense of the just-in-time criterion commonly used in industrial schedule practice.

## ARTICLE HISTORY

Received 4 March 2021
Accepted 12 June 2022

## KEYWORDS

High-throughput screening;
disturbance decoupling;
uncertain max-plus linear
systems; scheduling

## 1. Introduction

A successful drug discovery is an extremely time-consuming procedure, including initial target identification and validation, pre-clinical trials on animals, regulatory approval to start trials in humans, clinical trials, submission of marketing and manufacturing authorization, licensing review, product sale, and post-marketing surveillance (Major, 1998; Mayer et al., 2008; Mayr & Fuerst, 2008; Noah, 2010; Pereira & Williams, 2007). With the development of robotics and high-speed computing technology, it is feasible to develop automatic systems that can screen a large number of biochemical compounds in a short period. Such an automatic compound screening and analyzing process is called high-throughput screening (HTS) in drug discovery of pharmaceutical industries. HTS is a standard technology routinely employed in the pharmaceutical industry for drug discovery processes. It is used for initial screening in the process of drug discovery to reduce what is an almost infinite number of possible combinations of compounds to a reasonably few enough possibilities on which further testing can be carried out. Many HTS systems consist of several activities on several different resources. An HTS operation can incorporate multiple batches, each with hundreds of events, where a batch is a combination of all the operations to be performed on a set of substances for complete

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>
<figure><img src="image_12.png" /></figure>

analysis. HTS provides a practical and efficient method to test a large number of synthetic compounds in miniaturized in vitro assays to identify hit targets of interest. Then, the chemical compounds that have therapeutic and useful pharmacological or biological activities, called leads, are evaluated and undergo lead optimization to identify promising lead compounds. Followed by the initial synthesis and animal testing in preclinical trials and three phases of clinical trials on humans, a drug can be put on the market after the Food and Drug Administration (FDA) approval.

It is necessary to develop a schedule for an HTS system that ensures that analysis is completed in the shortest possible time, i.e. increase the throughput or reduce the makespan (makespan is the completion time of the last activity of a process). The optimal schedule must also ensure consistency in time spent on each activity for every batch. This means that the operation of the HTS system must be cyclic. There have been various methods developed to obtain an optimal schedule for an HTS system and other cyclic systems. The problem of obtaining a time-optimal schedule for an HTS system, as an example of a strictly cyclic system, was formulated and solved as mixed integer linear programming (MILP) problem by Mayer and Raisch (2004). This MILP approach is useful for obtaining a static schedule, which works well so long as there are no changes to the cyclic design. A time-optimal schedule was also obtained by Mayer and Raisch (2003) for cyclic systems with multicapacity resources. HTS processes with hierarchical nested cycles were considered by Mayer et al. (2008), where mixed integer optimization was adopted to minimize the overall cycle time. In these existing works, the schedules obtained are static as they do not tolerate any deviation from the cyclic processes.

Starting from a static schedule, a max-plus algebra model of an HTS system was obtained in Brunsch and Raisch (2009). In a max-plus linear system, the addition and multiplication operations in the classical control systems are replaced by max and addition, respectively. The real number is also replaced by a special algebra called max-plus algebra with the operations changed accordingly (Baccelli et al., 1992; Golan, 1999). The scheduling problem can be transformed into the problem of computing a feedback controller that generates possible control actions in the presence of disturbances for such a max-plus linear system model. The controller generates a control input that gives an optimal schedule in spite of system malfunctions and delays. Brunsch et al. (2012) presented an optimal feedback controller. There have been investigations in this vein into scheduling and control of dioid models of discrete event systems in the presence of disturbances (Shang et al. (2013, p. 2016a)), all of which can be applied to HTS systems. In Shang, et al. (2016b), the concept of disturbance decoupling problem (DDP) for uncertain max-plus linear systems whose system matrices vary between intervals was introduced.

In this paper, an HTS process is first modeled in max-plus algebra and a static schedule is developed for the process based on the model. Online scheduling was then considered by computing controllers that ensure an optimal schedule when there are disturbances. The result is extended to the uncertain case where system matrices lie in intervals. The rest of this paper is organized as follows. Section 2 gives the mathematical fundamentals of the max-plus algebra that is to be used in describing the HTS system. In Section 3, the specifications of the HTS system considered are described. Sections 4 and 5 detail online scheduling of certain and uncertain HTS systems, respectively. Simulation
<figure><img src="image_13.png" /></figure>

results are presented in Section 6 and a discussion of the result and contribution of the study is presented in Section 7. Section 8 concludes this paper with future research directions.

## 2. Mathematical preliminaries of max-plus algebra

### 2.1. Idempotent semirings

**Definition 1.** A semiring is a set $S$, equipped with two operations $\oplus$ and $\otimes$, such that $(S, \oplus)$ is a commutative monoid (the zero element will be denoted $\varepsilon$), $(S, \otimes)$ is a monoid (the unit element will be denoted $e$), operation $\otimes$ is right and left distributive over $\oplus$, and $\varepsilon$ is absorbing for the product (i.e. $\varepsilon \otimes a = a \otimes \varepsilon = \varepsilon, \forall a$).

A semiring $S$ is *idempotent* if $a \oplus a = a$ for all $a \in S$. In an idempotent semiring $S$, operation $\oplus$ induces a partial order relation:

$$ab \succeq \Leftrightarrow a = a \oplus b, \forall a, b \in S. \quad (1)$$

Then, $a \vee b = a \oplus b$. An idempotent semiring $S$ is *complete* if sums of infinite numbers of terms are always defined, and if multiplication distributes over infinite sums too. In particular, the sum of all the elements of the idempotent semiring is denoted $T$ (for ‘top’). In this paper, we denote $\overline{Z}_{\max} = (\mathbb{Z} \cup \{-\infty, +\infty\}, \max, -\infty, +, 0)$ as the integer max-plus semiring, where $\varepsilon = -\infty$ is the neutral (zero) element for $\max$ and $e = 0$ is the unit element for $+$. A non-empty subset $\mathcal{B}$ of a semiring $S$ is a sub-semiring of $S$ if for all $a, b \in \mathcal{B}$ we have $a \oplus b \in \mathcal{B}$ and $a \otimes b \in \mathcal{B}$.

**Definition 2.** A mapping $f: S \to S$, where $S$ is a complete idempotent semiring, is *residuated* if and only if $f(\varepsilon) = \varepsilon$ and $f$ is lower-semicontinuous, that is,

$$f\left(\bigoplus_{i \in I} a_i\right) = \bigoplus_{i \in I} f(a_i),$$

for any (finite or infinite) set $I$. The mapping $f$ is said to be *residuated* and $f^{\#}$ is called its *residual*.

**Theorem 1.** (Baccelli et al. 1992) When $f$ is residuated, $f^{\#}$ is the unique order preserving mapping such that

$$f \circ f^{\#} \preceq \mathrm{Id} \quad \text{and} \quad f^{\#} \circ f \succeq \mathrm{Id}, \qquad (2)$$

where $\mathrm{Id}$ is the identity mapping from $S$ to $S$.

The operator $f \circ f^{\#}$ denotes function composition, i.e. $f \circ f^{\#}$ means that $f$ is a function of $f^{\#}$.

It is straightforward that: $L_a : S \to S, x \mapsto ax$ and $R_a : S \to S, x \mapsto xa$ are lower semi-continuous. Therefore, these mappings are residuated i.e. $L_a(x) \preceq b$ (resp. $R_a(x) \preceq b$) admits a greatest solution, then the following notations are considered:
<figure><img src="image_14.png" /></figure>

$$L_a^\#(b) = a \circled{\vert} b = \oplus \{x \vert axb\} \text{ and}$$

$$R_a^\#(b) = b \circled{\vert} a = \oplus \{x \vert xab\}, \forall a, b \in S,$$

where $L_a^\#$ and $R_a^\#(b)$ are the residual mappings.

**Theorem 2.** (Baccelli et al. 1992) Over a complete max-plus algebra, the implicit equation $x = ax \oplus b$ admits $x = a^*b$ as the least solution, where $a^* = \oplus_{i \in \mathbb{N}} a^i ((\cdot)^*)$ is the Kleene Star operator) with $a^0 = e$.

**Definition 3** (Cohen et al. 1996, 1997, 2006). Let $S$ be a complete idempotent semiring and let $C$ be a $n \times p$ matrix with entries in $S$. We call null kernel of $C$ as the set of elements $x \in S^p$ such that $Cx = \varepsilon$, denoted as $ker C$. We call equivalence kernel of $L_C$ (denoted by $ker_{eq} C$), the subset of all pairs of elements of $S^p$ whose components are both mapped by $L_C$ to the same element in $S^n$, i.e. the following definition

$$ker_{eq} C := \{ (s, s') \in (S^p)^2 | Cs = Cs' \}. \quad (3)$$

Clearly $ker_{eq} C$, is an equivalence relation on $S^p$, i.e., $Cs = Cs' \Leftrightarrow s' \equiv s(\text{mod}ker_{eq}C)$ and furthermore it is a congruence and then we can define the quotient $S^p/ker_{eq}C$.

The subset of elements $s' \in S^p$ that are equivalent to $s$ modulo $ker_{eq} C$ is denoted as $[s]_C$, i.e.,

$$s_C = \{s' \in S^p | s' \equiv s(\text{mod } ker_{eq}C)\} \subset S^p. \quad (4)$$

**Definition 4** (Restricted map). Let $f: S^p \to S^n$ be a map and $A \subseteq S^p$. We will denote $f|_A: A \to S^n$ the map defined by $f|_A = f \circ \text{Id}|_A$ where $\text{Id}|_A: A \to S^p, x \mapsto x$ be the canonical injection. Identically, let $B \subseteq S^n$ with $\text{Imf} \subseteq B$. Map $\text{Bf}: S^p \to B$ is defined by $f = \text{Id}|_B \circ \text{Bf}$, where $\text{Id}|_B: B \to S^n, x \mapsto x$ be the canonical injection.

**Definition 5** (Isotone map). A map $f: S^p \to S^p$ is said to be order preserving or isotone if the following property holds: $a \preceq b \Rightarrow f(a) \preceq f(b)$.

**Definition 6** (Closure map). An isotone map $f: S^p \to S^p$ is a closure map iff $f \succeq \text{Id}_{S^p}$ and $f \circ f = f$.

**Proposition 1** (2003. A closure map $f: S^p \to S^p$ restricted to its image $\text{Imf}|f$ is a residuated map whose residual is the canonical injection $\text{Id}|_{\text{Imf}}: \text{Imf} \to S^p, s \mapsto s$.

**Corollary 1.** Let $K: S^p \to S^p, s \mapsto s^*$ be a map, where $s^* = \oplus_{i \in \mathbb{N}} s^i$. The map $\text{Imf}|_K K$ is a residuated map whose residual is $(\text{Imf}|_K K)^\# = \text{Id}|_{\text{Imf}}$. This means that $x = s^*$ is the greatest solution to inequality $x^* \preceq s^*$. Actually, the greatest solution achieves equality.

There are basic properties for star and residuation operations in the residuation theory (Baccelli et al. 1992), for example,

$$a^*(ba^*)^* = (a \oplus b)^* = (a^*b)^*a^*, \quad (5)$$

$$(ab)^*a = a(ba)^*, \quad (6)$$
<figure><img src="image_15.png" /></figure>

$$a(a \circled (ax)) = ax, \qquad (7)$$

$$(ab) \circled x = b \circled (a \circ \setminus x). \qquad (8)$$

## 2.2. Max-plus linear systems

A max-plus linear system is defined by the following equations:

$$x(k) = Ax(k-1) \oplus Bu(k) \\ y(k) = Cx(k), \qquad (9)$$

where $x(k) \in \bar{Z}_{\max}^n$, $u(k) \in \bar{Z}_{\max}^p$, $y(k) \in \bar{Z}_{\max}^q$ and $k \in \mathbb{Z}$.

This kind of system makes it possible to describe the behaviors of TEGs, by associating to each transition a firing date sequence $x_i(k) \in \bar{Z}_{\max}$ and predict the system evolution. For a state equation in (9), each increasing sequence $\{x(k)\}$, it is possible to define the transformation $X(y) = \bigoplus_{k \in \mathbb{Z}} x(k)y^k$ where $y$ is a backward shift operator in event domain (i.e. $Y(y) = \gamma X(y) \Leftrightarrow \{y(k)\} = \{x(k-1)\}$), (see, Baccelli et al. (1992), p. 228). This transformation is analogous to the z-transform used in discrete-time classical control theory and the formal series $X(y)$ is a synthetic representation of the trajectory $x(k)$. The set of the formal power series in $y$ is denoted by $\bar{Z}_{\max}[[y]]$ and constitutes an idempotent semiring. Therefore, the state equation in (9) becomes a polynomial equation or an event-domain representation,

$$X(y) = \bar{A}X(y) \oplus BU(y), \quad \text{where} \quad \bar{A} = \gamma A, \\ Y(y) = CX(y), \qquad (10)$$

where the state $X(y) \in (\bar{Z}_{\max}[[y]])^n$, the output $Y(y) \in (\bar{Z}_{\max}[[y]])^q$, and the input $U(y) \in (\bar{Z}_{\max}[[y]])^p$, and matrices $\bar{A} \triangleq \gamma A \in (\bar{Z}_{\max}[[y]])^{n \times n}$, $B \in (\bar{Z}_{\max}[[y]])^{n \times p}$ and $C \in (\bar{Z}_{\max}[[y]])^{q \times n}$ represent the link between transitions. According to the state equation (10), the evolution of the system is

$$X(y) = \bar{A}^* BU(y), \\ Y(y) = C\bar{A}^* BU(y). \qquad (11)$$

The trajectories $U(y)$ and $Y(y)$ can be related (Baccelli et al. (1992), p. 243) by the equation $Y(y) = H(y)U(y)$, where $H(y) = C\bar{A}^*B \in (\bar{Z}_{\max}[[y]])^{q \times p}$ is called the transfer matrix of the TEG. Entries of matrix $H$ are periodic series (Baccelli et al. (1992), p. 260) in the idempotent semiring, usually represented by $p(y) \oplus q(y)(\tau y^\nu)^*$, where $p(y)$ is a polynomial representing the transient behavior, $q(y)$ is a polynomial corresponding to a pattern which is repeated periodically, the period being given by the monomial $(\tau y^\nu)$. For example, consider a transfer matrix $H(y) = 0y^0 \oplus 1y \oplus 3y^4 \oplus (5y^5 \oplus 6y^7)(3y^4)^*$ represented in Figure 1. The transient behavior is given by $p(y) = 0y^0 \oplus 1y \oplus 3y^4$ while the periodic behavior is given by $q(y)(\tau y^\nu)^* = (5y^5 \oplus 6y^7)(3y^4)^*$.
<figure><img src="image_16.png" /></figure>

<figure><img src="image_17.png" /><figcaption>Figure 1. Diagram showing a series with transient and periodic behavior (Hardouin et al., 2010).</figcaption></figure>

### 3. Modeling of HTS systems

We consider an HTS system from the High-Throughput Screening Center (HTSC) at Washington University in St. Louis whose layout is shown in Figure 2. The system consists of different resources as identified in Figure 2 on which different activities are

<figure><img src="image_18.png" /><figcaption>Figure 2. Layout of an HTS system. Source: HTSC at Washington University in St. Louis.</figcaption></figure>
<figure><img src="image_19.png" /><figcaption>Table 1. Resources and their associated tasks at the HTS plant at Washington University in St. Louis.</figcaption></figure>

<table><thead><tr><td>Resource</td><td>Device Type</td><td>Tasks</td></tr></thead><tbody><tr><td>Cytomat2C</td><td>Microplate hotel</td><td>Feeds microplates into the system.<br/>Stores microplates after they are processed.</td></tr><tr><td>Teleshake</td><td>Incubator (4-shelf microplate shaker)</td><td>Allows the biochemical substances to incubate in a controlled environment.</td></tr><tr><td>Envision</td><td>Reader</td><td>Analyzes chemical changes occurring in biochemical substances.</td></tr><tr><td>SCARA</td><td>Robot arm (transportation device)</td><td>Transfers microplates between resources.</td></tr><tr><td>Multidrop-384</td><td>Liquid dispenser</td><td>Dispenses liquid substance into the microplates' wells.</td></tr></tbody></table>

carried out. Table 1 lists the resources used in the assay and their associated tasks. Cytomat2C is the microplate hotel where all the microplates are kept and where they are returned after all the processes are completed. The microplates have either 384 wells or 96 wells each. Teleshake is the station where the incubation occurs while reading is done in Envision. SCARA is the robot doing all of the transfer activities between resources while necessary biochemical substances are added to the microplates in the Multidrop.

An enzymatic assay process that consists of 11 activities was carried out on the HTS system of Figure 2. To develop a schedule for the process, we are interested in the start and release times of each activity of the process.

Figure 3 shows the start and release times of each of the activities. $o_i$ is the start time of activity $i$ and $r_i$ its release time. It should be noted that for this particular process, the start time of activity $i$ corresponds to the release time of activity ($i - 1$). It can be seen that the 11 activities of the process are carried out on 5 resources. Activities 1 and 11 are carried out on the Cytomat2C; activities 2, 4, 6, 8, and 10 are transfer activities carried out by the SCARA robot; activities 5 and 9 are done on the Envision. The Multidrop and the Teleshake are single activity resources responsible for activity 3 and activity 7, respectively. All the resources are single capacity resources, except for the Teleshake, which has a capacity of 4. The assay process utilized full 384 -well microplates. If the microplates were not full or 96 -well microplates were used, the duration of activities 5 and 9, which are done on the Envision, might be reduced. For this process, the duration of other activities that include transfer, liquid dispensing and incubation would remain unchanged because they are independent of the number of wells. The Gantt chart for one batch of the assay process is shown in Figure 4. We see that the makespan of a single process is 491 s. It can also be seen from the Gantt chart that there is some overlap between some

<figure><img src="image_20.png" /><figcaption>Figure 3. Diagram showing the start and release times of each activity of the assay process.</figcaption></figure>
<figure><img src="image_21.png" /></figure>

<figure><img src="image_22.png" /><figcaption>Figure 4. Gantt chart of a single batch of the assay process.</figcaption></figure>

activities during a batch in which case a batch occupies more than one resource at the same time. Such is the case during transfer activities 6 and 8 and during the incubation, activity 7.

## 4. Scheduling of HTS systems by max-plus algebra

A timed event graph (TEG), shown in Figure 5 was generated to capture the information provided in Figure 3. The TEG has 12 transitions representing the start and release events of the different activities. The 4 tokens in the TEG indicate that the Teleshake has a capacity of 4. Each transition can be described by max-plus algebra. For instance,

$$x_2(k) = x_1(k) \otimes 32 = 32x_1(k),$$

$$x_{12}(k) = 5 \otimes x_{11}(k) \oplus 0 \otimes x_9(k) \oplus 0 \otimes x_{13}(k - 4)$$

$$= 5x_{11}(k) \oplus x_9(k) \oplus x_{13}(k - 4).$$

In $\overline{\mathbb{Z}}_{\max}[\{y\}]$ domain, these equations become

$$X_2(y) = 32X_1(y),$$

$$X_{12}(y) = 5X_{11}(y) \oplus X_9(y) \oplus \gamma^4 X_{13}(y).$$

<figure><img src="image_23.png" /><figcaption>Figure 5. TEG of the assay process.</figcaption></figure>
<figure><img src="image_24.png" /></figure>

The HTS system can be modeled as a max-plus linear system as in Eq. (10). The elements of matrix A are the states of the system and are given by the minimal time of each activity as shown in Figure 5. Matrices B and C are the input and output matrices respectively. The output signal is the release time of the last activity. The elements of matrix A that are not ε are

$$ 
\begin{array}{l@{\quad}l@{\quad}l@{\quad}l@{\quad}l}
A_{1,22} = \gamma & A_{2,1} = 32 & A_{3,2} = 19 & A_{4,3} = 9 & A_{5,4} = 23 & A_{6,5} = 54 \\
A_{6,8} = e & A_{7,10} = \gamma & A_{8,6} = 18 & A_{8,7} = 1 & A_{9,8} = 18 & A_{9,12} = e \\
A_{10,9} = 1 & A_{11,14} = \gamma & A_{12,11} = 5 & A_{12,9} = e & A_{12,13} = \gamma^4 & A_{13,12} = 210 \\
A_{13,16} = e & A_{14,13} = 5 & A_{15,18} = \gamma & A_{16,15} = 1 & A_{16,13} = e & A_{17,16} = 18 \\
A_{18,17} = 1 & A_{19,17} = e & A_{20,19} = 54 & A_{21,20} = 20 & A_{22,21} = 34, & \\
\end{array}
$$

where $A_{ij}$ gives the element at row i and column j in matrix A.

The elements of matrices B and C that are not ε are

$$ B_{1,1} = B_{7,2} = B_{11,3} = B_{15,4} = C_{1,22} = e, $$

where $B_{ij}$ and $C_{ij}$ give the element at row i and column j in matrices B and C.

It is trivial to obtain an off-line schedule for the system by computing its output trajectory as given in Eq. (11). To obtain an on-line schedule in the presence of disturbances, we need to find an optimal controller that solves the disturbance decoupling problem as explained in Shang et al. (2013, 2016a).

## 4.1. Online scheduling for HTS system

Disturbances are unforeseen disruptions to the normal evolution of an HTS system where a disturbance can be a breakdown of the robot arm or an obstruction in its way. Disturbances are usually modelled as uncontrollable inputs in a system. For our running example, we assume that there is a disruption in the movement of the robot arm. These disturbances are shown as inputs $q_1$ to $q_5$ in Figure 6. The system equations become

$$ X(y) = \bar{A}X(y) \oplus BU(y) \oplus SQ(y), $$

$$ Y(y) = CX(y), $$

<figure><img src="image_25.png" /><figcaption>Figure 6. TEG of the assay process with disturbances.</figcaption></figure>
<figure><img src="image_26.png" /></figure>

where matrices A, B and C are as earlier defined and the non-ε elements of the
disturbance matrix are:

$$
S_{2,1} = S_{4,2} = S_{8,3} = S_{16,4} = S_{20,5} = e.
$$

**Definition 7.** (Shang, Hardouin, Lhommeau, Maia et al. 2016a) *The system (12) is called **disturbance decoupled** by a state-feedback control $u(k) = Fx(k-1) \oplus v(k)$ (or an open-loop control $u(k) = Pv(k)$) if and only if any disturbance signal will not affect the system output $y(k)$ for all $k \in \mathbb{Z}$ and for any initial condition $x_0$.*

Solving the DDP means that the control $U(y)$ has to achieve

$$
\overline{CA}^*BU(y) \oplus \overline{CA}^*SQ(y) = \overline{CA}^*BU(y) \quad (12)
$$

$$
\Leftrightarrow \overline{CA}^*SQ(y) \preceq \overline{CA}^*BU(y). \qquad (13)
$$

The control input is computed from an external input, V(y), and the disturbance, Q(y),
which is measurable, as U(y) = PQ(y) ⊕ V(y). For any external input V(y) and dis-
turbance Q(y), the inequality above is equivalent to

$$
\overline{CA}^*SQ(y)\overline{CA}^*BPQ(y) \oplus \overline{CA}^*BV(y). \quad (14)
$$

**Definition 8.** (Shang, Hardouin, Lhommeau, Maia et al. 2016a) *The max-plus linear system described in (12) is called **modified disturbance decoupled** by a state feedback control $u(k) = Fx(k-1) \oplus v(k)$ (or an open-loop control $u(k) = Pv(k)$) if and only if the system output signals will not be disturbed more than the output signals influenced by the disturbances.*

The objective of the modified disturbance decoupling problem (MDDP) is to find
the greatest open-loop or output feedback control U(y) such that the output
trajectories will not be delayed more than the disturbance signals have acted on
the system. This means finding the greatest control, U(y), such that the following
equation holds,

$$
\overline{CA}^*BU(y) \oplus \overline{CA}^*SQ(y) = \overline{CA}^*SQ(y) \quad (15)
$$

$$
\Leftrightarrow \overline{CA}^*BU(y)\overline{CA}^*SQ(y). \qquad (16)
$$

The scheduling problem can be transformed into a DDP or MDDP in max-plus linear
systems. This is done by finding a prefilter or feedback controller that generates a control
input to solve the DDP or MDDP. Such a control input generates an optimal schedule in
spite of system malfunctions and delays, based on the just-in-time criterion according to
the industry standard.

If the control $U(y) = PQ(y)$, where $P$ is a prefilter which generates the control by
taking the measured disturbances into account, then, solving the MDDP is equivalent to
finding a prefilter $P$ satisfying
<figure><img src="image_27.png" /></figure>

$$ \overline{CA}^* BPQ(y) \oplus \overline{CA}^* SQ(y) = \overline{CA}^* SQ(y), \forall Q(y). \quad (17) $$

The optimum prefilter matrix $P_{opt}$ that solves both the DDP and MDDP is given by

$$ P_{opt} = (\overline{CA}^* B) \circledast (\overline{CA}^* S). \quad (18) $$

According to Shang et al. (2016a), the optimal feedback controller $F_{opt}$ that generates a schedule according to the just-in-time criterion is given by

$$ F_{opt} = (\overline{CA}^* B) \circledast (\overline{CA}^* \tilde{B}) \circledast (\overline{CA}^* \tilde{B}) \text{ where } \tilde{B} = [B|S]. $$

When simulated in MinMaxGD, a C++ toolbox for handling periodic series introduced in Cottenceau et al. (2000), we obtain

$$ P_{opt} = \begin{bmatrix} 459\gamma(491\gamma)^* & 431\gamma(491\gamma)^* & 354\gamma(491\gamma)^* & 126\gamma(491\gamma)^* & 54\gamma(491\gamma)^* \\ 104(491\gamma)^* & 76(491\gamma)^* & 490\gamma(491\gamma)^* & 262\gamma(491\gamma)^* & 190\gamma(491\gamma)^* \\ 118(491\gamma)^* & 90(491\gamma)^* & 13(491\gamma)^* & 276\gamma(491\gamma)^* & 204\gamma(491\gamma)^* \\ 332(491\gamma)^* & 304(491\gamma)^* & 227(491\gamma)^* & 490\gamma(491\gamma)^* & 418\gamma(491\gamma)^* \end{bmatrix}, $$

and

$$ F_{opt} = \begin{bmatrix} \gamma(491\gamma)^* \\ 136\gamma(491\gamma)^* \\ 150\gamma(491\gamma)^* \\ 364\gamma(491\gamma)^* \end{bmatrix}. $$

Each element in $P_{opt}$ and $F_{opt}$ is periodic with a period of 491 time units. For instance,

$$ \begin{aligned} P_{opt}(2, 1) &= 104(491\gamma)^* \\ &= 104(e \oplus 491\gamma \oplus 982\gamma^2 \oplus 1473\gamma^3 \oplus 1964\gamma^4 \oplus \dots) \\ &= 104 \oplus 595\gamma \oplus 1086\gamma^2 \oplus 1577\gamma^3 \oplus 2068\gamma^4 \oplus \dots. \end{aligned} $$

A realization of the TEG together with the controllers is shown in Figure 7. It can be seen that the disturbances pass through the prefilter before going into the system. We see that

$$ V = P_{opt}Q, \text{ and} $$

$$ U = \hat{F}CX \oplus V = \hat{F}Y \oplus V. $$

For instance, the first external input, $V_1$, is gotten from the first row of $P_{opt}$ and the disturbance vector as follows:

$$ V_1 = \begin{bmatrix} 459\gamma(491\gamma)^* & 431\gamma(491\gamma)^* & 354\gamma(491\gamma)^* & 126\gamma(491\gamma)^* & 54\gamma(491\gamma)^* \end{bmatrix} \begin{bmatrix} q_1 \\ q_2 \\ q_3 \\ q_4 \\ q_5 \end{bmatrix}; $$

and the first control input $U_1$ is
<figure><img src="image_28.png" /><figcaption>Figure 7. TEG model of the HTS system with integrated prefilter and the output controllers.</figcaption></figure>

$$U_1 = \gamma(491\gamma)^* y \oplus V_1.$$

In $\gamma(491\gamma)^*$, the cycle time is 491, and $\gamma$ represents a token in the cyclic mode.

## 5. Online scheduling of uncertain HTS systems

For the HTS process we have been considering, the time delays indicated in Figure 5 and Figure 6 give the minimum amount of time to be spent on each activity. In some cases, certain activities have a lower limit and an upper limit of the time they take rather than a minimal time. In such a situation where the times for some activities are not fixed, we have what we call an uncertain HTS system, which can be described by uncertain max-plus linear system.

An uncertain max-plus linear system is defined as

$$x(k) = Ax(k-1) \oplus Bu(k) \oplus Sq(k),$$

$$y(k) = Cx(k), \qquad (19)$$

where system matrices lie in corresponding matrix intervals with known lower and upper bounds, specifically, $A \in A = [A_l, A_u] \in I(\bar{Z}_{\max})^{n \times n}$, $B \in B = [B_l, B_u] \in I(\bar{Z}_{\max})^{n \times p}$, $S \in S = [S_l, S_u] \in I(\bar{Z}_{\max})^{n \times r}$, and $C \in C = [C_l, C_u] \in I(\bar{Z}_{\max})^{q \times n}$. And the states are $x(k) \in X \cong \bar{Z}_{\max}^n$, the inputs are $u(k) \in U \cong \bar{Z}_{\max}^p$, the disturbances are $q(k) \in Q \cong \bar{Z}_{\max}^r$, and the outputs are $y(k) \in Y \cong \bar{Z}_{\max}^q$ and $k \in \mathbb{Z}$.

An uncertain HTS system with disturbances can be disturbance decoupled by different control inputs in much the same way as for systems without uncertainties. However, the controls have to be selected very restrictively in order to solve the DDP and the MDDP for all the variations of the system matrices, as well as for arbitrary disturbances.
<figure><img src="image_29.png" /></figure>

Sometimes, such non-trivial controls do not even exist for all variations of system matrices. Therefore, we relax the definitions of the DDP and MDDP in order to achieve more general results.

## 5.1. Weakly DDP and weakly MDDP

**Definition 9. (WDDP)** The uncertain max-plus linear system in (19) is called weakly disturbance decoupled if and only if there exists a control U belonging to an interval U such that the output trajectory interval, generated by the disturbances and the uncertain matrices, belongs to the output trajectory interval, generated by the control U. In other words, solving WDDP means that the following equality

$$ \overline{\mathbf{C}\mathbf{A}^*\mathbf{B}\mathbf{U}} \oplus \overline{\mathbf{C}\mathbf{A}^*\mathbf{S}\mathbf{Q}} = \overline{\mathbf{C}\mathbf{A}^*\mathbf{B}\mathbf{U}}, \quad (20) $$

holds for intervals of matrices A, B, C, S, and the disturbance interval Q.

**Definition 10. (WMDDP)** The max-plus linear system described in (19) is called weakly modified disturbance decoupled if and only if it exists a control U belonging to an interval U such that the output trajectories belong to the interval generated by any disturbances and uncertain matrices. In other words, solving WMDDP means that the following equality

$$ \overline{\mathbf{C}\mathbf{A}^*\mathbf{B}\mathbf{U}} \oplus \overline{\mathbf{C}\mathbf{A}^*\mathbf{S}\mathbf{Q}} = \overline{\mathbf{C}\mathbf{A}^*\mathbf{S}\mathbf{Q}}, \quad (21) $$

is satisfied for intervals of matrices A, B, C, S, and the disturbance interval Q.

The controllers solving the DDP and MDDP have to solve them for any arbitrary variations of system matrices, as well as the disturbances. On the other hand, the WDDP and WMDDP only requires that there exists an interval of control U which induces an output interval unchanged by the disturbance Q(y) in an interval Q.

## 5.2. Open-loop controller interval solving WMDDP

Solving the WMDDP by an open-loop controller requires finding an optimal prefilter interval $\mathbf{P}_{\text{opt}}^{\text{WMDDP}} = [P_{lopt}, P_{uopt}]$ such that the control interval $\mathbf{U} = \mathbf{PQ}$ solves the WMDDP. $\mathbf{U} = [U_l, U_u]$, $\mathbf{P} = [P_l, P_u]$ and $\mathbf{Q}$ is interval of the disturbances. $\mathbf{P}_{\text{opt}}^{\text{WMDDP}}$ ensures the following equality holds

$$ \overline{\mathbf{C}\mathbf{A}^*\mathbf{B}} \mathbf{P}_{\text{opt}}^{\text{WMDDP}} \oplus \overline{\mathbf{C}\mathbf{A}^*\mathbf{S}} = \overline{\mathbf{C}\mathbf{A}^*\mathbf{S}}, \quad (22) $$

$$ \Leftrightarrow \overline{\mathbf{C}\mathbf{A}^*\mathbf{B}} \mathbf{P}_{\text{opt}}^{\text{WMDDP}} \preceq \overline{\mathbf{C}\mathbf{A}^*\mathbf{S}}. $$

The optimal prefilter interval was found in Shang, Hardouin, Lhommeau, Maia et al. (2016a) as

$$ \mathbf{P}_{\text{opt}}^{\text{WMDDP}} \triangleq \overline{\mathbf{C}\mathbf{A}^*\mathbf{B}} \overline{\mathbf{B}^T} \overline{\mathbf{C}\mathbf{A}^*\mathbf{S}} = [P_{lopt}^{\text{WMDDP}}, P_{uopt}^{\text{WMDDP}}], \quad (23) $$

where
<figure><img src="image_30.png" /></figure>

<figure><img src="image_31.png" /><figcaption>Figure 8. The output feedback structure.</figcaption></figure>

$$P_{lopt}^{WMDDP} = C_l \overline{A}_l^* B_l \circled{C} C_l \overline{A}_l^* S_l \wedge C_u \overline{A}_u^* B_u \circled{C} C_u \overline{A}_u^* S_u,$$

$$P_{uopt}^{WMDDP} = C_u \overline{A}_u^* B_u \circled{C} C_u \overline{A}_u^* S_u. \qquad (24)$$

### 5.3. Feedback controller interval solving WMDDP

As introduced in Shang, Hardouin, Lhommeau, Maia et al. (2016a), the output feedback controller structure solving the WMDDP is shown in Figure 8. The output feedback control interval that solves the WMDDP is $\mathbf{U} = \mathbf{FY} \oplus \mathbf{V} = \mathbf{FCX} \oplus \mathbf{V}$. If we define $\tilde{\mathbf{H}} = C\overline{\mathbf{A}}^*\tilde{\mathbf{B}}$, the controller $\mathbf{F}$ has to satisfy

$$\left(\tilde{\mathbf{H}}\left(\begin{matrix} e \\ \varepsilon \end{matrix}\right)\right)^* \tilde{\mathbf{H}} \preceq \tilde{\mathbf{H}}, \qquad (25)$$

for all $\tilde{\mathbf{H}} \in \tilde{\mathbf{H}} = [C_l \overline{A}_l^* \tilde{B}_l, C_u \overline{A}_u^* \tilde{B}_u]$, $\tilde{B}_l = [B_l | S_l]$ and $\tilde{B}_l = [B_u | S_u]$.

The interval of output feedback controls $\mathbf{F} = [F_l, F_u]$ is

$$\begin{align}
\widehat{\mathbf{F}} &= \mathbf{H} \widehat{\circled{\overline{\mathbf{H}}}} \widehat{\circled{\overline{\mathbf{H}}}} \widehat{\phi} \widehat{\mathbf{H}} = [\hat{F}_{lopt}, \hat{F}_{uopt}] \nonumber \\
&= [H_l \circled{\overline{\mathbf{H}}}_l \circled{\overline{\mathbf{H}}}_l \wedge H_u \circled{\overline{\mathbf{H}}}_u \circled{\overline{\mathbf{H}}}_u, H_u \circled{\overline{\mathbf{H}}}_u \circled{\overline{\mathbf{H}}}_u \wedge \widehat{\mathbf{H}}_l] \tag{26}
\end{align}$$

## 6. Simulation results

The uncertain case of the HTS assay process is considered to simulate controllers that preserve open-loop behavior in the presence of disturbance. The incubation time in the uncertain HTS process has an upper and a lower bound; it takes up to 210 s and no more than 235 s as shown in Figure 9.

The controllers given in Eqs (24) and (26) are computed for the HTS system. The optimal lower and upper bounds of the prefilter are
<figure><img src="image_32.png" /></figure>

<figure><img src="image_33.png" /><figcaption>Figure 9. TEG of the assay process of the uncertain case with disturbances.</figcaption></figure>

$$P_{lopt}^{WMDDP} = \begin{bmatrix} 459\gamma(491\gamma)^* & 431\gamma(491\gamma)^* & 354\gamma(491\gamma)^* & 126\gamma(491\gamma)^* & 54\gamma(491\gamma)^* \\ 104(491\gamma)^* & 76(491\gamma)^* & 490\gamma(491\gamma)^* & 262\gamma(491\gamma)^* & 190\gamma(491\gamma)^* \\ 118(491\gamma)^* & 90(491\gamma)^* & 13(491\gamma)^* & 276\gamma(491\gamma)^* & 204\gamma(491\gamma)^* \\ 332(491\gamma)^* & 304(491\gamma)^* & 227(491\gamma)^* & 490\gamma(491\gamma)^* & 418\gamma(491\gamma)^* \end{bmatrix},$$

and

$$P_{uopt}^{WMDDP} = \begin{bmatrix} 484\gamma(516\gamma)^* & 456\gamma(516\gamma)^* & 379\gamma(516\gamma)^* & 126\gamma(516\gamma)^* & 54\gamma(516\gamma)^* \\ 104(516\gamma)^* & 76(516\gamma)^* & 515\gamma(516\gamma)^* & 262\gamma(516\gamma)^* & 190\gamma(516\gamma)^* \\ 118(516\gamma)^* & 90(516\gamma)^* & 13(516\gamma)^* & 276\gamma(516\gamma)^* & 204\gamma(516\gamma)^* \\ 357(516\gamma)^* & 329(516\gamma)^* & 252(516\gamma)^* & 515\gamma(516\gamma)^* & 443\gamma(516\gamma)^* \end{bmatrix}.$$

The optimal lower and upper bounds of the output feedback controller are

$$\hat{F}_{lopt} = \begin{bmatrix} \gamma(491\gamma)^* \\ 136\gamma(491\gamma)^* \\ 150\gamma(491\gamma)^* \\ 364\gamma(491\gamma)^* \end{bmatrix},$$

and

$$\hat{F}_{uopt} = \begin{bmatrix} \gamma(516\gamma)^* \\ 136\gamma(516\gamma)^* \\ 150\gamma(516\gamma)^* \\ 389\gamma(516\gamma)^* \end{bmatrix}.$$

In Figure 10, the upper-bound prefilter, $P_{uopt}^{WMDDP}$ and the upper-bound output feedback controller, $\hat{F}_{uopt}$ are integrated into the TEG of the HTS process. The Gantt chart of the process is shown in Figure 11. The states corresponding to the activities marked in the Gantt chart are given in Table 2. It is seen from the Gantt chart that the prefilter eliminates unnecessary wait time caused by the disturbance. With the prefilter, the completion time of each activity satisfies the just-in-time criterion.
<figure><img src="image_34.png" /><figcaption>Figure 10. TEG model of the HTS system with integrated upper bound prefilter and the output controllers.</figcaption></figure>

This is better seen in Tables 3 and 4. Table 3 shows the start times, release times and durations of each activity with and without disturbance and when the prefilter is included. We see, for example, that activity 2 has a duration of 516 seconds without prefilter. This is reduced to just 20 s, which is the minimum required time, by the prefilter. Table 4 shows that this is a 96% improvement. Therefore, the prefilter ensures a just-in-time completion time.

## 7. Discussion

This study has successfully obtained a max-plus algebraic model for an HTS system that is operated cyclically. We take further the method introduced in Brunsch and Raisch (2009) for obtaining an online schedule for cyclically operated HTS systems with disturbances by obtaining a schedule for an HTS system with disturbances and uncertainties. This study also presented optimal controllers that solve the DDP, and thereby solving the scheduling problem for an HTS system with disturbances and uncertainties. Whereas in Shang et al. (2016b), controllers that solve DDP for uncertain max-plus linear systems whose system matrices vary between intervals was introduced, we show that the controllers solving DDP also solve the scheduling problem in cyclically operated discrete event systems. The supervisory control proposed in this study will be useful for generating static schedules for HTS systems and other discrete event systems. More importantly, it can be used to generate online schedule, which is important for when there are disturbances or deviations from a predetermined schedule during operation.

In summation, an HTS process is first modeled in max-plus algebra and a static schedule is developed for the process based on the model. Online scheduling was then considered by computing controllers that ensure an optimal schedule when there are disturbances. The result is extended to the uncertain case where system matrices lie in intervals. This paper
<figure><img src="image_35.png" /></figure>

<figure><img src="image_36.png" /><figcaption>Figure 11. Gantt charts of the process with disturbance and with disturbance and prefilter.</figcaption></figure>
<figure><img src="image_37.png" /><figcaption>**Table 2.** Activities and corresponding states.</figcaption></figure>

<table><caption>**Table 3.** Comparison of activity duration without prefilter and with prefilter for $k = 1$.</caption><thead><tr><td>Activity</td><td>States</td></tr></thead><tbody><tr><td>1</td><td>x<sub>1</sub> - x<sub>6</sub></td></tr><tr><td>2</td><td>x<sub>7</sub> - x<sub>10</sub></td></tr><tr><td>3</td><td>x<sub>11</sub> - x<sub>14</sub></td></tr><tr><td>4</td><td>x<sub>15</sub> - x<sub>18</sub></td></tr><tr><td>5</td><td>x<sub>19</sub> - x<sub>22</sub></td></tr></tbody></table>

<table><caption>**Table 4.** Waiting time eliminated with prefilter.</caption><thead><tr><th rowspan="2">Activity</th><th colspan="3">Without prefilter</th><th colspan="3">With prefilter</th></tr><tr><th>Start time (s)</th><th>Release time (s)</th><th>Duration (s)</th><th>Start time (s)</th><th>Release time (s)</th><th>Duration (s)</th></tr></thead><tbody><tr><td>1</td><td>484</td><td>621</td><td>137</td><td>484</td><td>621</td><td>137</td></tr><tr><td>2</td><td>124</td><td>640</td><td>516</td><td>620</td><td>640</td><td>20</td></tr><tr><td>3</td><td>363</td><td>879</td><td>516</td><td>634</td><td>879</td><td>245</td></tr><tr><td>4</td><td>377</td><td>893</td><td>516</td><td>873</td><td>893</td><td>20</td></tr><tr><td>5</td><td>892</td><td>1000</td><td>108</td><td>892</td><td>1000</td><td>108</td></tr></tbody></table>

<table><thead><tr><th>Activity</th><th>Waiting time eliminated</th><th>% Time saved</th></tr></thead><tbody><tr><td>1</td><td>0</td><td>0</td></tr><tr><td>2</td><td>496</td><td>96</td></tr><tr><td>3</td><td>275</td><td>52</td></tr><tr><td>4</td><td>496</td><td>96</td></tr><tr><td>5</td><td>0</td><td>0</td></tr></tbody></table>

presents the solution of the scheduling problem for HTS systems, with and without uncertainties, by solving different types of disturbance decoupling problems of max-plus linear systems. Open-loop controls and the output feedback controls solving the DDP, which ultimately solves an online scheduling problem, were constructed using the residualtion the ory and proved to be optimal schedules for the just-in-time control criterion.

## 8. Future research

In future research, the method employed here will be extended to other high-throughput processes such as High-Throughput *in vitro* Absorption, Distribution, Metabolism and Excretion (HT-ADME) and High-Throughput Mass Spectrometry (HT-MS) by modeling these processes as a max-plus linear systems. Of particular interest in future work will be HTS systems with nested cycles. The authors are working on other optimization models such as mixed integer programming and queueing theory to obtain additional results and to compare with the results obtained here.
<figure><img src="image_38.png" /></figure>

## Acknowledgments

The authors would like to thank Dr Ma Xenia Ilagan, Associate Professor in the Department of Biochemistry and Molecular Biophysics, for providing access to the data from the High Throughput Screening Center at Washington University in St. Louis.

## Disclosure statement

No potential conflict of interest was reported by the author(s).

## References

* Baccelli, F., Cohen, G., Olsder, G. J., & Quadrat, J. P. (1992). *Synchronization and linearity: An algebra for discrete event systems*. John Wiley and Sons. http://www-rocq.inria.fr/metalau/cohen/SED/SED1-book.html

* Brunsch, T., & Raisch, J. (2009). Max-plus algebraic modeling and control of high-throughput screening systems. *The 2nd IFAC workshop on dependable control of discrete systems*.

* Brunsch, T., Raisch, J., & Hardouin, L. (2012). Modeling and control of high-throughput screening systems. *Control Engineering Practice, 20*(1), 14–23. https://doi.org/10.1016/j.conengprac.2010.12.006

* Cohen, G., Gaubert, S., & Quadrat, J. P. (1996). Kernels, images, and projections in dioids. *The IEEE International Workshops on Discrete Event Systems*, Edinburgh, Scotland, UK. http://www-rocq.inria.fr/metalau/quadrat/kernel.pdf

* Cohen, G., Gaubert, S., & Quadrat, J. P. (1997). Linear projectors in the max-plus algebra. *The 5th IEEE-Mediterranean Conference*, Paphos, Cyprus. http://www-rocq.inria.fr/metalau/quadrat/projector.pdf

* Cohen, G., Gaubert, S., & Quadrat, J. P. (2006). Projection and aggregation in maxplus algebra. In *Current Trends in Nonlinear Systems and Control, in Honor of Peter Kokotovic and Turi Nicosia*. Birkhäuser. 422–442. http://www-rocq.inria.fr/metalau/quadrat/Rome.pdf

* Cottenceau, B., Hardouin, L., Lhommeau, M., & Boimond, J.-L. (2000). Data processing tool for calculation in dioid. *The IEEE International Workshop on Discrete Event Systems*, Ghent, Belgium. http://www.istia.univ-angers.fr/hardouin/outils.html

* Cottenceau, B., Lhommeau, M., Hardouin, L., & Boimond, J.-L. (2003). On timed event graph stabilization by output feedback in dioid. *Kybernetika*, 39(2), 165–176. https://www.kybernetika.cz/content/2003/2/165/paper.pdf

* Golan, J. S. (1999). *Semirings and Their Applications*. Kluwer Academic Publishers.

* Hardouin, L., Maia, C. A., Cottenceau, B., & Lhommeau, M. (2010). Observer design for (max, plus)-linear systems. *The IEEE Transactions on Automatic Control*, 5(2), 538–543. https://doi.org/10.1109/TAC.2009.2037477

* Hardouin, L., Lhommeau, M., & Shang, Y. (2011). Towards geometric control of max-plus linear systems with applications to manufacturing systems. *The 50th IEEE Conference on Decision and Control and European Control Conference*, Orlando, FL, USA, 1149–1154.

* Major, J. (1998). Challenges and opportunities in high throughput screening: Implications for new technologies. *Journal of Biomolecular Screening*, 3(1), 13–17. https://doi.org/10.1177/108705719800300102

* Mayer, E., & Raisch, J. (2003). Throughput-optimal scheduling for cyclically repeated processes. *Proceedings of the 9th IEEE International Conference on Methods and Models in Automation and Robotics*. Miedzyzdroje, Poland, 871–876.

* Mayer, E., & Raisch, J. (2004). Time-optimal scheduling for high throughput screening processes using cyclic discrete event models. *Mathematics and Computers in Simulation*, 66(2–3), 181–191. https://doi.org/10.1016/j.matcom.2003.11.004
Mayer, E., Haus, -U.-U., Raisch, J., & Weismantel, R. (2008). Throughput-optimal sequences for cyclically operated plants. *Discrete Event Dynamic Systems, 18*(3), 355-383. https://doi.org/10.1007/s10626-008-0038-3

Mayr, L. M., & Fuerst, P. (2008). The future of high-throughput screening. *Journal of Biomolecular Screening, 13*(6), 443-448. https://doi.org/10.1177/1087057108319644

Noah, J. W. (2010). New developments and emerging trends in high-throughput screening methods for lead compound identification. *International Journal of High Throughput Screening, 1*(1), 141-149. https://doi.org/10.2147/IJHTS.S8683

Pereira, D. A., & Williams, J. A. (2007). Origin and evolution of high throughput screening. *British Journal of Pharmacology, 153*(1), 53-61. https://doi.org/10.1038/sj.bjp.0707373

Shang, Y., Hardouin, L., Lhommeau, M., & Maia, C. A. (2013). Open-loop controllers to solve the disturbance decoupling problem for max-plus linear systems. *European Control Conference, Zurich, Switzerland*, 125-130.

Shang, Y., Hardouin, L., Lhommeau, M., & Maia, C. A. (2016a). An integrated control strategy to solve the disturbance decoupling problem for max-plus linear systems with applications to a high throughput screening system. *Automatica, 63*, 338-348. https://doi.org/10.1016/j.auto-matica.2015.10.030

Shang, Y., Hardouin, L., Lhommeau, M., & Maia, C. A. (2016b). Robust controllers in disturbance decoupling of uncertain max-plus linear systems: An application to a high throughput screening system for drug discovery. *The 13th International Workshop on Discrete Event Systems (WODES)*, 404-409. https://doi.org/10.1109/WODES.2016.7497880