A two-agent scheduling problem in a two-machine flowshop

Article history: Received March 18 2017 Received in Revised Format August 2


Introduction
The scheduling problems in which several agents compete for using the shared resources have recently attracted the attention of many researchers.In these problems, contrary to the general multi-objective problems in which the optimality criteria reflect data from the consideration of all jobs, each agent has its own set of jobs and follows its own optimality criterion only.Baker and Smith (2003) and Allesandro Agnetis et al. (2004) introduced the two-agent scheduling problems for the first time.They introduced the positive combination of each agent's objective, constrained optimization, and Pareto optimization case of the problem.In the constrained optimization case, the goal is to optimize the objective function of one agent under the condition that the objective functions of the other agents are smaller than the predetermined values.Also, in the Pareto optimization problem, the goal is to find a set of non-dominant solutions and their related sequences.
In this paper, a two-agent constrained optimization scheduling problem in the two-machine flowshop is studied.The goal of this problem is to minimize the total tardiness of the first agent's jobs subject to the condition that the makespan of the second agent's jobs is allowed to be lower than an upper bound.In the following, a number of related studies dealing with the constrained optimization case of the multiagent scheduling problems are reviewed.Kellerer and Strusevich (2010) provided a Fully Polynomial Time Approximation Scheme (FPTAS) for the scheduling problem in a single-machine environment with the objective of minimizing the total weighted completion times of the first agent while the makespan of the second agent was bounded.Wu et al. (2011) offered a branch and bound algorithm and some heuristic procedures for the problem in a single-machine environment with learning effects to minimize the total tardiness of the first agent and a limited number of tardy jobs of the second agent.Wu (2014) investigated a two-agent single-machine scheduling problem with learning effects via an objective function which minimizes the weighted completion time of all the jobs subject to a constraint that one agent's makespan cannot exceed a prescribed upper bound.They proposed a branch and bound algorithm along with three simulated annealing procedures for the problem.Gajpal et al. (2014) investigated a single machine scheduling problem with the aim of minimizing the total weighted completion times of the first agent's jobs and bounded completion times of the second agent, and presented three heuristics for this problem.Yin et al. (2012a) addressed a two-agent scheduling problem in a single-machine environment with arbitrary release dates, where the objective was to minimize the total tardiness of one agent, while keeping the total lateness of the other agent below or at a fixed level.They developed a mixed integer programming (MIP), a branch and bound procedure, and a marriage in honey-bees optimization algorithm (MBO) for this problem.Yin et al. (2012b) considered a single machine scheduling problem with deteriorating jobs where the objective of each agent was the maximum amount of a regular function of completion times of jobs.For this problem, they provided an optimal procedure with polynomial time complexity.Wu et al. (2013) deliberated upon a two-agent single-machine scheduling problem with deteriorating jobs.The objective of this problem was to minimize the total weighted number of tardy jobs of the first agent's jobs subject to the condition that the maximum lateness of the second agent is allowed to be lower than an upper bound.They proposed a branch and bound algorithm and a tabu search algorithm for this problem.Liu et al. (2013) proposed a new scheduling model with two-agent and sum-of-processingtimes-based deterioration.The objective is to minimize the total completion time of the first agent with the restriction that the makespan of the second agent cannot exceed a given upper bound.They proposed optimal properties of the problems and the optimal polynomial time algorithms to solve it.Cheng et al. (2011) considered a two-agent single-machine scheduling problem involving deteriorating jobs and learning effects, simultaneously.The objective is to minimize the total weighted completion time of the jobs of the first agent with the restriction that no tardy job is allowed for the second agent.They developed a branch-and-bound and several simulated annealing algorithms to solve the problem.For more studies on this line of research, the reader may refer to Alessandro Agnetis et al. (2007) and Leung et al. (2010).
As seen above, all the reviewed studies are concerned with the single machine environment.However, it is of note that usually more than one operation must be executed on every job in many manufacturing systems.Often, these operations have to be performed on all the jobs in the same order, implying that the jobs have to follow the same route.This environment is referred to as a flowshop (Pinedo, 2002).According to the literature, few studies have been conducted on multi-agent scheduling in a flowshop environment.Alessandro Agnetis et al. (2004) showed that the two-agent constrained optimization scheduling problem in a two machine flowshop intended to minimize the makespan of the first agent's jobs with the limited makespan of the second agent is NP-hard.For this problem, Luo et al. (2012) provided a dynamic programing algorithm with pseudo-polynomial time complexity and a FPTAS.Lee et al. (2010) considered a two-agent scheduling problem in the two-machine flowshop, where the objective was to minimize the total tardiness of the first agent's jobs with the restriction that the number of tardy jobs of the second agent is zero.They proposed a branch and bound algorithm and a simulated annealing algorithm and claimed that the branch and bound algorithm can solve problems involving up to 15 jobs in size in a reasonable amount of time.Lee et al. (2011) considered a two-agent scheduling problem in the two-machine flowshop, where the objective was to minimize the total completion time of the first agent's jobs with the restriction that the number of tardy jobs of the second agent is zero.They presented a branch and bound and a simulated annealing procedure to solve the problem, and showed that the branch and bound algorithm is capable of solving instances containing up to 20 jobs in size in a reasonable amount of time.Fan and Cheng (2016) studied two-agent scheduling in a two-machine flowshop.The cost function was the weighted sum of some common regular functions, including the makespan and the total completion time.For the first problem, they proposed an ordinary NP-hardness proof and a pseudo-polynomial-time algorithm.For the second problem, they proposed an approximation algorithm based on linear programming relaxation of the problem.Shiau et al. (2015) studied a two-agent scheduling problem in a two machine flowshop with learning effects.The objective is to minimize the total completion time of the jobs from one agent, given that the maximum tardiness of the jobs from the other agent cannot exceed a bound.They also provided a branch-and-bound algorithm for the problem.In addition, They presented several genetic algorithms to obtain near-optimal solutions.Lei (2015) studied flow shop scheduling problem with two agents and considered its feasibility model.In this problem the goal was to minimize the makespan of the first agent and the total tardiness of the second agent simultaneously under the given upper bounds.They proposed a simple variable neighborhood search (VNS) algorithm for this problem.Recently, Perez-Gonzalez and Framinan (2014) and Alessandro Agnetis et al. (2014) have reviewed the multi agent scheduling problems, which can be referred to for more insight into such problems.
As single agents (one objective) in the two machine flowshop are related to the problem considered in this paper, some relevant studies are reviewed in the following.At this point, it is noteworthy to say that the essential conclusion ever reached for minimizing makespan in a two-machine flowshop in the scheduling theory is known as Johnson's rule (Johnson, 1954).Johnson (1954) provided the optimal schedule using a simple procedure with polynomial time complexity for this problem.Although finding the optimal solution for the problem of minimizing the total tardiness in this environment is limited to search in the permutation schedule, the problem of minimizing total tardiness in the case of considering zero value of the due date values reduces to the problem of minimizing the total completion time, which is a problem with strongly NP-hard computational complexity (Graham, 1979).Thus, in most of such studies, implicit enumeration approaches have been used.Sen et al. (1989), Kim (1993), Pan and Fan (1997), Pan et al. (2002) and Schaller (2005) provided a branch and bound algorithm for minimizing the total tardiness in the two-machine flowshop.For this problem, Kharbeche and Haouari (2012) presented simple MIP models based on the position of jobs in the sequence.Their results showed that the efficiency of mathematical models is higher than that of the branch and bound algorithm.Hoogeveen and Velde (1995), Della Croce et al. (1996), Federico Della Croce et al. (2002) and Akkan and Karabatı (2004) developed a branch and bound algorithm for the two-machine flowshop scheduling problem with the objective of minimizing the total completion times.Haouari and Kharbeche (2013) presented an assignment-based lower bound for the two-machine flowshop scheduling problems with a regular additive performance criterion.Their computational results showed that the quality of this lower bound is better than the lower bounds in the literature focusing on the special cases of minimizing both the total completion time and the total tardiness in the two-machine flowshop.
The paper is organized as follows: the problem definition and used symbols are given in the next section.The theorems and properties of the optimal solutions are presented in section 3.In section 4, the mathematical model is introduced.Furthermore, tabu search method in section 5, computational results in section 0, and in the last section, conclusions and suggestions for future research are presented.

