# Toward Robust Manufacturing Scheduling: Stochastic Job-Shop Scheduling

Mikhail A. Bragin, , Matthew E. Wilhelm, Nanpeng Yu, and Matthew D. Stuber
Mikhail A. Bragin is with the Institute for Advanced Systems Engineering,
University of Connecticut, Storrs, CT, 06269, USA; Matthew E. Wilhelm and
Matthew D. Stuber are with the Process Systems and Operations Research
Laboratory, Department of Chemical and Biomolecular Engineering, University of
Connecticut, Storrs, CT, 06269, USA; Nanpeng Yu is with the Department of
Electrical and Computer Engineering, University of California Riverside,
Riverside, CA, 92521, USA;

###### Abstract

Manufacturing plays a significant role in economic development, production,
exports, and job creation, which ultimately contribute to improving the
quality of life. The presence of manufacturing defects is, however, inevitable
leading to products being discarded, i.e. scrapped. In some cases, defective
products can be repaired through rework. Scrap and rework cause a longer
completion time, which can contribute to orders being shipped late. Moreover,
the presence of uncertainties and combinatorial complexity significantly
increases the difficulty of complex manufacturing scheduling. This paper
tackles this challenge, exemplified by a case study on stochastic job-shop
scheduling in low-volume, high-variety manufacturing contexts. To ensure on-
time delivery, high-quality solutions are required, and near-optimal solutions
must be obtained within strict time constraints to ensure smooth operations on
the job-shop floor. To efficiently solve the stochastic job-shop scheduling
(JSS) problem, a recently-developed Surrogate “Level-Based” Lagrangian
Relaxation is used to reduce computational effort while efficiently exploiting
the geometric convergence potential inherent to Polyak’s step-sizing formula
thereby leading to fast convergence. Numerical testing demonstrates that the
new method is two orders of magnitude faster as compared to commercial
solvers.

Note to Practitioners—Manufacturing defects leading to scrap or rework create
significant challenges, from both computational as well as on-time delivery
standpoints. To assist practitioners in overcoming these challenges, this
paper presents a case study of stochastic job-shop scheduling as well as
Surrogate “Level-Based” Lagrangian Relaxation as an efficient solution
methodology, which not only reduces the time required to obtain high-quality,
near-optimal solutions for complex scheduling problems but also minimizes the
dependence on hyperparameters. As a result, the method becomes more user-
friendly, reducing the need for domain knowledge or instance-specific
heuristical adjustments. By implementing the proposed method, practitioners in
the manufacturing industry can expect a substantial increase in the speed of
obtaining solutions to their scheduling problems, expecting to significantly
outperform commercial solvers. This improvement will help organizations meet
delivery deadlines, maintain high-quality production standards, and ultimately
enhance their competitiveness in the market.

###### Index Terms:

Advanced Manufacturing, Stochastic Job-Shop Scheduling, Lagrangian Relaxation,
Markov Processes

##  I Introduction

Mmanufacturers’ ultimate goal is customer satisfaction, achieved through
delivering high-quality products on time to remain competitive in the
marketplace [1, 2, 3]. Late deliveries, however, can diminish the value of
products in the eyes of customers. Agile and flexible tailored production
greatly influence the value created. The trend towards improving flexibility
is through make-to-order production, allowing the ordering of bespoke products
[4, 5, 6]. Manufacturing operations must account for tardiness to avoid late
shipments and subsequent tarnishing reputations in the marketplace. Correcting
quality problems at the factory level is more cost-effective than resolving
quality problems after the product has reached the end customer [7]. To this
end, machining-level actions (e.g., feed rate and depth of cut) [8], as well
as tardiness and reactive re-scheduling, have been optimized [9].

Scope: This paper addresses tardiness-based job-shop scheduling (JSS), which
is crucial to avoid late shipments and maintain customer satisfaction. Job
shops are a specific type of manufacturing environment designed for low-
volume/high-variety production, making them ideal for producing bespoke
products. However, the JSS problem is difficult due to the combinatorial
complexity caused by discrete scheduling decision variables. As a result,
finding even near-optimal solutions can be challenging, leading to inefficient
operations and excessive costs, or long solving times that contribute to late
deliveries and associated consequences.

Additionally, the paper considers the possibility of damaged parts or orders
resulting in scrapped or reworked parts. The problem is modeled
stochastically, leading to a probabilistic routing of jobs that takes
uncertainties into account. Inaccurate decisions in the presence of
uncertainty can be particularly costly, as scrapping one operation of one job
may require remanufacturing the entire part, causing unexpected delays in
production that can lead to customer dissatisfaction. This paper aims to
obtain near-optimal solutions quickly to mitigate the impact of these
uncertainties and delays.

The rest of this section is organized as follows: In Subsection I-A, we review
JSS problems with the consideration of tardiness. In Subsection I-B, we review
stochastic JSS. In Subsection I-C, we review recent formulations of JSS
problems. Finally, in Subsection I-D, we provide a review of solution
methodologies for solving JSS problems and general discrete programming
problems.

###  I-A Reducing Tardiness: On-Time Delivery in Job Shops

The goal of a JSS problem is to assign a set of jobs, each consisting of a
sequence of consecutive non-preemptive (uninterruptible) operations, to be
processed on a set of machines. Generally, JSS problems have been recognized
as difficult combinatorial optimization problems [10, 11].

As previously reviewed, on-time delivery has long been recognized as an
important aspect of production; the minimization of tardiness has been used as
an objective within JSS problems as well. The seminal works on “tardy
problems” begin with the research of [12] and [13] with subsequent research
of, for example, [14], [15], [16], [17], [18], [19], and [20]. The problem is
known to be NP-hard (e.g., [17]), even single-machine scheduling [21, 22, 23,
17] belongs to a class of NP-hard problems. Yet, the job-shop scheduling
problems are also subject to strict computational requirements; long-solving
times are not allowed since they may contribute to orders shipped late.

For this reason, several rules have been developed in early implementations
such as “first-come, first-serve”, “shortest processing time”, or “earliest
operation due date” rules; these and several other rules are described and
compared by [19]. The minimization of tardiness, however, requires the
coordination of multiple jobs and multiple machines, while dispatching rules
typically ignore other resources. In fact, as pointed out by [24], “no single
priority rule performs well.” Therefore, most of the research has been heavily
invested in the development of solution methodologies beyond rule-based ones.

###  I-B Robustness to Uncertainty: Stochastic Job-Shop Scheduling

Uncertainty is a concomitant process within any production system and must be
accounted for. Therefore, it is important to account for potential causes of
disruptions to mitigate the consequences should such unforeseen events occur
in actual operations. Without the appropriate incorporation of uncertainty
into problem formulations, deterministic scheduling decisions may not be
optimal during operations, if feasible at all. In terms of JSS, neglecting
uncertainties within the optimization may result in unexpectedly high
tardiness in actual operations, thereby compromising the on-time delivery of
products.

In the extant research, uncertain processing times have been considered [25,
26, 27, 28, 29, 30, 31, 32, 33, 34]111Exceptions are, perhaps, [35, 36] where
machine unavailability is considered, though, both papers consider
minimization of makespan, not tardiness.. Within these papers, while the
processing time is uncertain, a job is assumed to be processed with certainty.
However, in actual operations, a product may be damaged (scrapped) after being
processed and must be remade anew thereby potentially leading to higher
tardiness. This addition of tardiness cannot be captured through stochastic
modeling of processing times, since a job may require several operations, all
of which may need to be repeated. Another salient feature of the above
existing methods is their meta-heuristic nature (with the exception of the
work by [25]). Heuristic-based algorithms generally do not provide a lower
bound to quantify the solution quality. Some notable exceptions are algorithms
such as those proposed by [37].

A fast resolution of JSS problems rests upon two pillars: formulation and
solution methodology, both are reviewed next.

###  I-C “Tight” Job-Shop Scheduling ILP/MILP Formulations

An aspect of problem formulation that has a direct impact on the solution
process is the formulation “tightening,” the gist of which is to delineate
facets of the convex hull. If successful, then the combinatorial problem
reduces to a much easier-to-solve LP problem; even if a problem is “partially”
tightened, the CPU time is greatly reduced [38, 39]. A recent advancement
reduces decision variables and constraints in ILP formulations, proving to be
tighter than previous ones [40]: only the beginning times of operations are
decision variables thereby reducing the number of decision variables and
constraints compared to those within the previous formulations. Since the
formulation has few decision variables and constraints, it leads to several
orders of magnitude improvement over previously-used approaches in terms of
CPU time and to optimality of solution obtained (from 3,600 s and 3.72% gap
down to 3.31 s and 0% gap) for a problem instance from [41], when solved by
using standard B&C [40].

###  I-D Previous Solution Methodologies

In this subsection, methods used to solve general scheduling problems are
briefly reviewed. Then, methods specifically developed for the JSS problem are
reviewed. Since the methodology of this paper is based on Lagrangian
Relaxation, the LR-based methods are reviewed as related to JSS problems.
Finally, the most recent advancements of LR are presented.

Methodologies for General Scheduling Problems. To solve scheduling problems
(including JSS problems minimizing makespan, and flow shop scheduling), the
following methods have been used: ant colony optimization [42], conflict-
directed search [43], constraint programming [44, 45], decomposition methods
[46], various version of evolutionary algorithms [47, 48], greedy algorithms
[49, 50], various versions of genetic algorithms [51, 52, 11, 53, 54],
heuristics/metaheuristics [37, 50, 55], Lagrangian Relaxation [56, 57], local
search [58, 59, 44], particle swarm optimization [60], Petri nets [61, 50],
reinforcement learning [62, 63, 64, 65] and taboo search [66].

Methodologies for Job-Shop Scheduling Problems with the Consideration of
Tardiness. Because of the complexity, JSS problems require a more
sophisticated methodology beyond the simple rules reviewed in subsection I-A.
Toward this goal, several methods have been developed such as Lagrangian
Relaxation [41, 25], genetic algorithms [67, 28], heuristics [26], neural
networks with simulated annealing [68], evolutionary algorithms [27, 29, 30],
ordinal optimization [31], and evolutionary learning-based simulation
optimization [69].

Lagrangian Relaxation for Job-Shop Scheduling Problems. To address the
combinatorial complexity inherent in Job Shop Scheduling (JSS) problems,
Lagrangian Relaxation (LR) has been employed to “reverse” the complexity,
leading to an exponential reduction in the effort required to solve
subproblems. For JSS problems, relaxing machine capacity constraints allows
the subproblems associated with individual jobs to be solved with
significantly reduced complexity. The application of standard Lagrangian
Relaxation in JSS problems can be traced back to the work of Hoitomt et al.
[41], where Lagrange multipliers were updated using the Polyak stepsize [70]
and sub-gradient directions. While there were other advancements in improving
the coordination of LR, in particular, to solve JSS problems [71, 38, 72, 40],
Polyak’s results have the potential to achieve geometric convergence (fastest
possible) and deserve a separate discussion.