Problem definition and notations
In this paper, the two-agent flowshop scheduling problem is dealt with.As delay penalties and efficient resource utilization are more practical criteria in the scheduling problems, in the considered problem the objective is to find an optimal schedule to minimize the total tardiness of the first agent, under the situation that the makespan of second agent is bounded.All jobs are available at the zero time and have to follow the same route i.e. they must be processed first on machine 1 and then on machine 2. Let JA and JB (JA∩JB=∅) denote the job sets belonging to agent A and agent B, respectively, where JA consists of nA jobs and JB has nB jobs and therefore n (n=nA + nB) jobs must be scheduled.For each job j belonging to agent x, a processing time on the first machine, a processing time on the second machine, and a due date are defined.Considering a schedule S, the completion time of job j of the agent x on the second machine is denoted by and the tardiness of this job is defined as 0, . Also, the makespan of agent x is defined as .Following the notations offered by Graham et al. (1979) as well as the classification of Allesandro Agnetis et al. (2004), this problem is shown as 2|| ∑ : . When the number of jobs belonging to agent B is zero or Q value is large enough, the problem 2|| ∑ : is reduced to minimizing the total tardiness in a two-machine flowshop scheduling problem, which is strongly NPhard (Lenstra et al., 1977).Therefore, 2|| ∑ : is also strongly NP-hard.
Because of the limited makespan of the second agent's jobs and influence of the Q value on determining the solution space, different aspects of the problem with respect to different values of the Q are investigated as follows.
As the makespan of the second agent is limited to Q, if Q is less than the minimum time required to complete the second agent's jobs, the problem will be infeasible.Also, if it is large enough, the problem reduces to the problem of minimizing the total tardiness in the two-machine flowshop scheduling.This is why the Q value has to be especially addressed to study the problem.Suppose that * and # indicate the minimum and maximum time required to complete all the jobs belonging to agent x, respectively.# is obtained by applying Johnson's rule (Johnson, 1954) and # is calculated by the inverse of the sequence obtained by Johnson's rule (Johnson, 1954) By considering these values, there exist the following three cases: Second case: If * = Q, the second agent's jobs are arranged at the beginning of the sequence based on Johnson's rule (Johnson, 1954) and the first agent's jobs are arranged after them at the end of the sequence.In this case, the problem is to find the optimal sequence of the first agent's jobs arranged after the second agent's jobs according to the total tardiness as the objective function, Third case: Consider the situation in which the Q value is so large that if the first agent's jobs are arranged at the beginning of the sequence with the minimum total tardiness (the optimal sequence), the second agent's jobs can be placed after them at the end of the sequence while maintaining the feasibility condition.First, the maximum possible required time to arrange the first agent's jobs at the beginning of the sequence is calculated.Then, with regard to this calculated value as a virtual job (indicator of the first agent's jobs) the second agent's jobs are assigned after the first agent's jobs.In other words, in this case we need to suppose that L is the lower bound of Q (L ≤ Q).To determine the value of L, at first, # is calculated which is equal to the upper bound of the completion time of the first agent's jobs.Then, the virtual job that has zero processing time on the first machine and processing time # ∑ on the second machine is considered.Establish ∪ i.e. the set of the second agent's jobs and virtual job and calculate ∪ * value as the minimum time required to complete the second agent's jobs arranged after the first agent's jobs.Therefore, with the definition of ∑ ∪ *

, If
, the second agent's jobs are arranged at the end of the sequence according to the Johnson's rule (Johnson, 1954).The first agent's jobs are placed at the beginning of the sequence, and the two-agent scheduling problem reduces to finding the optimal sequence of the first agent's jobs at the beginning of the sequence according to the total tardiness as the objective function.According to the above cases, the problem is investigated further for the values * ∑ ∪ * , where the Q value shows its effect on the optimal sequence apart from the three mentioned cases.

Theorems and optimal properties of the problem
In this section, some efficient theorems and optimal properties of the problem used in the developed methods of the next sections are introduced and some of them are given in the appendix.The presented propositions are classified based on the jobs belonging to each agent.In these theorems and propositions, is a schedule in which job j is processed before job i and ′ is a schedule that is identical to S, except for the fact that jobs j and i are interchanged.Also, when both job i and job j belong to one agent, for brevity the indexes of agent are removed.
For the cases in which both job i and job j belong to the first agent with the objective of minimizing the total tardiness of jobs, the following theorems are discussed.
Theorem 1 (Pan & Fan, 1997): For any two jobs i and j belonging to the first agent, if conditions , and are met, then there exists an optimal schedule such that job j is processed after job i.
Theorem 2 (Pan et al., 2002): For job i, if there exists a job j satisfying conditions , and , then there exists an optimal schedule such that job i is not the first job of the sequence.
When both job i and job j belong to the second agent with the bounded makespan, the following propositions are developed in this paper.
Proposition 1: For any two adjacent jobs i and j belonging to the second agent, the Johnson's rule (Johnson, 1954) is established between them.
Proof: Based on the Johnson's rule (Johnson, 1954), if job i is assigned before job j, it must be shown that by interchanging jobs i and j in S, . Let Ii and Ij denote idle times in S occurring on the second machine immediately prior to the processing of jobs i and j, respectively.Also, let Ii′ and Ij′ denote idle times in S′ occurring on the second machine immediately prior to the processing of jobs i and j, respectively.Further, let F and C denote completion times of jobs that are completed before these jobs on the first machine and the second one, respectively.According to Fig. 1 .With regard to these cases, by interchanging job i and job j in S, the makespan of these two jobs does not increase.■