Surrogate Sub-Gradient Method. Within the Surrogate Sub-Gradient Method [73],
multipliers are updated after solving one subproblem at a time rather than
solving all the subproblems as within the standard sub-gradient method. This
significantly reduces the computational effort. The step-sizes are updated by
using the following variation of the Polyak’s formula (originally developed in
[70]):

$\displaystyle
0<s^{k}<\gamma\cdot\frac{q(\lambda^{*})-L^{k}}{\left\|\tilde{g}^{k}\right\|^{2}_{2}},\gamma<1,$
(1)

where $L^{k}$ is the Lagrangian value and $\tilde{g}^{k}$ are the levels of
constraint violations (“surrogate” sub-gradients). These values are used in
place of dual values $q(\lambda^{k})$ and sub-gradients $g^{k}$. The
convergence to $\lambda^{*}$ is guaranteed [73]222Unlike that in Polyak’s
formula [70], parameter $\gamma$ is less than 1 rather than 2.. The
concomitant reduction of multiplier zigzagging has been also observed. Even
though the geometric convergence rate is expected, e.g., multipliers strictly
move closer to $\lambda^{*}$, the implementation of the method is problematic
because of the unavailability of the knowledge about the optimal dual value
$q(\lambda^{*})$.

Surrogate “Level-Based” Lagrangian Relaxation (SLBLR). To overcome the above
issue, the SLBLR method has been recently developed by [74] to obtain “level”
(over)estimates of the optimal dual value $q(\lambda^{*})$ within the Polyak’s
step-sizing formula (1) to exploit the geometric convergence potential.

The main idea of SLBLR is: if multipliers do not approach $\lambda^{*}$, then
stepsizes violate (1) as proved in [74]. Specifically, if the following
“multiplier-divergence-detection” problem (with $\lambda$ being a continuous
decision variable: $\lambda\in\mathbb{R}^{M}$)

$\begin{cases}\|\lambda-\lambda^{k_{j}+1}\|_{2}\leq\|\lambda-\lambda^{k_{j}}\|_{2},\\\
\|\lambda-\lambda^{k_{j}+2}\|_{2}\leq\|\lambda-\lambda^{k_{j}+1}\|_{2},\\\
\qquad\qquad\vdots\\\
\|\lambda-\lambda^{k_{j}+n_{j}}\|_{2}\leq\|\lambda-\lambda^{k_{j}+n_{j}-1}\|_{2},\end{cases}$
(2)

admits no feasible solution with respect to $\lambda$ for some $k_{j}$ and
$n_{j}$, then there exists $i\in[k_{j},k_{j}+n_{j}]$ such that

$s^{i}\geq\gamma\cdot\frac{q(\lambda^{*})-L^{k}}{\|\tilde{g}^{k}\|^{2}_{2}}.$
(3)

From (3) it follows that there exists an overestimate of $q(\lambda^{*})$, and
the overestimate (the “level value”) has been derived to be:

$\overline{q}_{j}=\max_{i\in[k_{j},k_{j}+n_{j}]}\left(\frac{1}{\gamma}\cdot
s^{i}\cdot\|\tilde{g}^{k}\|^{2}_{2}+L^{k}\right)>q(\lambda^{*}).$ (4)

As a result, the step-sizes are set by using a series of decreasing (not
necessarily monotonically) overestimates $\\{\overline{q}_{k}\\}$ as:

$s^{k}=\zeta\cdot\gamma\cdot\frac{\overline{q}_{k}-L^{k}}{\left\|\tilde{g}^{k}\right\|^{2}_{2}},\gamma<1,\zeta<1.$
(5)

The factor of $\gamma$ in the above formula is inherited from the formula (1),
and the extra factor of $\zeta$ is introduced to counter the fact that the
stepsizes are set by using overestimates of the optimal dual value
$q(\lambda^{*})$.

The rest of the paper is organized as follows. After presenting stochastic JSS
formulation in Section II, the SLBLR method is used to solve the problem in
Section III. The numerical case studies and solution results are presented in
Section IV demonstrates the advantages of the SLBLR method. The conclusion is
provided in Section V.

##  II Job-Shop Scheduling Problem Description and Formulation

Within a typical job shop, machines are grouped by their functionality (e.g.,
milling, lathing, drilling, etc.), and the number of such groups is denoted by
$M$. Within each group, the number of machines referred to as the
“$\underline{\text{c}}$apacity” is denoted as $C_{m},$ where $m=1,\dots,M.$
For every order $i=1,\dots,I$ (e.g., a part to be processed) arriving at a job
shop with priority $w_{i}$ a specific sequence of $J_{i}$ operations,
processing times $p_{i,j,m}$ (for each part $i$, operation $j$333Operation $j$
is denoted as $(i,j)$ whenever appropriate to indicate specific part $i$
requiring operation $j$. and machine type $m$), as well as an integer-valued
due date $d_{i}$ is specified. Processing times are assumed to be an integer
and each operation can occupy several contiguous time blocks, that is, every
operation is assumed to be non-preemptive. Accordingly, the time horizon is
assumed to consist of $T$ discrete equidistant time blocks. A part $i$ may
need to go through a sequence of operations, each is eligible to be performed
by one or a few machine groups. Let $E_{i,j}$ be a group of machines
$\underline{\text{e}}$ligible to process operation $(i,j)$, and $O_{m}$ be a
set of $\underline{\text{o}}$perations that can be processed on machine $m$.

Figure 1: A schematic of job flow within a typical job shop. Circles represent
inputs for each job and squares represent the corresponding products as
outputs. Machine groups may contain 1 or more machines.

Figure 1 illustrates possible, albeit simplified, job flow within a job shop
omitting the time aspect. As shown in Figure 1, Job 1 requires 3 operations to
be performed by Machine Groups 1, 6, and 3, sequentially. Likewise, Job 2
requires 3 operations to be performed by Machine Groups 2, 3, and 5,
sequentially; and Job 3 requires 4 operations to be performed by Machine
Groups 1, 2, 5, and 4, sequentially. If Machine Group 1 contains only one
machine available and Jobs 1 and 3 arrive at the same time, only one job can
be processed during the next time period — this is to be decided through
optimization. A similar argument holds for other machine groups, except
Machine Groups 4 and 6, which only operates on Job 3 and Job 1, respectively.
In terms of the job shop shown in Figure 1, $E_{1,1}=\\{1\\}$,
$E_{1,2}=\\{6\\}$, etc., and
$O_{1}=\\{(1,1),(3,1)\\},O_{2}=\\{(2,1),(3,2)\\},$ etc.

There is a nonzero probability $p_{(i,j)}^{s}$ of scrap after operation $j$ of
part $i$, and production must restart from Operation $1$. Additionally, there
is a nonzero probability $p_{(i,j)}^{r}$ that a scrapped part $i$ can be
recovered through rework and its production must restart from operation $j$.
Transitions among the states of a part follow a Markov process schematically
illustrated in Figure 2. For simplicity of exposition, probabilities of scrap
and rework are chosen to be 20% and 50%, respectively. Transitions from one
operation to the next thus happen with a probability of 80%, and in the case
of a defective part, with equal probabilities of 10%, a part will either be
reworked or completely discarded, triggering the restart of processing from
Operation 1.

Figure 2: A Markov process is illustrated and used in this work to capture
transitions among states (operations) of a part. In this illustration, with a
probability of 10% a part is either scrapped or reworked; at Operation 1,
scrap and rework are lumped together since a defective part needs to begin
getting processed at Operation 1.

###### Assumption 1

The models are developed in this paper under the following assumptions:

  1. 1.

Each machine can process at most one operation at a time;

  2. 2.

Preemption is not allowed, so no operating can be interrupted during
processing;

  3. 3.

No two operations of a given job can be processed at the same time, and
moreover, the predefined sequence of the operations should be kept: no two
operations can be swapped;

  4. 4.

All jobs and machines are available at Time 1;

  5. 5.

All the non-processing related times such as setup times and transportation
times are included in the processing time.

The job-shop scheduler’s task is to proactively assign each operation of each
order to specifically designated machines in anticipation of scrap and rework
while minimizing the overall expected tardiness. Accordingly, the overall
layout of the proactive scheduling consists of two attempts — if scrap happens
on a first attempt, the second attempt is initiated. These two attempts are
captured probabilistically and are modeled within one problem formulation. The
formulation developed by [40] is used as a basis, and stochastic elements due
to scrap and rework are added, as explained in the following sections. Once a
feasible schedule is obtained, the schedule corresponding to the first attempt
is generated at the beginning of a shift: which job is to be assigned to what
machine group and at what time. The integer-valued beginning and completion
times for all operations within both attempts ((1) and (2)) are collectively
denoted as vectors $b=\\{\\{b^{1}\\},\\{b^{2}\\}\\}$ and
$c=\\{\\{c^{1}\\},\\{c^{2}\\}\\}$. After a schedule is obtained, a sequence of
operations to be processed by each machine group is inferred from the
beginning times of the “first attempt” $\\{b^{1}\\}$ and the resulting
sequence is passed on to the job-shop floor. After each shift, if there is no
scrap, job-shop operations follow the aforementioned sequence. Otherwise, a
rescheduling is triggered and a new sequence of operations is inferred. The
short-solving times and high-quality schedules are thus especially important
in the presence of uncertainties.

###  II-A Objective Function

To minimize the delays, it is important to ensure the on-time shipment of
processed orders. While it may not be possible to completely avoid delays,
minimizing weighted tardiness can help mitigate their impact. Weighted
tardiness is used as the objective function in the optimization problem, with
the goal of minimizing the expected tardiness of all orders

$\displaystyle
o(c)=\sum_{i=1}^{I}w_{i}\cdot\Bigg{\\{}\prod_{j=1}^{J_{i}}(1-p^{s}_{(i,j)})\Bigg{\\}}\cdot\max(c_{i,J_{i}}^{1}-d_{i},0)+$
(6)
$\displaystyle\sum_{i=1}^{I}w_{i}\\!\cdot\\!\sum_{j^{\prime}=1}^{J_{i}}\Bigg{\\{}\\!\Bigg{\\{}\\!\prod_{j=1}^{j^{\prime}-1}(1-p^{s}_{(i,j)})-\prod_{j=1}^{j^{\prime}}(1-p^{s}_{(i,j)})\Bigg{\\}}\\!\times\\!$
$\displaystyle(1-p^{r}_{(i,j^{\prime})})\\!\cdot\\!\max(c_{i,J_{i},j^{\prime},\bf{0}}^{2}-d_{i},0)\\!\Bigg{\\}}+$
$\displaystyle\sum_{i=1}^{I}w_{i}\\!\cdot\\!\sum_{j^{\prime}=1}^{J_{i}}\Bigg{\\{}\\!\Bigg{\\{}\\!\prod_{j=1}^{j^{\prime}-1}(1-p^{s}_{(i,j)})-\prod_{j=1}^{j^{\prime}}(1-p^{s}_{(i,j)})\Bigg{\\}}\\!\times\\!$
$\displaystyle(p^{r}_{(i,j^{\prime})})\\!\cdot\\!\max(c_{i,J_{i},j^{\prime},\bf{1}}^{2}-d_{i},0)\\!\Bigg{\\}}.$

The linearization of max operators follows the standard special ordered set
procedure [38, 39]. Despite the complex appearance of the objective function
within (6), the important salient feature of it is additivity in terms of jobs
$i:o(c)=\sum_{i=1}^{I}o_{i}(c_{i}),$ which is efficiently exploited in Section
III.

###  II-B Constraints

To ensure the feasibility of the schedule, processing-time-requirement
constraints (how long an operation needs to be processed), operation
precedence constraints (to keep the order of operations as specified within a
job), and machine capacity constraints (the expected number of operations to
be processed on machines of a certain type cannot exceed the number of
machines) are imposed at a modeling stage.

Beginning/Completion Time Constraints: During Attempt $1,$ completion time
$c_{i,j}^{1}$ of operation $(i,j)$ is expressed as:444Due to the discrete
nature of time, one unit is subtracted. For instance, if an operation starts
at 4 (e.g., 4:00 pm) with a processing time of 2, it finishes at 5 (e.g.,
5:59:59 pm). If the next operation starts after the previous one’s completion,
it begins at time 6, leaving 2 units of time between the start times of
subsequent operations.

$c_{i,j}^{1}=\sum_{m\in E_{i,j}}\sum_{t=1}^{T}(t+p_{i,j,m})\cdot
x_{i,j,m,t}^{1}-1,\forall(i,j),$ (7)

where $x_{i,j,m,t}^{1}$ is an indicator: $x_{i,j,m,\hat{t}}^{1}=1,$ if an
operation $(i,j)$ begins at time $\hat{t}$ 555Here $\hat{t}$ is used to
distinguish a specific beginning time from $t,$ which is a dummy variable
within summations. on a machine group $m,$ and $0,$ otherwise; the above
indicator variable can take the value of 1 for only one combination of
$\\{i,j,m,t\\},$ therefore, the following constraint is imposed to satisfy
Assumptions 1.1 and 1.3:

$\sum_{m\in E_{i,j}}\sum_{t=1}^{T}x_{i,j,m,t}^{1}=1,\forall(i,j).$ (8)

Therefore, all the terms within (7) with the exception of the $\hat{t}^{\rm
th}$ term are zero, and the part is completed at a time $\hat{t}+p_{i,j,m}.$
Accordingly, during Attempt 1, the beginning time $b_{i,j}^{1}$ of operation
$(i,j)$ is expressed as:

$b_{i,j}^{1}=\sum_{m\in E_{i,j}}\sum_{t=1}^{T}t\cdot
x_{i,j,m,t}^{1},\forall(i,j).$ (9)

The completion time within the second attempt is modeled in a similar way:

$c_{i,j,j^{\prime},r}^{2}\\!=\\!\\!\\!\\!\sum_{m\in
E_{i,j}}\\!\\!\sum_{t=1}^{T}(t+p_{i,j,m})\\!\cdot\\!x_{i,j,j^{\prime},m,t,r}^{2}\\!-\\!1,\\!\forall(i,j,j^{\prime},r).$
(10)

The additional subscripts $j^{\prime}$ and $r$ indicate that the completion
time is computed for scenarios whereby a part is scrapped after operation
$j^{\prime}$ within the first attempt, and then the part is either reworked
($r=1$) or completely discarded ($r=0$). Analogously to (8) and (9), the
beginning times of the second attempt are formulated as follows:

$\displaystyle\sum_{m\in
E_{i,j}}\sum_{t=1}^{T}x_{i,j,j^{\prime},m,t,r}^{2}=1,\forall(i,j,j^{\prime},r),$
(11) $\displaystyle b_{i,j,j^{\prime},r}^{2}=\sum_{m\in
E_{i,j}}\sum_{t=1}^{T}t\cdot
x_{i,j,j^{\prime},m,t,r}^{2},\forall(i,j,j^{\prime},r).$ (12)

Operation Precedence: Within each attempt, operation precedence constraints
ensure that subsequent operation $j+1$ can only begin after the previous
operation $j$ is completed:

$b_{i,j+1}^{1}\geq c_{i,j}^{1}+1,\;\;b_{i,j+1,j^{\prime},r}^{2}\geq
c_{i,j,j^{\prime},r}^{2}+1,\forall(i,j,j^{\prime},r).$ (13)

Scrap/Rework Constraints: Scrap/rework constraints ensure that the second
attempt begins at Operation $1$ but only after operation $(i,j)$ leads to
scrap ($r=0$) of part $i$:

$b_{i,1,j^{\prime},0}^{2}\geq c_{i,j^{\prime}}^{1}+1,\forall(i,j^{\prime}).$
(14)

If a part can be reworked ($r=1$), then the operation $j^{\prime}$ (second
subscript) is repeated:

$b_{i,j^{\prime},j^{\prime},1}^{2}\geq
c_{i,j^{\prime}}^{1}+1,\forall(i,j^{\prime}).$ (15)

Beginning of Shift Constraints: It is a common practice to reschedule
operations after the end of a shift to generate the schedule for the upcoming
shifts. After part $i$ is scrapped after operation $j^{\prime},$ the second
attempt needs to be initiated at or after the beginning of future shifts. This
is captured through the following constraints:

$\displaystyle
S\cdot\Big{\lceil}{\frac{c_{i,j^{\prime}}^{1}}{S}\Big{\rceil}}+1$
$\displaystyle\leq
b_{i,j^{\prime},j^{\prime},0}^{2},S\cdot\Big{\lceil}{\frac{c_{i,j^{\prime}}^{1}}{S}\Big{\rceil}}+1\leq
b_{i,1,j^{\prime},1}^{2},\forall(i,j^{\prime}),$ (16)

where $S$ is the length of a shift. The ceiling operator is linearized after
introducing integer variables $y_{i,j^{\prime}}$ as:

$\frac{c_{i,j^{\prime}}^{1}}{S}\leq
y_{i,j^{\prime}}\leq\frac{c_{i,j^{\prime}}^{1}}{S}+1-\varepsilon,$ (17)

where $\varepsilon$ is a small positive number.

Expected Machine Capacity Constraints: Machine capacity constraints ensure
that at any point in time, the expected number of operations processed does
not exceed the total number of machines within a group eligible to perform the
corresponding operation:

$\displaystyle\sum_{(i,j)\in
O_{m}}\prod_{j^{\prime}=1}^{j-1}(1-p^{s}_{{(i,j^{\prime})}})\cdot\left\\{\sum_{t^{\prime}=t-p_{i,j,m,1}:t^{\prime}\geq
1}^{t}x_{i,j,m,t^{\prime}}^{1}\right\\}+$ (18)
$\displaystyle\sum_{j^{\prime},(i,j)\in
O_{m}}\\!\\!\\!\\!\Bigg{(}\\!\prod_{j=1}^{j^{\prime}-1}(1\\!-\\!p^{s}_{{(i,j^{\prime})}})-\prod_{j=1}^{j^{\prime}}(1\\!-\\!p^{s}_{(i,j^{\prime})})\\!\\!\Bigg{)}\\!\times\\!$
$\displaystyle(1-p^{r}_{{(i,j^{\prime})}})\\!\cdot\\!\left\\{\sum_{t^{\prime}=t-p_{i,j,m,1}:t^{\prime}\geq
1}^{t}\\!\\!\\!\\!\\!x_{i,j,j^{\prime},m,t^{\prime},0}^{2}\\!\right\\}+$
$\displaystyle\sum_{j^{\prime},(i,j)\in
O_{m}}\\!\\!\\!\\!\Bigg{(}\prod_{j=1}^{j^{\prime}-1}(1\\!-\\!p^{s}_{{(i,j^{\prime})}})-\prod_{j=1}^{j^{\prime}}(1\\!-\\!p^{s}_{(i,j^{\prime})})\\!\\!\Bigg{)}\\!\times\\!$
$\displaystyle
p^{r}_{{(i,j^{\prime})}}\\!\cdot\\!\left\\{\sum_{t^{\prime}=t-p_{i,j,m,1}:t^{\prime}\geq
1}^{t}\\!\\!\\!\\!\\!x_{i,j,j^{\prime},m,t^{\prime},1}^{2}\\!\right\\}\\!\leq\\!C_{m},\forall(m,t).$

Equation (18) ensures that at any time $t$ and for any machine group $m$, the
expected number of parts processed does not exceed the number of machines
$C_{m}$, and the expectation within (18) follows the similar logic as within
the objective function (6). The exception is the first term, thereby the upper
limit of the summation is $j-1.$ For example, the probability of a part
surviving up to the last operation is
$\prod_{j^{\prime}=1}^{J_{i}-1}(1-p^{s}_{{(i,j^{\prime})}})$, while the
probability of a part surviving after the last operation is complete is
$\prod_{j^{\prime}=1}^{J_{i}}(1-p^{s}_{{(i,j^{\prime})}}),$ which is
appropriate for calculating the tardiness as in (6).

For further convenience and compactness of notation, the objective function is
denoted as $o(c)=\sum_{i=1}^{I}o_{i}(c_{i})$ due to the additivity of (6) and
the optimization problem becomes:

$\min_{b,c,x}\left\\{\sum_{i=1}^{I}o_{i}(c_{i}),\;\;{\rm
s.t.}\;\;\eqref{eq2}-\eqref{eq11}\right\\}.$ (19)

The constraints (18) are additive in terms of parts $i$ and are converted to
equality constraints through the use of nonnegative slack decision variables
$z_{m,t}$, and are expressed as:

$\sum_{i=1}^{I}g_{i}(x_{i})+z_{m,t}-C_{m}=0.$ (20)

Here, $\sum_{i=1}^{I}g_{i}(x_{i})$ compactly denotes the left-hand side of
(18) with
$x_{i}=\\{\\{x_{i,j,m,t}^{1}\\},\\{x_{i,j,j^{\prime},m,t,r}^{2}\\}\\}$
denoting a vector containing all decisions with respect to operations and
attempts. While all the variables within (6)-(18) are integer, because of the
probabilities involved within (18), the slack variables are modeled as
continuous variables. In the following section, the solution methodology is
presented.

##  III Solution Methodology

Lagrangian Relaxation relies on the optimization of the dual function, which,
in a general form can be written as:

$\displaystyle\max_{\lambda}\\{q_{\rho}(\lambda):\lambda\in\mathbb{R}^{M}\times\mathbb{R}^{T}\\},$
(21)

with

$\displaystyle q_{\rho}(\lambda)=$ (22)
$\displaystyle\min_{b,c,x,z}\Bigg{\\{}L_{\rho}(c,x,z,\lambda),(b_{i},c_{i},x_{i})\in\mathcal{F}_{i},i=1,\dots,I,z\geq
0.\Bigg{\\}},$

where

$\displaystyle L_{\rho}(c,x,z,\lambda)\equiv$
$\displaystyle\sum_{i=1}^{I}o_{i}(c_{i})\\!+\\!\lambda^{\rm
T}\\!\cdot\\!\left(\sum_{i=1}^{I}g_{i}(x_{i})+z_{m,t}-C_{m}\right)$
$\displaystyle\\!+\\!\rho\\!\cdot\\!\left\|\sum_{i=1}^{I}g_{i}(x_{i})+z_{m,t}-C_{m}\right\|_{1}\\!,$