Fig. 1. Interchanging two adjacent jobs
Corollary 1: Due to the transitivity property in the Johnson's rule (Johnson 1954), if the adjacent jobs belonging to the second agent are more than two jobs, then there exists an optimal schedule such that the Johnson's rule (Johnson, 1954) is established among all of them.
Proposition 2: For two jobs i and j belonging to the second agent, if conditions ai≤aj and bi = bj are met, then there exists an optimal schedule such that job j is processed after job i.
Proof: According to Fig. 2, let M1 and M1' be the total processing times of jobs between job i and job j on the first machine in S and S', respectively.Also, let M2 and M2' be the total processing times plus idle times on the second machine in S and S', respectively.The completion times of job i and job j in schedule S and S' are calculated by relations and , and therefore we have .Because and , by interchanging these jobs in S, the completion times of jobs between job i and job j will not increase.So, if 0, we have and if 0, we have .Therefore, , and the completion times of subsequent jobs will not increase.□ Fig. 2. Interchanging two non-adjacent jobs Proposition 3: For two jobs i and j belonging to the second agent, if conditions , and are met, then there exists an optimal schedule such that job i is not the first job of the sequence.
Proof: According to Fig. 3, if the first and third conditions are met, by interchanging job i and job j in S, the completion time of jobs between them will not increase.Thus, because of no increase in the completion time of jobs between job i and job j, similar to the proof of Proposition 2, the relation is established.Therefore, according to the second condition, the completion times of subsequent jobs will not increase.■

Mathematical programming model
Due to the complexity of the studied problem, a branch and bound algorithm is proposed to achieve the optimum solution.In this algorithm, the theorems, optimal properties (proposed in section 3 and appendix), and some lower bounds are used.However, computational results show that the branch and bound algorithm is only able to solve all the problem instances up to 16 jobs in size.Due to the low efficiency of the branch and bound algorithm compared with the mathematical model of the problem, the details of this algorithm are not presented in this paper, and instead the mathematical model is provided as the efficient exact method to optimally solve the problem.Kharbeche and Haouari (2012) presented three mathematical models for the two-machine flowshop problem to minimize the total tardiness.These models are based on the jobs' positions in the sequence and formulated by defining waiting time decision variables, completion time decision variables, and machine's idle time decision variables.Also, Chandra et al. (2009) proposed precedence-based formulations for the permutation flowshop with a common due date.However, according to the computational results of Kharbeche and Haouari (2012) and M'Hallah (2014) as well as preliminary tests which were conducted in this study, the mathematical programming models based on the position of jobs in the sequence showed higher efficiency to achieve optimal solutions.Also, among these models, the mathematical model based on the completion times of jobs on the machines showed the best performance.According to the mentioned issues, the mathematical programming model of 2|| ∑ : is presented below, which is denoted as ACT (Agent's jobs Completion Times) model.This model is based on the jobs' positions in the sequence and obtained by defining the completion times on machines as decision variables.This model is formulated as relations (1) to ( 12).In this model, decision variable T[k] is equal to the tardiness of job that is scheduled in the position k, and is a binary variable which takes value 1 if job j is assigned to position k in the sequence and zero otherwise.Decision variable and decision variable are completion times of the scheduled job at position k on the first machine and on the second machine, respectively.