is the “absolute-value” Lagrangian function [72]. The feasible set
$\mathcal{F}_{i}$ for each job $i$ is delineated by constraints (7)-(15).
Lagrangian multipliers $\lambda$ (“dual” variables) are the decision variables
with respect to the dual problem (21) and it is assumed that the set of
optimal solutions
$\Omega=\\{\lambda\in\mathbb{R}^{M}\times\mathbb{R}^{T}\;\;|\;\;q(\lambda)=q(\lambda^{*})\\}$
is not empty. The minimization within (22) with respect to $\\{b,c,x,z\\}$ is
referred to as the “relaxed problem”. Due to integer variables
$\\{b,c,x,z\\}$, the dual function (21) is non-smooth with facets (each
corresponding to a particular solution to the relaxed problem within (22))
intersecting at ridges whereby derivatives of $q(\lambda)$ exhibit
discontinuities; in particular, the dual function is not differentiable at
$\lambda^{*}$.

To maximize the dual function, “Surrogate Level-Based Lagrangian Relaxation”
[74] is chosen. The dual function (21) is maximized by updating Lagrange
multipliers $\lambda$ by taking a series of steps $s^{k}$ along “surrogate”
sub-gradient directions
$\sum_{i=1}^{I}g_{i}(\tilde{x}_{i}^{k})+\tilde{z}_{m,t}^{k}-C_{m}$ (levels of
constraint violation) as

$\lambda^{k+1}=\left[\lambda^{k}+s^{k}\cdot\left(\sum_{i=1}^{I}g_{i}(\tilde{x}_{i}^{k})+\tilde{z}_{m,t}^{k}-C_{m}\right)\right]^{+},$
(23)

where $[\cdot]^{+}$ is a projection operator onto a positive orthant
$\left\\{\lambda|\lambda\geq 0\right\\}$. Following (5), the step-sizes are
set as:

$s^{k}=\\!\zeta\\!\cdot\\!\gamma\\!\cdot\\!\frac{\overline{q}_{k}-L(\tilde{b}^{k},\lambda^{k})}{\left\|\sum_{i=1}^{I}g_{i}(\tilde{b}_{i}^{k})\\!+\\!\tilde{s}_{m,t}^{k}\\!-\\!M_{m}\right\|^{2}_{2}},\gamma<1,\zeta<1.$
(24)

Here the “$\sim$” is used to distinguish solutions $\\{b^{k},x^{k},z^{k}\\}$
that are obtained by solving optimally the relaxed problem (minimization
within (22)) from subproblem solution
$\\{\tilde{b}^{k},\tilde{x}^{k},\tilde{z}^{k}\\}$. Subproblems are formulated
by forming a group of parts $\mathcal{I}$ and optimizing the relaxed problem
with respect to the associated variables while keeping decision variables not
belonging to a subset $\mathcal{I}$ fixed as:

$\displaystyle\min_{b_{i},c_{i},x_{i},z;i\in\mathcal{I}}$
$\displaystyle\left\\{\sum_{i\in\mathcal{I}}o_{i}(c_{i})+\lambda^{\rm
T}\cdot\left(\sum_{i\in\mathcal{I}}g_{i}(x_{i})+z_{m,t}\right)+\right.$
$\displaystyle\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!+\left.\rho\cdot\left\|\sum_{i\in\mathcal{I}}g_{i}(x_{i})+\sum_{i\notin\mathcal{I}}g_{i}(x_{i}^{k-1})+z_{m,t}-C_{m}\right\|_{1}\right\\}.$
(25)

The above minimization involves piecewise linear penalties ($l_{1}$ norms),
that efficiently penalize constraint violations and are exactly linearizable
through the use of special ordered sets, thereby enabling the use of MILP
solvers.

Heuristics. Heuristics are generally the last, yet, necessary step required to
obtain feasible solutions. Feasible solutions are obtained by repairing
subproblem solutions and the feasible solution search is initiated after the
norm of constraint violations is “sufficiently small.” For the JSS problem,
the following simple rule along the lines of a “local search” is used: adjust
subproblem beginning times of operations $(i,j)$ by no more than $\delta_{j}$
time periods. This rule is operationalized by solving the original problem (to
enforce feasibility) subject to extra restrictions to model the elements of
the “local search” as follows:

$\min_{b,c,x,z}\;\;\eqref{eq1},\;\;{\rm
s.t.}\;\;\eqref{eq2}-\eqref{eq11},\;\;-\delta_{j}\leq
b_{i,j}-b_{i,j}^{k}\leq\delta_{j}.$ (26)

The above problem is then solved by a MILP solver. The steps of the SLBLR
method are summarized in the algorithm below.

Surrogate “Level-Based” Lagrangian Relaxation

  * Step 1: Initialization. Initialize $\lambda^{0},$ $s^{0},$ $c^{0}.$

  * Step 2: Subproblem Solving. Solve a subproblem (25).

  * Step 3: Stepsize Update. Use (24) to update stepsizes.

  * Step 4: Multiplier Update. Use (23) to update multipliers.

  * Step 5: Penalty Coefficient Update. Update $\rho^{k+1}=\rho^{k}+\beta,\beta>0.$666When solving continuous problems, penalty coefficients are typically updated by using a multiplicative constant $\beta$ as $\rho^{k+1}=\rho^{k}\cdot\beta,\beta>1$, for example, within the Method of Multipliers [75]. However, an $additive$ constant $\beta$ is adopted here.

  * Step 6: Constraint Violation Check. Check the Levels of Constraint Violation. If $\left\|\sum_{i=1}^{I}g_{i}(x_{i})+z_{m,t}-C_{m}\right\|\leq\epsilon$, go to Step 7, otherwise, go to Step 2.

  * Step 7: Feasible Solution Search. Solve (26). If the CPU time limit is not reached, go to 2, otherwise, terminate with the best feasible solution obtained.

The efficiency of the overall approach is discussed with respect to numerical
case studies in the following section.

##  IV Numerical Testing Results

The solution methodology was implemented in an open-source Julia package,
Jobshop.jl available at https://github.com/PSORLab/Jobshop.jl [76]. Each
subproblem in the numerical examples was solved using a 64-bit CPLEX 12.10.0
optimizer and tested using a CPU (48 threads) of Intel(R) Xeon(R) E-2286M CPU
@ 2.99 Hz, 192 GB of RAM, and Windows 10. Example 1 is to demonstrate the
computational efficiency of the method, it is an instance with 20 jobs, 5
operations per job, and 5 machines, each specializing in a specific operation.
Example 2 is to demonstrate the scalability of the SLBLR method with respect
to the increase in the number of jobs as well as in duration of the processing
times.

###  IV-A Example 1. Instance with 20 Jobs and Short Processing Times.

Within this example, a base-case instance with 20 jobs, 5 operations per job,
and 5 machines, each designated to perform one specific operation is
considered. Processing times are generated randomly based on discrete uniform
distribution $U[1,5]$. The data used in this example are shown in Table I in
the Appendix.

The due dates are $d=\\{$15, 25, 32, 36, 21, 27,26, 13, 29, 12, 35, 31, 19,
24, 33, 23, 18, 21, 17, 22$\\}.$ The machine capacities are
$M=\\{2,3,2,2,3\\}.$ The scrap rate is 5% and the rework rate is 20%. The
priorities are all chosen to be equal $w_{i}=1,i=1,\dots,I$. The results of
the SLBLR method and standard B&C are shown in Figure 3.

Figure 3: The results for Example 1 are illustrated. The SLBLR method performs
approximately 250 times faster than the B&C algorithm.

As illustrated in Figure 3, the SLBLR method is more than two orders of
magnitude faster than standard B&C.

To test robustness, the base case is modified by re-generating processing
times based on the discrete uniform distribution $U[1,5]$ as well as by
generating due dates based on discrete uniform distribution $U[10,40].$
Accordingly, 5 new instances are generated. Because of the complexity of the
problem, the stopping criterion is 10% of the duality gap. The results are
demonstrated in Table I.

TABLE I: Robustness results for Example 1.

| B&C | SLBLR |
---|---|---|---
Case | Feas. | Lower | Gap (%) | CPU | Feas. | Gap (%) | CPU | Improvement
| Cost | Bound |  | time (sec) | Cost |  | time (sec) |
1 | 133.50 | 120.27 | 9.91% | 12189.88 | 131.45 | 8.51% | 401.20 | 30.38
2 | 158.89 | 143.25 | 9.84% | 78690.8 | 158.82 | 9.80% | 313.17 | 251.26
3 | 137.05 | 123.34 | 10.00% | 78788.47 | 135.31 | 8.85% | 843.47 | 93.40
4 | 150.31 | 135.31 | 9.98% | 6610.47 | 149.83 | 9.69% | 322.04 | 20.52
5 | 118.52 | 106.98 | 9.74% | 41301.77 | 118.28 | 9.56% | 404.07 | 102.21
|  |  | Avg. time: | 43516.27 |  | Avg. time: | 456.79 | 95.26
|  |  | St. Dev.: | 34747.87 |  | St. Dev.: | 220.32 |

Table I demonstrates that the CPU time improvement ranges from 20.52 times to
251.26 times, for an average of 95.26 times, which is almost two orders of
magnitude improvement.

###  IV-B Example 2. Instances with 20 and 100 Jobs and Long Processing Times.

Results for a 20-Job Case. The base case study from Example 1 is modified by
increasing due dates by a factor of 10 and by generating processing times by
using uniform discrete distribution $U[1,50].$ The results are shown in Figure
4.

Figure 4: The results for Example 2 are illustrated. The SLBLR method performs
approximately 4 times faster than the B&C algorithm and 2.59 times more
accurately.

Because of the associated increase of the time periods required to accommodate
much-increased processing times, the CPLEX is out of memory after 10,255
seconds. In contrast, the SLBLR method is 4 times faster and 2.59 times more
accurate.

Results for a 100-Job Case. The instance is created by generating processing
times based on the discrete uniform distribution $U[1,50]$ as well as by
generating due dates based on discrete uniform distribution $U[100,400].$ For
this case, no feasible solution was obtained by CPLEX after 15,129.42 sec. In
contrast, a feasible solution with a corresponding cost of 9697.59 is found by
SLBLR after 3618.53 seconds.

##  V Conclusions, Limitations, and Future Work.

This study presented a proactive scheduling approach that minimizes the
expected tardiness in the presence of scrap and rework in job-shop
manufacturing environments. The major contributions and takeaways of this
study can be summarized as follows:

  1. 1.

Formulation. The use of Markov processes to model probabilistic defects in
manufacturing production, and the formulation of an optimization problem that
minimizes expected tardiness.

  2. 2.

Online Capabilities of SLBLR. Drawing on the geometric potential inherited
from the Polyak formula, the SLBLR method can effectively operate online and
accommodate unexpected or urgent orders by continuously updating multipliers.
This capability enhances the resilience of manufacturing operations.

  3. 3.

Consistent Improvement Capabilities. The SLBLR method’s ability to
systematically improve solutions through the update of multipliers.

  4. 4.

Optimization Enhancement through QC and ML Integration. The SLBLR method’s
adaptability to different optimization techniques facilitates seamless
integration with quantum computing (QC) and machine learning (ML). This
combination holds the potential to substantially improve subproblem-solving
and feasible solution quality. Recent advancements have shown promising
results in tackling job-shop scheduling challenges within the automation
community, highlighting the value of merging Lagrangian Relaxation and Machine
Learning methods [77].

  5. 5.

Broader Impact. The potential for significant reductions in tardiness and
improved downstream supply chain operations management.

While the results and future directions are promising, some limitations should
be acknowledged. First, the computational effort can be considerable for
large-scale problems, as seen in the 100-job case. Second, the model assumes
all jobs are available at the beginning of each shift, which may not always be
accurate. Lastly, the model does not account for sequence-dependent setup
times, which can significantly affect manufacturing operations.

Future research can focus on overcoming these limitations by developing more
efficient algorithms for large-scale problems, incorporating real-time job
availability, and considering setup times in the optimization formulation.
Additionally, exploring the integration of the SLBLR method with emerging
technologies like artificial intelligence and Industry 4.0 could lead to more
comprehensive and effective manufacturing operations management.

## Disclaimer

This paper was prepared as an account of work sponsored by an agency of the
United States Government. Neither the United States Government nor any agency
thereof, nor any of their employees, makes any warranty, express or implied,
or assumes any legal liability or responsibility for the accuracy,
completeness, or usefulness of any information, apparatus, product, or process
disclosed or represents that its use would not infringe privately owned
rights. Reference herein to any specific commercial product, process, or
service by trade name, trademark, manufacturer, or otherwise does not
necessarily constitute or imply its endorsement, recommendation, or favoring
by the United States Government or any agency thereof. The views, opinions,
and/or findings contained in this paper are those of the authors and should
not be interpreted as representing the official views or policies, either
expressed or implied, of the United States Government or any agency thereof,
including the Air Force Research Laboratory, the United States Air Force, and
the Department of Defense.

## Funding

This material is based upon work supported by the U.S. Department of Energy’s
Office of Energy Efficiency and Renewable Energy (EERE) under the Advanced
Manufacturing Office Award No. DE-EE0007613. We also gratefully acknowledge
the Air Force Research Laboratory, Materials and Manufacturing Directorate
(AFRL/RXMS) for support via Contract No. FA8650-20-C-5206. This work is also
supported in part by the US NSF under award ECCS-1810108.