∑
(1) In this model, the objective function (1) minimizes the total tardiness for jobs which are assigned in different positions of the sequence.Since tardiness is calculated in this model only for the first agent's jobs, this equation is equal to minimizing the total tardiness of the first agent's jobs.Eq. ( 2) and Eq. ( 3) are assignment constraints and respectively state that at each position there is sequenced only one job and each job must be assigned to exactly one position in the sequence.Constraints (4) to (8) determine the completion times of jobs at different positions of the sequence on the first machine and the second one.Constraint (9) states that tardiness of a position is calculated if its sequenced job belongs to the first agent.Similarly, Constraint (10) states that the completion time of a position is limited to the Q value if the first agent's job is placed in this position.The formulation of ACT includes O(n 2 ) binary variables, O(n) continuous variables, and O(n) constraints.

Optimal properties in MIP
In section 3, Theorem 1 (Pan & Fan, 1997) and Proposition 2 was proposed for the problem 2|| ∑ : . Set Γ including the pairs of jobs i and j (in which job i is processed before job j) is formed by using the mentioned properties.Since ∑ states the position index of job i, the following inequality is established for members of the set Γ: Also, by using Theorem 2 (Pan et al., 2002) and Proposition 3, the values of variables associated with the jobs which cannot be placed at the beginning and at the end of the sequence, are set to zero.Also, since in the members of Γ, job i is processed before job j, job j and job i cannot be placed at the beginning and at the end of the sequence, respectively.
Finally, by calculating the lower bound of completion times of jobs in different positions of the sequence, the procedure which was provided by Haouari and Kharbeche (2013), some variables corresponding to the second agent's jobs can be set to zero because of their calculated completions times' lower bound which is greater than Q.To handle this case, the following procedure is presented: Step 1: Sort all jobs (including the first and second agent's jobs) according to the Johnson's rule (Johnson, 1954).
Step 2: Calculate the lower bound of completion times of jobs in different positions according to the procedure proposed by Haouari and Kharbeche (2013).
Step 3: Because the lower bound of the completion times of jobs in different positions has been identified, for each one of the second agent's jobs, specify the first position in which the completion time's lower bound of the job is larger than Q.
Step 4: For each one of the second agent's jobs, set to zero the value of corresponding variables in the first position determined in Step 3 up to the last position of the sequence.
Among the offered optimal properties, inequality (13) was introduced by Kharbeche and Haouari (2012) to improve the efficiency of mathematical programming model of minimizing the total tardiness in a twomachine flowshop.They applied Theorem 1 (Pan & Fan, 1997) to form set Γ.As mentioned before, in addition to utilizing this inequality in this study to form set Γ, Proposition 2 has been used.Also, Kharbeche and Haouari (2012) used Theorem 2 (Pan et al., 2002) to determine those jobs that are not placed at the beginning of the sequence.The other used properties in the mathematical programming model were developed in this study.

Tabu search algorithm
Tabu search is one of the well-known meta-heuristic algorithms which guide local search procedure by exploring the solution space of a problem beyond local optimality.This algorithm and its principles were introduced by Glover (1989Glover ( , 1990) ) for the first time.In this section, a Multi Start Tabu Search (MSTS) algorithm is proposed to solve 2|| ∑ : . The main components of this algorithm include different methods to create the neighborhoods, the neighborhood searching type including the number of neighbors and used properties to determine the neighbors, memory structures including tabu list length and aspiration criterion, and finally the condition of terminate search.MSTS creates several initial solutions from solution space at first.For each of these solutions named 'starts', a separate tabu list is created.Tabu search continues to search until reaching the stopping criterion.In each iteration, the algorithm chooses a 'start' probabilistically based on its cost to carry out searching.The best sequence corresponding to this start is replaced with this and the corresponding tabu list is updated.When tabu search is stopped, the best solution among the latest found solutions stored for each 'start' is chosen for a local search.

Initial solutions
The success of MSTS strictly depends on the initial solutions, and the aim of heuristic algorithms is to disperse these solutions widely on the solution space.To generate initial solutions, at first the second agent's jobs are sorted according to the Johnson's rule (1954).For the first agent's jobs, three orders are considered including the earliest due date (EDD) rule, the shortest processing time (SPT) rule of sum of job's processing time on both machines, and the random order of jobs (Józefowska et al., 1994).According to these orders, algorithms H1 and H2 are proposed as follows: Algorithm H1: To generate the initial solution in this method, the sequence is divided into two parts.In the first part, the largest possible number of the first agent's jobs is placed with respect to the feasibility condition of the second agent's jobs.Then, the second agent's jobs are placed after them.The remainder of the first agent's jobs is placed at the end of sequence in the second part.This procedure is performed according to the following steps: Step 1: Sort the second agent's jobs based on the Johnson's rule (Johnson, 1954) and sort the first agent's jobs according to each of the EDD, SPT and Rnd rules.
Step 2: Select the maximum number of the not sequenced first agent's jobs, such that the not sequenced second agent's jobs can be placed at their last authorized position before Q.Place these selected jobs in the sequence and remove them from the not sequenced first agent's jobs.
Step 3: If the entire first agent's jobs are sequenced, place the entire second agent's jobs after the last sequenced job and go to step 5; otherwise, place the first not sequenced second agent's job after the sequenced jobs in the sequence and remove them from the not sequenced second agent's jobs.
Step 4: If all the second agent's jobs are sequenced, place the not sequenced first agent's job after the last sequenced jobs in the sequence and go to step 5; otherwise go to step 2.

Algorithm H2:
In this algorithm, the sequence is divided into two parts.First, the second agent's jobs are placed in the first part of the sequence and then the first agent's jobs are placed after them.
In the MSTS algorithm, by using H1 and H2 algorithms and considering the three EDD, SPT and Rnd rules as the orders of the first agent's jobs, 6 initial solutions are produced.Also, by applying a local search on any of these 6 solutions, 6 other initial solutions are generated.In the applied local search, the generated initial sequences are improved through interchanging each pair of jobs and using swapping and insertion procedures to be described in the next section.First, all the possible interchanges are checked and then a move with the greatest amount of improvement is performed.This process is continued until no further interchange leads to improvement.Thus, according to the 6 obtained initial solutions, 12 initial solutions are produced for MSTS algorithms.As noted, the MSTS algorithm selects a 'start' with the corresponding tabu list based on its cost to search in each iteration.Before starting the searching in the algorithm, the initial solutions (starts) are sorted in non-ascending order of their costs.In this order of 'starts', the probability of selecting the ith 'start' from among all L 'starts' for searching in current iteration is (2*i)/L(L+1).After selecting a 'start' and completion of the search process, this "start" and corresponding tabu list are updated, and with regard to its cost it is placed in the proper position between all other "starts".Searching algorithm is similarly continued until the stop condition is seen.

Neighbor generation
The search space is a space of feasible solutions which can be visited during the searching procedure.Also, the neighborhood structure definition depends on the search space.At each iteration of the algorithm, a local moving is performed on the current solution, and neighbor solutions in the search space are generated.In the MSTS, the three methods including swapping, insertion and reversion are used to create neighbors.In these methods, two positions of the sequence are randomly chosen.Suppose that job i and job j are corresponding jobs to these positions in the sequence of (σ i π j ω).In this sequence, σ is a partial sequence of jobs before job i, π is a partial sequence of jobs between job i and job j and ω is a partial sequence of jobs after job j.In the swapping, the sequence is changed to (σ j π i ω).There are two cases in the insertion; in the first case, job i is transmitted to after job j and the sequence is changed to (σ π j i ω).In the other case, job j is transmitted to before job i and the sequence is changed to (σ j i π ω).In the reversion, jobs i and j and their between jobs are reversed and the sequence is changed to (σ j π' i ω) where π' is a reverse of π.

Neighbor search
In each iteration of the MSTS, a number of neighbors are generated for searching.Some optimal properties of the problem are used to generate these neighbors.As stated before, to create new sequences, two positions are randomly chosen.Then, one of the three methods of swapping, insertion and reversion are selected randomly.After generating a new neighbor, Corollary 1 (local property of Johnson's rule (Johnson, 1954)) is applied to the second agent's jobs.Also according to the proposed properties, neighbors are not produced in the three cases below: First: If the selected jobs for moving and jobs which are placed between them belong to the second agent, according to Corollary 1, none of the three methods are applied and the neighbor is not created.Second: If the selected jobs for moving belong to the second agent, and the conditions of Proposition 2 are established for them, swapping method is not applied.Third: If any of the interchange methods violate the feasibility condition of the second agent, the moving will not be applied and the neighbor not created.

Memory structure and search terminate condition
From among the tabu search elements which distinguish tabu search from local search, its memory and tabu structures are more important.Tabu structures are used to prevent the cycle and escape of local optimums.As noted, there can be set some conditions for the repeal of tabu structures which are known as aspiration criteria.In MSTS, the aspiration criteria allow to accept a tabu movement whose obtained solution is the best solution encountered by the algorithm up to the current iteration.Finally, after a stop condition is met, the best found solution is improved by a local search.Three stop conditions are used in MSTS, which can stop the algorithm by establishing each of them.These conditions include the number of iterations without improving the best found solution, reaching a certain maximum iteration, and reaching a zero objective value as the optimal solution during the search process.The pseudo-code of the MSTS algorithm is shown in Fig. 3.

Computational results
In this section, the results of solving some problem instances are provided to evaluate the efficiency of the mathematical programming model and tabu search procedures.The proposed procedures were coded in Visual studio 2013 programming environment and in C# programming language.Also, CPLEX 12.6 was used to solve the mathematical programming model.The time limit to optimally solve the instances was set to 1800 seconds.The proposed procedures have been implemented on a computer with Intel® Core ™ i7-2600M CPU @ 3.4 GHz and 4 GB RAM in Windows 7 with 32-bit operating system.

Instances
To evaluate the effectiveness of the provided procedures for solving the problem, the number of instances was considered based on the study of Lee et al. (2010), Lee et al. (2011) andYin et al. (2012ab).For all instances, the processing times of the first agent's jobs and second agent's jobs were randomly generated from a uniform distribution on the interval [1,10].The due dates of the first agent's jobs follow a discrete uniform distribution on the interval [δ (1-τ-R/2), δ (1-τ+R/2)], where R and τ are measures of dispersion and tardiness factor, respectively.Also, the value δ is the lower bound of completion time of all jobs and it is equal to the summation of the processing times on the second machine and the minimum processing time of all jobs on the first machine.For parameters R and τ, the values {0.25, 0.50} and {0.25, 0.50, 0.75} are selected, respectively.Also, three values of {0.25, 0.50, 0.75} for the proportions of the first agent jobs to the total jobs (ρ) are considered.According to the study of Yin et al. (2012ab), if L is the lower bound of Q and U is its upper bound, the value of L + q (U-L) is considered in the generated instances as the Q value, where q is one of the three values {0.25, 0.50, 0. 75}.Given the negligible impact of various R values on the performance of the presented methods, the results are not provided for various R values.Therefore, according to the values of τ, ρ and q, the instances were generated in 18(2×3×3) different groups.These groups are named as G01 to G18.For each one of the three intended values of R, 10 instances and therefore 30 instances were generated in each group.

Evaluation of mathematical models
In this section, the results of solving instances by the mathematical programming model are presented.As noted, some optimum properties were used in this model in order to increase its efficiency.Computational results are presented in Table 1.In this table, the column 'ACT' shows the results of the basic mathematical model without using optimum properties, and the column 'MACT' (Modified ACT) indicates the results of the mathematical model with the consideration of optimal properties.In this table, the column 'No of Opt' shows the number of instances which were solved optimally.In this column, if for one of the three values for R at least one instance could not reach the optimal solution in 1800 seconds as solution time, solving other instances with this value of R and results of the 10 instances corresponding to this value of R are not reported.Thus, according to the three values of R, values of this field can be 30, 20 and 10.The column 'Avg.CPU time' represents the mean solution times of instances.Also, the column 'Improve' shows the percentage value of improvement in solution time by using optimum properties.These values were calculated by the following equation:

Evaluation of the performance of Tabu search algorithm
In order to evaluate the performance of the tabu search algorithm, the error of its solutions according to the optimum solutions was calculated.In the following, the tuned parameters of the algorithm along with the results are presented.

The MSTS tuned parameters
In this section, parameters of the MSTS are tuned by using the Taguchi method.The tuned parameters include the number of examined neighbours in each iteration, the tabu list length and the number of iterations without any improvement as a stop condition.As mentioned, there are 4 parameters, and each one has 3 levels that showed in Table 2. Thus, the proper orthogonal array was selected by Minitab Software and the L9 orthogonal array (OA) was used.The average mean and S/N ratio at each level for relative deviation (RD) values are shown in Fig. 4 and Fig. 5, respectively.As indicated in these figures, better robustness of the algorithm is achieved when the number of examined neighbours in each iteration set as n and the maximum number of iteration set as 10000.For parameter "the number of iterations without any improvement", whether 40n or 60n can be selected and for increasing effectiveness of MSTS, this parameter set as 60n.Finally, for parameter "the tabu list length", two level n/2 and n can be selected, too.To increase efficiency and saving time, the length of the tabu list was set to n/2.In Table 3, the selected values for the parameters of the algorithm are presented.Fig. 6 shows the AEP of the MSTS algorithm for the proportion of the first agent's jobs (ρ).The results indicate that this value decreases by increasing the number of the second agent's jobs.This is because of applying Corollary 1 (local property of Johnson's rule (Johnson, 1954)) in searching procedure.When the number of second agent's jobs increases, according to Corollary 1, the optimal order of the second agent's adjacent jobs is Johnson's order that reduces the search space.However, for the first agent's adjacent jobs a predefined sequence does not exist.So, when the first agent's job is increased, the solution space and the number of its points, which the search algorithm fails to investigate, are increased and this leads to the increase in the value of AEP.

Conclusion
In this paper, the two-agent scheduling problem in a two-machine flowshop with the aim of minimizing the total tardiness of first agent's jobs under the bounded makespan of the second agent's jobs was studied.A mathematical programming model with the theorems and useful properties were used to optimally solve the problem.Computational results showed that the presented mathematical programming model is able to solve all the instances in 84.56% of groups with up to 60 jobs in size.Also, the computational results indicated that the proposed tabu search algorithm had high efficiency according to the average absolute error value which is less than 0.18%.Future research line related to this problem could involve feasibility optimization, Pareto optimization and optimization of the weighted sum of the agent's objective.Also, because few studies have been done on the multi-agent scheduling in the flowshop environment, the important problems not solved in this field of scheduling are recommended as further areas of inquiry.

Fig. 3 .
Fig. 3. Interchanging any job of the sequence with the first one

Fig. 4 .
Fig. 4. Main effects plot for means Fig. 5. Main effects plot for SN ratios

Fig. 6 .
Fig. 6.The AEP of the MSTS for ρ values , the completion times of jobs i and j in S and S′ are calculated by relations and , respectively.So, it is enough to show that .
The best neighbor of S which has the minimum cost among all Ss neighbors encountered by the algorithm.Code 1. Begin 2. Set nbr ← Null Sequence; nbr-cost← ∞; 3. Sort the agent A jobs and agent B jobs in EDD and JR, respectively; 4. Generate L initial sequences for the L starts; 5. Sort the L sequences generated in Step 4 in non-decreasing order of their costs; 6. for iter from 1 to k do 7. Choose a start (S1) probabilistically to carry out a tabu search iteration, for this start; 8. Obtain the best tabu and non-tabu neighbors of S1 by using neighborhood search procedure; 9. Select the best neighbor Snbr keeping a check on the aspiration criterion, and update the tabu list using tabu tenure t; 10.Replace S1 by Snbr as the latest sequence for that start.If stop criteria encountered, go to 12; 11.End for; 12. Perform a neighborhood local search on the best sequence among all starts; 13.Obtain the best sequence S * encountered during the search; 14.Set nbr ← S * ; nbr-cost ← cost of S * ; 15.Output nbr and nbr-cost; 16.End : An instance of size n, parameters L, t, and k.Output:

Table 1
The MIP models resultsResults indicate that the MACT is more efficient than ACT in most cases.The MACT is able to solve all the instances containing up to 40 jobs in size, in 92.56% of groups and instances having up to 60 jobs in size, in 84.92% of groups.These results for ACT are equivalent to 86.11% and 77.77%, respectively.Also, mean percentage of improvement in solution time is equal to 30.40%.According to the results, when the number of first agent's jobs increases, more time is required to solve instances.

Table 2
Factors levels

Table 3
The MSTS parameters tuning resultsThe absolute error percentage (AEP) of the MSTS algorithm is calculated according to Eq. (15).The values of average absolute error percentage (AEP) of the MSTS algorithm are presented in Table4.These values are obtained by comparing the solution of MSTS algorithm with the optimum solutions of MACT model.In Table4, 'No of Inst.' column shows the number of instances which were solved by MSTS algorithm.The results show that AEP and the average solution time of the MSTS algorithm for solved problems are 0.08% and 26.16 seconds, respectively.

Table 4
The AEP for the MSTS algorithm