## References

  * [1] B.-Y. Cheng, J. Y.-T. Leung, and K. Li, “Integrated scheduling of production and distribution to minimize total cost using an improved ant colony optimization method,” _Computers & Industrial Engineering_, vol. 83, pp. 217–225, 2015. [Online]. Available: https://doi.org/10.1016/j.cie.2015.02.017
  * [2] T.-M. Choi, “Supply chain systems coordination with multiple risk sensitive retail buyers,” _IEEE Transactions on Systems, Man, and Cybernetics: Systems_ , vol. 46, no. 5, pp. 636–645, 2015. [Online]. Available: https://doi.org/10.1109/TSMC.2015.2452894
  * [3] H. F. Rahman, M. N. Janardhanan, L. P. Chuen, and S. G. Ponnambalam, “Flowshop scheduling with sequence dependent setup times and batch delivery in supply chain,” _Computers & Industrial Engineering_, vol. 158, p. 107378, 2021\. [Online]. Available: https://doi.org/10.1016/j.cie.2021.107378
  * [4] J. Chen, G. Q. Huang, H. Luo, and J. Wang, “Synchronisation of production scheduling and shipment in an assembly flowshop,” _International Journal of Production Research_ , vol. 53, no. 9, pp. 2787–2802, 2015. [Online]. Available: https://doi.org/10.1080/00207543.2014.994075
  * [5] H. F. Rahman, R. Sarker, and D. Essam, “A genetic algorithm for permutation flowshop scheduling under practical make-to-order production system,” _AI EDAM_ , vol. 31, no. 1, pp. 87–103, 2017. [Online]. Available: https://doi.org/10.1016/j.cie.2015.08.006
  * [6] K. Wang, H. Luo, F. Liu, and X. Yue, “Permutation flow shop scheduling with batch delivery to multiple customers in supply chains,” _IEEE Transactions on Systems, Man, and Cybernetics: Systems_ , vol. 48, no. 10, pp. 1826–1837, 2017. [Online]. Available: https://doi.org/10.1109/TSMC.2017.2720178
  * [7] O. Ergun, W. J. Hopp, and P. Keskinocak, “A structured overview of insights and opportunities for enhancing supply chain resilience,” _IISE Transactions_ , vol. 55, no. 1, pp. 57–74, 2023. [Online]. Available: https://doi.org/10.1080/24725854.2022.2080892
  * [8] J. P. Wilson, Z. Shen, U. Awasthi, G. M. Bollas, and S. Gupta, “Multi-objective optimization for cost-efficient and resilient machining under tool wear,” _Journal of Advanced Manufacturing and Processing_ , vol. 4, no. 4, p. e10140, 2022.
  * [9] Y. Sun, J. Tu, M. Bragin, and L. Zhang, “A simulation-based integrated virtual testbed for dynamic optimization in smart manufacturing systems,” _Journal of Advanced Manufacturing and Processing_ , vol. 4, no. 4, p. e10141, 2022.
  * [10] M. R. Garey, D. S. Johnson, and R. Sethi, “The complexity of flowshop and jobshop scheduling,” _Mathematics of operations research_ , vol. 1, no. 2, pp. 117–129, 1976. [Online]. Available: https://doi.org/10.1287/moor.1.2.117
  * [11] X. Zhu, W. Wang, X. Guo, and L. Shi, “A genetic programming-based evolutionary approach for flexible job shop scheduling with multiple process plans,” in _2020 IEEE 16th International Conference on Automation Science and Engineering (CASE)_. IEEE, 2020, pp. 49–54. [Online]. Available: https://doi.org/10.1109/CASE48305.2020.9216783
  * [12] J. B. Sidney, “Optimal single-machine scheduling with earliness and tardiness penalties,” _Operations Research_ , vol. 25, no. 1, pp. 62–69, 1977. [Online]. Available: https://doi.org/10.1287/opre.25.1.62
  * [13] S. Lakshminarayan, R. Lakshmanan, R. L. Papineau, and R. Rochette, “Optimal single-machine scheduling with earliness and tardiness penalties,” _Operations Research_ , vol. 26, no. 6, pp. 1079–1082, 1978. [Online]. Available: https://doi.org/10.1287/opre.26.6.1079
  * [14] A. P. J. Vepsalainen and T. E. Morton, “Priority rules for job shops with weighted tardiness costs,” _Management Science_ , vol. 33, no. 8, pp. 1035–1047, 1987. [Online]. Available: https://doi.org/10.1287/mnsc.33.8.1035
  * [15] G. Chryssolouris, K. Wright, J. Pierce, and W. Cobb, “Manufacturing systems operation: dispatch rules versus intelligent control,” _Robotics and Computer-Integrated Manufacturing_ , vol. 4, no. 3-4, pp. 531–544, 1988. [Online]. Available: https://doi.org/10.1016/0736-5845(88)90026-9
  * [16] E. J. Anderson and J. C. Nyirenda, “Two new rules to minimize tardiness in a job shop,” _The International Journal of Production Research_ , vol. 28, no. 12, pp. 2277–2292, 1990. [Online]. Available: https://doi.org/10.1080/00207549008942866
  * [17] J. Du and J. Y.-T. Leung, “Minimizing total tardiness on one machine is np-hard,” _Mathematics of Operations Research_ , vol. 15, no. 3, pp. 483–495, 1990. [Online]. Available: https://doi.org/10.1287/moor.15.3.483
  * [18] G. Chryssolouris, J. Pierce, and K. Dicke, “An approach for allocating manufacturing resources to production tasks,” _Journal of Manufacturing Systems_ , vol. 10, no. 5, pp. 368–382, 1991. [Online]. Available: https://doi.org/10.1016/0278-6125(91)90055-7
  * [19] J. J. Kanet and Z. Zhou, “A decision theory approach to priority dispatching for job shop scheduling,” _Production and Operations Management_ , vol. 2, no. 1, pp. 2–14, 1993. [Online]. Available: https://doi.org/10.1111/j.1937-5956.1993.tb00035.x
  * [20] T. R. Rohleder and G. Scudder, “A comparison of order-release and dispatch rules for the dynamic weighted early/tardy problem,” _Production and Operations Management_ , vol. 2, no. 3, pp. 221–238, 1993. [Online]. Available: https://doi.org/10.1111/j.1937-5956.1993.tb00099.x
  * [21] A. Seidmann, S. S. Panwalkar, and M. L. Smith, “Optimal assignment of due-dates for a single processor scheduling problem,” _The International Journal Of Production Research_ , vol. 19, no. 4, pp. 393–399, 1981\. [Online]. Available: https://doi.org/10.1080/00207548108956667
  * [22] S. S. Panwalkar, M. L. Smith, and A. Seidmann, “Common due date assignment to minimize total penalty for the one machine scheduling problem,” _Operations Research_ , vol. 30, no. 2, pp. 391–399, 1982. [Online]. Available: https://doi.org/10.1287/opre.30.2.391
  * [23] P. S. Ow and T. E. Morton, “The single machine early/tardy problem,” _Management science_ , vol. 35, no. 2, pp. 177–191, 1989. [Online]. Available: https://doi.org/10.1287/mnsc.35.2.177
  * [24] S. Barman and R. L. Laforge, “The impact of simple priority rule combinations on delivery speed and delivery reliability in a hybrid shop,” _Production and Operations Management_ , vol. 7, no. 4, pp. 402–416, 1998\. [Online]. Available: https://doi.org/10.1111/j.1937-5956.1998.tb00132.x
  * [25] P. B. Luh, D. Chen, and L. S. Thakur, “An effective approach for job-shop scheduling with uncertain processing requirements,” _IEEE Transactions on Robotics and Automation_ , vol. 15, no. 2, pp. 328–339, 1999. [Online]. Available: https://doi.org/10.1109/70.760354
  * [26] D. Golenko-Ginzburg and A. Gonik, “Optimal job-shop scheduling with random operations and cost objectives,” _International Journal of Production Economics_ , vol. 76, no. 2, pp. 147–157, 2002. [Online]. Available: https://doi.org/10.1016/S0925-5273(01)00140-2
  * [27] D.-M. Lei and H.-J. Xiong, “An efficient evolutionary algorithm for multi-objective stochastic job shop scheduling,” in _2007 International Conference on Machine Learning and Cybernetics_ , vol. 2, 2007, pp. 867–872. [Online]. Available: https://doi.org/10.1109/ICMLC.2007.4370264
  * [28] D. Lei, “Simplified multi-objective genetic algorithms for stochastic job shop scheduling,” _Applied Soft Computing_ , vol. 11, no. 8, pp. 4991–4996, 2011\. [Online]. Available: https://doi.org/10.1016/j.asoc.2011.06.001
  * [29] S.-C. Horng, S.-S. Lin, and F.-Y. Yang, “Evolutionary algorithm for stochastic job shop scheduling with random processing time,” _Expert Systems with Applications_ , vol. 39, no. 3, pp. 3603–3610, 2012. [Online]. Available: https://doi.org/10.1016/j.eswa.2011.09.050
  * [30] R. Zhang, S. Song, and C. Wu, “A hybrid differential evolution algorithm for job shop scheduling problems with expected total tardiness criterion,” _Applied Soft Computing_ , vol. 13, no. 3, pp. 1448–1458, 2013. [Online]. Available: https://doi.org/10.1016/j.asoc.2012.02.024
  * [31] H.-a. Yang, Y. Lv, C. Xia, S. Sun, and H. Wang, “Optimal computing budget allocation for ordinal optimization in solving stochastic job shop scheduling problems,” _Mathematical Problems in Engineering_ , vol. 2014, 2014. [Online]. Available: https://doi.org/10.1155/2014/619254
  * [32] J. Shen and Y. Zhu, “Chance-constrained model for uncertain job shop scheduling problem,” _Soft Computing_ , vol. 20, no. 6, pp. 2383–2391, 2016\. [Online]. Available: https://doi.org/10.1007/s00500-015-1647-z
  * [33] A. Jamili, “Job shop scheduling with consideration of floating breaking times under uncertainty,” _Engineering applications of artificial intelligence_ , vol. 78, pp. 28–36, 2019. [Online]. Available: https://doi.org/10.1016/j.engappai.2018.10.007
  * [34] S.-C. Horng and S.-S. Lin, “Apply ordinal optimization to optimize the job-shop scheduling under uncertain processing times,” _Arabian Journal for Science and Engineering_ , pp. 1–13, 2021. [Online]. Available: https://doi.org/10.1007/s13369-021-06317-9
  * [35] S. K. Hasan, R. Sarker, and D. Essam, “Genetic algorithm for job-shop scheduling with machine unavailability and breakdowns,” _International Journal of Production Research_ , vol. 49, no. 16, pp. 4999–5015, 2011. [Online]. Available: https://doi.org/10.1080/00207543.2010.495088
  * [36] Y.-K. Lin, P.-C. Chang, L. C.-L. Yeng, and S.-F. Huang, “Bi-objective optimization for a multistate job-shop production network using nsga-ii and topsis,” _Journal of Manufacturing Systems_ , vol. 52, pp. 43–54, 2019. [Online]. Available: https://doi.org/10.1016/j.jmsy.2019.05.004
  * [37] M. Dawande, S. Gavirneni, and R. Rachamadugu, “Scheduling a two-stage flowshop under makespan constraint,” _Mathematical and computer modelling_ , vol. 44, no. 1-2, pp. 73–84, 2006. [Online]. Available: https://doi.org/10.1016/j.mcm.2004.12.016
  * [38] B. Yan, M. A. Bragin, and P. B. Luh, “Novel formulation and resolution of job-shop scheduling problems,” _IEEE Robotics and Automation Letters_ , vol. 3, no. 4, pp. 3387–3393, 2018. [Online]. Available: https://doi.org/10.1109/LRA.2018.2850056
  * [39] ——, “An innovative formulation tightening approach for job-shop scheduling,” _IEEE Transactions on Automation Science and Engineering_ , 2021\. [Online]. Available: https://doi.org/10.1109/TASE.2021.3088047
  * [40] A. Liu, P. Luh, B. Yan, and M. Bragin, “A novel integer linear programming formulation for job-shop scheduling problems,” _IEEE Robotics and Automation Letters_ , 2021. [Online]. Available: https://doi.org/10.1109/LRA.2021.3086422
  * [41] D. J. Hoitomt, P. B. Luh, and K. R. Pattipati, “A practical approach to job-shop scheduling problems,” _IEEE Transactions on Robotics and Automation_ , vol. 9, no. 1, pp. 1–13, 1993. [Online]. Available: https://doi.org/10.1109/70.210791
  * [42] C. W. Leung, T. N. Wong, K.-L. Mak, and R. Y. K. Fung, “Integrated process planning and scheduling by an agent-based ant colony optimization,” _Computers & Industrial Engineering_, vol. 59, no. 1, pp. 166–180, 2010\. [Online]. Available: https://doi.org/10.1016/j.cie.2009.09.003
  * [43] D. Grimes and E. Hebrard, “Solving variants of the job shop scheduling problem through conflict-directed search,” _INFORMS Journal on Computing_ , vol. 27, no. 2, pp. 268–284, 2015.
  * [44] J. C. Beck, T. Feng, and J.-P. Watson, “Combining constraint programming and local search for job-shop scheduling,” _INFORMS Journal on Computing_ , vol. 23, no. 1, pp. 1–14, 2011.
  * [45] A. Malapert, H. Cambazard, C. Guéret, N. Jussien, A. Langevin, and L.-M. Rousseau, “An optimal constraint programming approach to the open-shop problem,” _INFORMS Journal on Computing_ , vol. 24, no. 2, pp. 228–244, 2012\.
  * [46] Z. Zhao, S. Liu, M. Zhou, X. Guo, and L. Qi, “Decomposition method for new single-machine scheduling problems from steel production systems,” _IEEE Transactions on Automation Science and Engineering_ , vol. 17, no. 3, pp. 1376–1387, 2019.
  * [47] Y. K. Kim, K. Park, and J. Ko, “A symbiotic evolutionary algorithm for the integration of process planning and job shop scheduling,” _Computers & operations research_, vol. 30, no. 8, pp. 1151–1171, 2003. [Online]. Available: https://doi.org/10.1016/S0305-0548(02)00063-1
  * [48] Q. Wu, M. Zhou, Q. Zhu, Y. Xia, and J. Wen, “Moels: Multiobjective evolutionary list scheduling for cloud workflows,” _IEEE Transactions on Automation Science and Engineering_ , vol. 17, no. 1, pp. 166–176, 2019.
  * [49] Z. Zhao, M. Zhou, and S. Liu, “Iterated greedy algorithms for flow-shop scheduling problems: A tutorial,” _IEEE Transactions on Automation Science and Engineering_ , 2021.
  * [50] Z. Zhao, S. Liu, M. Zhou, D. You, and X. Guo, “Heuristic scheduling of batch production processes based on petri nets and iterated greedy algorithms,” _IEEE Transactions on Automation Science and Engineering_ , vol. 19, no. 1, pp. 251–261, 2020.
  * [51] X. Li and L. Gao, “An effective hybrid genetic algorithm and tabu search for flexible job shop scheduling problem,” _International Journal of Production Economics_ , vol. 174, pp. 93–110, 2016. [Online]. Available: https://doi.org/10.1016/j.ijpe.2016.01.016
  * [52] X. Li, L. Gao, Q. Pan, L. Wan, and K.-M. Chao, “An effective hybrid genetic algorithm and variable neighborhood search for integrated process planning and scheduling in a packaging machine workshop,” _IEEE Transactions on Systems, Man, and Cybernetics: Systems_ , vol. 49, no. 10, pp. 1933–1945, 2018\. [Online]. Available: https://doi.org/10.1109/TSMC.2018.2881686
  * [53] H. Li, B. Wang, Y. Yuan, M. Zhou, Y. Fan, and Y. Xia, “Scoring and dynamic hierarchy-based nsga-ii for multiobjective workflow scheduling in the cloud,” _IEEE Transactions on Automation Science and Engineering_ , vol. 19, no. 2, pp. 982–993, 2021.
  * [54] Z. Cao, C. Lin, M. Zhou, C. Zhou, and K. Sedraoui, “Two-stage genetic algorithm for scheduling stochastic unrelated parallel machines in a just-in-time manufacturing context,” _IEEE Transactions on Automation Science and Engineering_ , 2022.
  * [55] Y. Pan, K. Gao, Z. Li, and N. Wu, “Improved meta-heuristics for solving distributed lot-streaming permutation flow shop scheduling problems,” _IEEE Transactions on Automation Science and Engineering_ , 2022.
  * [56] S. Iloglu and L. A. Albert, “An integrated network design and scheduling problem for network recovery and emergency response,” _Operations Research Perspectives_ , vol. 5, pp. 218–231, 2018. [Online]. Available: https://doi.org/10.1016/j.orp.2018.08.001
  * [57] T.-H. Tsai, B. Yan, P. B. Luh, H.-C. Yang, M. A. Bragin, and F.-T. Cheng, “Near-optimal scheduling for ic packaging operations considering processing-time variations and factory practices,” _IEEE Robotics and Automation Letters_ , pp. 1–8, 2023.
  * [58] E. H. Aarts, P. J. van Laarhoven, J. K. Lenstra, and N. L. Ulder, “A computational study of local search algorithms for job shop scheduling,” _ORSA Journal on Computing_ , vol. 6, no. 2, pp. 118–125, 1994.
  * [59] R. J. M. Vaessens, E. H. Aarts, and J. K. Lenstra, “Job shop scheduling by local search,” _Informs Journal on computing_ , vol. 8, no. 3, pp. 302–317, 1996.
  * [60] M. Nouiri, A. Bekrar, A. Jemai, S. Niar, and A. C. Ammari, “An effective and distributed particle swarm optimization algorithm for flexible job-shop scheduling problem,” _Journal of Intelligent Manufacturing_ , vol. 29, no. 3, pp. 603–615, 2018. [Online]. Available: https://doi.org/10.1007/s10845-015-1039-3
  * [61] A. Casalino, A. M. Zanchettin, L. Piroddi, and P. Rocco, “Optimal scheduling of human–robot collaborative assembly operations with time petri nets,” _IEEE Transactions on Automation Science and Engineering_ , vol. 18, no. 1, pp. 70–84, 2019.
  * [62] I.-B. Park, J. Huh, J. Kim, and J. Park, “A reinforcement learning approach to robust scheduling of semiconductor manufacturing facilities,” _IEEE Transactions on Automation Science and Engineering_ , vol. 17, no. 3, pp. 1420–1431, 2019.
  * [63] S. Luo, L. Zhang, and Y. Fan, “Real-time scheduling for dynamic partial-no-wait multiobjective flexible job shop by deep reinforcement learning,” _IEEE Transactions on Automation Science and Engineering_ , vol. 19, no. 4, pp. 3020–3038, 2021.
  * [64] W. Song, N. Mi, Q. Li, J. Zhuang, and Z. Cao, “Stochastic economic lot scheduling via self-attention based deep reinforcement learning,” _IEEE Transactions on Automation Science and Engineering_ , 2023.
  * [65] T. Li, Y. Meng, and L. Tang, “Scheduling of continuous annealing with a multi-objective differential evolution algorithm based on deep reinforcement learning,” _IEEE Transactions on Automation Science and Engineering_ , 2023\.
  * [66] E. D. Taillard, “Parallel taboo search techniques for the job shop scheduling problem,” _ORSA journal on Computing_ , vol. 6, no. 2, pp. 108–117, 1994\.
  * [67] J. W. Herrmann, C.-Y. Lee, and J. Hinchman, “Global job shop scheduling with a genetic algorithm,” _Production and Operations Management_ , vol. 4, no. 1, pp. 30–45, 1995. [Online]. Available: https://doi.org/10.1111/j.1937-5956.1995.tb00039.x
  * [68] R. Tavakkoli-Moghaddam, F. Jolai, F. Vaziri, P. Ahmed, and A. Azaron, “A hybrid method for solving stochastic job shop scheduling problems,” _Applied Mathematics and Computation_ , vol. 170, no. 1, pp. 185–206, 2005\. [Online]. Available: https://doi.org/10.1016/j.amc.2004.11.036
  * [69] A. Ghasemi, A. Ashoori, and C. Heavey, “Evolutionary learning based simulation optimization for stochastic job shop scheduling problems,” _Applied Soft Computing_ , vol. 106, p. 107309, 2021. [Online]. Available: https://doi.org/10.1016/j.asoc.2021.107309
  * [70] B. T. Polyak, “Minimization of unsmooth functionals (in russian),” _USSR Computational Mathematics and Mathematical Physics_ , vol. 9, no. 3, pp. 14–29, 1969. [Online]. Available: https://doi.org/10.1016/0041-5553(69)90061-5
  * [71] M. A. Bragin, P. B. Luh, J. H. Yan, N. Yu, and G. A. Stern, “Convergence of the surrogate lagrangian relaxation method,” _Journal of Optimization Theory and applications_ , vol. 164, no. 1, pp. 173–201, 2015. [Online]. Available: https://doi.org/10.1007/s10957-014-0561-3
  * [72] M. A. Bragin, P. B. Luh, B. Yan, and X. Sun, “A scalable solution methodology for mixed-integer linear programming problems arising in automation,” _IEEE Transactions on Automation Science and Engineering_ , vol. 16, no. 2, pp. 531–541, 2019. [Online]. Available: https://doi.org/10.1109/tase.2018.2835298
  * [73] X. Zhao, P. B. Luh, and J. Wang, “Surrogate gradient algorithm for lagrangian relaxation,” _Journal of optimization Theory and Applications_ , vol. 100, no. 3, pp. 699–712, 1999. [Online]. Available: https://doi.org/10.1023/A:1022646725208
  * [74] M. A. Bragin and E. L. Tucker, “Surrogate “level-based” lagrangian relaxation for mixed-integer linear programming,” _Scientific Reports_ , vol. 12, no. 1, pp. 1–12, 2022. [Online]. Available: https://www.nature.com/articles/s41598-022-26264-1
  * [75] D. P. Bertsekas, _Nonlinear Programming: Third Edition_. Belmont, MA: Athena Scientific, 2016.
  * [76] M. E. Wilhelm, M. Bragin, and M. D. Stuber, “Jobshop.jl,” 2022. [Online]. Available: https://github.com/PSORLab/Jobshop.jl
  * [77] A. Liu, P. B. Luh, K. Sun, M. A. Bragin, and B. Yan, “Integrating machine learning and mathematical optimization for job shop scheduling,” _Conditionally accepted to IEEE Transactions on Automation Science and Engineering_ , 2023.

##  Appendix A Example Data

TABLE I: Processing time $p_{i,j}$ data for job $i$ and operation $j$ for Example 1 are contained in this table. $i$ | $j$ | $p_{i,j}$ | $i$ | $j$ | $p_{i,j}$ | $i$ | $j$ | $p_{i,j}$
---|---|---|---|---|---|---|---|---
1 | 1 | 1 | 2 | 1 | 2 | 3 | 1 | 3
1 | 2 | 2 | 2 | 2 | 1 | 3 | 2 | 2
1 | 3 | 2 | 2 | 3 | 2 | 3 | 3 | 3
1 | 4 | 1 | 2 | 4 | 2 | 3 | 4 | 2
1 | 5 | 1 | 2 | 5 | 4 | 3 | 5 | 2
4 | 1 | 1 | 5 | 1 | 3 | 6 | 1 | 1
4 | 2 | 3 | 5 | 2 | 2 | 6 | 2 | 3
4 | 3 | 2 | 5 | 3 | 5 | 6 | 3 | 2
4 | 4 | 2 | 5 | 4 | 1 | 6 | 4 | 3
4 | 5 | 1 | 5 | 5 | 4 | 6 | 5 | 2
7 | 1 | 2 | 8 | 1 | 1 | 9 | 1 | 1
7 | 2 | 5 | 8 | 2 | 5 | 9 | 1 | 2
7 | 3 | 4 | 8 | 3 | 1 | 9 | 1 | 2
7 | 4 | 1 | 8 | 4 | 5 | 9 | 1 | 5
7 | 5 | 1 | 8 | 5 | 1 | 9 | 1 | 3
10 | 1 | 1 | 11 | 1 | 5 | 12 | 1 | 3
10 | 1 | 3 | 11 | 2 | 2 | 12 | 2 | 1
10 | 1 | 2 | 11 | 3 | 3 | 12 | 3 | 3
10 | 1 | 3 | 11 | 4 | 4 | 12 | 4 | 1
10 | 1 | 1 | 11 | 5 | 5 | 12 | 5 | 1
13 | 1 | 5 | 14 | 1 | 2 | 15 | 1 | 3
13 | 2 | 1 | 14 | 2 | 4 | 15 | 2 | 1
13 | 3 | 3 | 14 | 3 | 2 | 15 | 3 | 1
13 | 4 | 1 | 14 | 4 | 1 | 15 | 4 | 3
13 | 5 | 5 | 14 | 5 | 1 | 15 | 5 | 5
16 | 1 | 3 | 17 | 1 | 3 | 18 | 1 | 1
16 | 2 | 3 | 17 | 2 | 1 | 18 | 2 | 2
16 | 3 | 2 | 17 | 3 | 2 | 18 | 3 | 5
16 | 4 | 5 | 17 | 4 | 3 | 18 | 4 | 2
16 | 5 | 3 | 17 | 5 | 1 | 18 | 5 | 3
19 | 1 | 1 | 20 | 1 | 2 |  |  |
19 | 1 | 4 | 20 | 1 | 4 |  |  |
19 | 1 | 4 | 20 | 1 | 4 |  |  |
19 | 1 | 2 | 20 | 1 | 2 |  |  |
19 | 1 | 5 | 20 | 1 | 1 |  |  |


# Forward-backward correlations between intensive observables

Vladimir Kovalenko  Vladimir Vechernin  Saint Petersburg State University, 7/9
Universitetskaya nab., St. Petersburg, 199034, Russia<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

We demonstrate that the investigations of the forward-backward correlations
between intensive observables enable to obtain more clear signal about the
initial stage of hadronic interaction, e.g. about the process of string
fusion, compared to usual forward-backward multiplicity correlations. As an
example, the correlation between mean-event transverse momenta of charged
particles in separated rapidity intervals is considered. We performed
calculations in the framework of dipole-based Monte Carlo string fusion model.
We obtained the dependence of the correlation strength on the collision
centrality for different initial energies and colliding systems. It is shown
that the dependence reveals the decline of the correlation coefficient for
most central Pb-Pb collisions at LHC energy. We compare the results both with
the ones obtained in alternative models and with the ones obtained by us using
various MC generators.

##  1 Forward-Backward Rapidity Correlations

The study of the correlations between observables in two separated rapidity
windows (the so-called long-range forward-backward correlations) has been
proposed [1] as a signature of the string fusion and percolation phenomenon
[2, 3], which is one of the collectivity effects in ultrarelativistic heavy
ion collisions. Later it was realized [4, 5, 6, 7] that the investigations of
the forward-backward correlations between intensive observables, such e.g. as
mean-event transverse momenta, enable to obtain more clear signal about the
initial stage of hadronic interaction, including the process of string fusion,
compared to usual forward-backward multiplicity correlations. As an example,
the correlation between mean-event transverse momenta of charged particles in
separated rapidity intervals can be considered.

Numerically, the long-range correlations are studied in terms of correlation
functions and correlation coefficients. Correlation function is defined as the
mean value of variable B in the backward window as a function of another
variable F in the forward rapidity window:

$\langle{B}\rangle_{F}=f(F).$

The correlation function can be approximated by the linear function:
$\langle{B}\rangle_{F}=a+b_{BF}F$ using the linear regression. Here $b_{BF}$
is the correlation coefficient, a slope of correlation function

Alternative definition of the correlation coefficient is

$b_{BF}=\frac{\langle FB\rangle-\langle F\rangle\langle B\rangle}{\langle
F^{2}\rangle-\langle F\rangle^{2}}=\frac{\mathrm{cov}(F,B)}{D_{F}}.$ (1)

For a correlation between relative variables, $F/\langle F\rangle$ and
$B/\langle B\rangle$:

$b^{rel}_{FB}=\frac{\langle F\rangle}{\langle B\rangle}b_{FB}\ .$ (2)

For the observables, $B$ and $F$, in backward and forward rapidity windows,
multiplicity ($n$) or mean event transverse momentum
($\text{p}_{\text{t}}$=$\frac{1}{n}\sum\limits_{i=1}^{n}{\text{p}_{\text{t}}}_{i}$)
can be taken [8, 9, 10].

Note that the multiplicities in backward and forward rapidity windows $n_{B}$,
$n_{F}$ are extensive variables, which means that they are proportional to the
size of the system. In contrast, the transverse momentum of particles in
backward and forward rapidity windows $p_{tB}$, $p_{tF}$ are intensive
variables, so they don’t depend on the number of sources which produce
particles.

The long-range multiplicity correlation $b_{nn}$ at large rapidity gap is
dominated by event-by-event variance in the number of sources (cut pomerons,
strings)[11]. In particular, in nuclear collisions, it also reflects the
fluctuation in the number of participant nucleons. The influence of possible
string fusion processes on this type of correlation is rather small [12]. It
also has been shown, that $b_{nn}$ depend strongly on the method of the
centrality determination and the on the width of centrality class [13, 14,
15].

However, the event-by-event fluctuations in the number of cut pomerons
(strings) (the “volume” fluctuations) do not lead to the correlation between
the intensive variables, e.g. the $p_{tB}$-$p_{tF}$ correlation
($b_{p_{t}\text{-}p_{t}}$). This makes
${\text{p}_{\text{t}}\text{-}\text{p}_{\text{t}}}$ correlations robust against
the volume fluctuations and the details of the centrality selection. Rather
the long-range $p_{tB}$-$p_{tF}$ correlation indicates the fluctuations in
“quality” of sources.

##  2 String Fusion Model

Due to the fact that quark-gluon strings have the finite transverse size they
can overlap in the impact parameter plane. The interaction of color strings
can be performed in the framework of local string fusion model [12].

According to this model, mean multiplicity of charged particles and mean
$p_{t}$ originated from the area $S_{k}$, where k strings are overlapping are
the following:

$\left\langle\mu\right\rangle_{k}=\mu_{0}\sqrt{k}\frac{S_{k}}{\sigma_{0}},\hskip
28.45274pt\left\langle p_{t}\right\rangle_{k}=p_{0}\sqrt[4]{k},$

where $\sigma_{0}=\pi r_{str}^{2}$ – string transverse area.

The calculations are simplified in a cellular variant of the string fusion
model by the introduction of the lattice in the impact parameter plane [7, 6,
8, 16, 10]. In the discrete model, the transverse plane is considered as a
grid with the cell area equals to string transverse area, and strings are
fused if their centers occupy the same cell. The cellular variant of string
fusion has been shown to provide the results very close to the ones in the
original string fusion model [17].

##  3 Dipole-based MC SFM model

The Monte Carlo model [16, 18, 13] used in this paper is based on dipole
picture of elementary partonic interactions. The model preserves energy and
angular momentum conservation in the initial state of a nucleon [18, 19]. The
probability of dipoles interaction depends on their transverse coordinates
$\vec{r}_{1},\vec{r}_{2},\vec{r}_{1}^{\prime},\vec{r}_{2}^{\prime}$ with
effective coupling:

$f=\dfrac{\alpha_{s}^{2}}{2}\ln^{2}\dfrac{|\vec{r}_{1}-\vec{r}_{1}^{\prime}||\vec{r}_{2}-\vec{r}_{2}^{\prime}|}{|\vec{r}_{1}-\vec{r}_{2}^{\prime}||\vec{r}_{2}-\vec{r}_{1}^{\prime}|},$
(3)

where $(\vec{r}_{1},\vec{r}_{2}),(\vec{r}_{1}^{\prime},\vec{r}_{2}^{\prime})$
are transverse coordinates of the projectile and target dipoles.
With confinement effects taking into account, [20, 21], the probability
amplitude is:

$f=\frac{\alpha_{S}^{2}}{2}\Big{[}K_{0}\left(\frac{|\vec{r}_{1}-\vec{r}_{1}^{\prime}|}{r_{\text{max}}}\right)+K_{0}\left(\frac{|\vec{r}_{2}-\vec{r}_{2}^{\prime}|}{r_{\text{max}}}\right)-K_{0}\left(\frac{|\vec{r}_{1}-\vec{r}_{2}^{\prime}|}{r_{\text{max}}}\right)-K_{0}\left(\frac{|\vec{r}_{2}-\vec{r}_{1}^{\prime}|}{r_{\text{max}}}\right)\Big{]}^{2}.$
(4)

Here $r_{max}$ is characteristic confinement scale.

Multiplicity and transverse momentum in the model are calculated in the
approach of color strings, stretched between projectile and target partons,
taking into account their finite rapidity width and local cellular string
fusion. Strings are assumed to emit particles independently according to
Poisson distribution (no short-range effects are included). Parameters of the
model are constrained from the p-p, p-Pb and Pb-Pb data on total inelastic
cross section and multiplicity [22].

##  4 Results in dipole-based MC SFM model

In dipole-based Monte Carlo string fusion model, we calculated the mean
transverse momentum correlation coefficient ($b_{pt-pt}$) for Au-Au collisions
RHIC energy and p-Pb and Pb-Pb collisions at the LHC energies. The results are
shown in figure 1.

Figure 1:  Centrality dependence of $b_{p_{t}-p_{t}}$ for Au-Au collisinos at
$\sqrt{s_{NN}}$=200 GeV (top left), p-Pb collisinos at $\sqrt{s_{NN}}$=5.02
TeV (top right) and Pb-Pb collisions at $\sqrt{s_{NN}}$=2.76 TeV (bottom). MC
simulations at $r_{\mathrm{str}}$=0.2 fm.

The results show a general growth of $b_{p_{t}\text{-}p_{t}}$ with centrality
reflecting the increase of string density and, hence, the role of string
fusion from peripheral to central collisions. We note that in Pb-Pb collisions
at the LHC energy the correlation coefficient has a maximum in mid-central
collisions and further decreases with centrality. This regime reflects the
attenuation of color field fluctuations due to the string fusion at large
string density.

We note that the non-monotonic behavior of $b_{pt-pt}$ with centrality is
achieved only in heavy ion collisions at the LHC, while in Au-Au collisions at
RHIC and p-Pb at LHC the maximal string density is not enough to provide a
decline of the correlation coefficient for most central collisions.

##  5 Mean ${\text{p}_{\text{t}}\text{-}\text{p}_{\text{t}}}$ correlations in
different models

We studied the mean ${\text{p}_{\text{t}}\text{-}\text{p}_{\text{t}}}$
correlations LHC energy in different Monte Carlo generators. The following
Monte Carlo model were considered:

  * •

MC model with string fusion [18, 13] (see above).

  * •

THERMINATOR 2 (THERMal heavy IoN generATOR 2) [23]. Based on parametrized
freeze-out hypersurface, Cooper-Frye particlization and decays.

  * •

HIJING (Heavy Ion Jet INteraction Generator) [24]. It includes gluon shadowing
effects and jet quenching.

  * •

DPMJET, two-component Dual Parton Model, based on the Gribov-Glauber approach
[25]. It models soft and hard components and performs fragmentation of partons
by the Lund model.

  * •

AMPT (A Multi-Phase Transport Model for Relativistic Heavy Ion Collisions)
[26]. It includes gluon shadowing, Zhang’s Parton Cascade, string melting,
relativistic transport.

The results for $\text{p}_{\text{t}}\text{-}\text{p}_{\text{t}}$ correlations
in Pb-Pb collisions at LHC energy are shown in figure 2.

\begin{overpic}[width=429.28616pt]{all_s_5}
\put(4.0,50.0){\large$\displaystyle b_{p_{t}\text{-}p_{t}}$} \end{overpic}
Figure 2:  Centrality dependence of $b_{p_{t}\text{-}p_{t}}$ for Pb-Pb
collisions at $\sqrt{s_{NN}}$=2.76 TeV in different Monte Carlo models.

The results demonstrate that the model with parameterized initial states shows
zero correlation coefficient. In the models, which include only initial state
shadowing and include soft+hard components, the forward-backward
${\text{p}_{\text{t}}\text{-}\text{p}_{\text{t}}}$ correlation is small and
independent on centrality.

The relativistic transport model and string fusion model demonstrate
significant non-trivial centrality dependence of $b_{p_{t}\text{-}p_{t}}$ as
well as its non-monotonic behavior.

The comparison clearly shows that
${\text{p}_{\text{t}}\text{-}\text{p}_{\text{t}}}$ forward-backward
correlation and its centrality dependence is sensitive to the initial stages
of heavy ion collisions.

##  6 Summary and conclusions

The dependence of the correlation strength between mean-event transverse
momenta on the collision centrality is obtained for different collision
energies. It is shown that above RHIC energy the dependence reveals the
decline of the correlation coefficient for most central collisions, reflecting
the attenuation of color field fluctuations due to the string fusion at large
string density.

The long-range correlation between intensive observables, being robust against
the volume fluctuations and the details of the centrality determination,
enables to obtain the signatures of string fusion at the initial stage of
hadronic interaction in relativistic heavy ion collisions at LHC energy.

The sensitivity of the long-range
${\text{p}_{\text{t}}\text{-}\text{p}_{\text{t}}}$ correlations to the
properties of the initial stages of heavy ion collisions can be further
studied in fully event-by-event hydrodynamical models, such as iEBE-VISHNU
[27] or EKRT [28].

The authors acknowledge Saint-Petersburg State University for the research
grant 11.38.242.2015.

## References

## References

  * [1] Amelin N _et al._ 1994 Phys. Rev. Lett. 73 2813–2816
  * [2] Braun M and Pajares C 1992 Phys. Lett. B 287 154–158
  * [3] Braun M and Pajares C 1993 Nucl. Phys. B 390 542–558
  * [4] Braun M A and Pajares C 2000 Phys. Rev. Lett. 85 4864–4867 (arXiv: hep-ph/0007201)
  * [5] Alessandro B et al. (ALICE Collaboration) 2006 J. Phys. G 32 1295--2040
  * [6] Vechernin V and Kolevatov R 2007 Phys. Atom. Nucl. 70 1797--1808
  * [7] Vechernin V and Kolevatov R 2007 Phys. Atom. Nucl. 70 1809--1818
  * [8] Braun M, Kolevatov R, Pajares C and Vechernin V 2004 Eur. Phys. J. C 32 535--546 (arXiv: hep-ph/0307056)
  * [9] V V Vechernin and R S Kolevatov 2004 Vestn. SPbU. Ser. 4. Fiz. Khim (4) 11--27 (arXiv: hep-ph/0305136)
  * [10] V V Vechernin and R S Kolevatov 2004 Vestn. SPbU. Ser. 4. Fiz. Khim (2) 12--23 (arXiv: hep-ph/0304295)
  * [11] Capella A and Krzywicki A 1978 Phys. Rev. D 18 3357
  * [12] Braun M and Pajares C 2000 Eur. Phys. J. C 16 349--359 (arXiv: hep-ph/9907332)
  * [13] Kovalenko V and Vechernin V 2014 EPJ Web Conf. 66 04015 (arXiv: 1308.6618 [nucl-th])
  * [14] Sputowska I and Rybicki A 2016 Correlations in Particle Production in Nuclear Collisions at LHC Energies Ph.D. thesis Cracow, INP presented 10 Oct 2016 URL http://cds.cern.ch/record/2223886
  * [15] Drozhzhova T A _et al._ 2016 Phys. Atom. Nucl. 79 737--748 (original Russian text: 2016 Yad. Fiz. 79(5) 508)
  * [16] Kovalenko V and Vechernin V 2012 PoS Baldin ISHEPP XXI 077 (arXiv: 1212.2590 [nucl-th])
  * [17] Vechernin V and Lakomov I 2012 PoS Baldin ISHEPP XXI 072
  * [18] Kovalenko V 2013 Phys. Atom. Nucl. 76 1189--1195 (arXiv: 1211.6209 [hep-ph])
  * [19] Kovalenko V 2013 arXiv: 1308.1932 [hep-ph]
  * [20] Gustafson G 2009 Acta Phys. Polon. B 40 1981--1996 (arXiv: 0905.2492 [hep-ph])
  * [21] Flensburg C, Gustafson G and Lonnblad L 2011 JHEP 1108 103 (arXiv: 1103.4321 [hep-ph])
  * [22] Kovalenko V 2013 PoS QFTHEP2013 052
  * [23] Chojnacki M _et al._ 2012 Comput. Phys. Commun. 183 746--773 (arXiv: 1102.0273 [nucl-th])
  * [24] Wang X N and Gyulassy M 1991 Phys. Rev. D 44 3501--3516
  * [25] Roesler S, Engel R and Ranft J Advanced Monte Carlo for radiation physics, particle transport simulation and applications. Proceedings of the Monte Carlo 2000 Conference, Lisbon, 23--26 October 2000
  * [26] Lin Z W, Ko C M, Li B A, Zhang B and Pal S 2005 Phys. Rev. C 72 064901 (arXiv: nucl-th/0411110)
  * [27] Shen C _et al._ 2016 Comput. Phys. Commun. 199 61--85 (arXiv: 1409.8164 [nucl-th])
  * [28] Niemi H, Eskola K J and Paatelainen R 2016 Phys. Rev. C 93 024907 (arXiv: 1505.02677 [hep-ph])