Time–Cost–Quality Trade ‐ Off in a Broiler Production Project Using Meta ‐ Heuristic Algorithms: A Case Study

: The global production of broiler meat was forecasted to be 97.8 MT in 2019. The cost fluctuations create risks in production. In order to have an effective management system, process uncertainty must be taken into account. This approach considers the process as an interval with fuzzy numbers and, for managing the risks, it uses the variable α , a parameter determined by the manager in an interval between 0 and 1. Then two algorithms, namely the multi ‐ objective imperialist competitive algorithm (MOICA) and multi ‐ objective particle swarm optimization (MOPSO), were compared and applied. Since the process of production has many activities and each activity has possible options, the process does not have a unique solution. Therefore, the objective function and its assigned weights in terms of time, cost, and quality can be applied to select the best solution from those obtained. A vast amount of uncertainty can be observed, and effective management necessitates dealing with these uncertainty issues. The MOPSO algorithm showed better performance than the MOICA algorithm in this problem. Based on fuzzy logic and influenced by the uncertainty condition ( α = 0), time, cost, and quality in the MOPSO and the MOICA algorithms were 1793.8 h, $260,571.7, and 46.66%, and 1792.5 h, $260,585.7, and 51.19%, respectively.


Introduction
The poultry sector, the fastest expanding livestock sector, is the leader in global meat production growth in response to rising demand for this healthy, high protein, and low-fat type of meat [1]. Poultry is the second most consumed meat after pork in the world (14.99 and 16.02 kg/capita/year, respectively) and is expected to take over the leading position in the future. High costs of production, high food exchange coefficient in Iran in comparison with developed countries, high mortality of chickens during the growing period, and fluctuations in the price of the animals' feed make, as for other processes, an effective control and management system necessary. Furthermore, to accelerate a process, some of the process activities can be completed faster than normal by either spending less time on them or more money to cover the expenses of additional resources allotted to them. Thus, accelerating a process increases its cost.
Generally, several methods of time and cost optimization are categorized as heuristic, mathematical, and meta-heuristic. Added to dimensions and complexity matters, the possibility of solving common optimum or fast calculating methods in a suitable time has been reduced, and as a result, receiving the whole set of optimum solutions in such a situation is very difficult. The In this study, the life span of thirty thousand chickens from eggs of broiler breeder hens to the slaughterhouse are summarized in fourteen activities. Each activity has some options for resource utilization and the goal is finding the optimal/near optimal ways of project completion in the search space of whole possible combinations of these resource utilization options to activities. Considering the situation of the implementation of activities, on average, each activity has 4 different executive approaches. Each approach itself contains time, cost, and quality in the shape of the individual triangular fuzzy number ( Figure 2).

Fuzzy Theory
In this section, some basic definitions of fuzzy theory that was defined by Kaufmann, Gupta and Qi, and Passino and Yurkovich [28][29][30], are presented. The characteristic function γA of a crisp set A ⊆ X assigns a value of either 0 or 1 to each member in X. A fuzzy function can be generalized to a function μA that is called the membership function, and the set A = {(A, μA(x): x ∈ X} defined by μA (x) for each x ∈ X is called a fuzzy set. A fuzzy set A defined on the set of real numbers R is said to be a fuzzy number. A fuzzy number is a generalization of a regular, real number. It refers to a connected set of possible values, where each possible value has its own weight between 0 and 1. A fuzzy number is thus a special case of a convex, normalized fuzzy set of the real line.
A triangular fuzzy number is one of the fuzzy numbers that is represented with three points as follows: A = (I, M, U): this representation is interpreted as a membership function. These parameters represent the smallest possible value (I), the most promising value (M), and the largest possible value (U).
Different α cuts are defined by the manager between 0 and 1, where when α cut decreases to zero, uncertainty and risk increase, whereas when it increases to one, uncertainty and risk decrease. When defining α cut, it is crucial to consider also the fluctuations and the risk prediction for the production cost. For each different α, different solutions for cost, time, and quality are produced and, consequently, the objective function is calculated.
In this study, the number of completion conditions of the broiler production process was the act of rearranging the members of a set into a sequence or order; it was represented by m n , where 'm' represents the number of different approaches for each activity and 'n' represents the number of activities in the production process (Table 1), and each activity could be done by four different approaches, as explained in Section 2.1. It was therefore the aim to find the best condition (smallest amount of cost and time and highest quality) among 4 14 possibilities. Table 2 illustrates the sample of approaches for activity 7 that explains that activity 7 has four approaches. As mentioned, each approach has three fuzzy numbers for each of the total time, total cost, and total quality.

Name of Activity Number of Activity
Buy and transport eggs from broiler breeders to incubation house and put them into the trays and fumigate them 1 Put trays into the setter and hatcher for hatching 2 Count broiler chicks and put them into the chick boxes 3 Wash and fumigate tray setter and hatcher 4 Wash and fumigate poultry house 5 Buy and transport broiler chicks from incubation house to poultry house 6 Grow up broiler chicks in the poultry house 7 When broilers reach market weight, they are caught and put into transport crates 8 Buy and transport poultry from poultry house to slaughterhouses 9 Slaughter, pick feathers, and eviscerate 10 Wash poultry and put them into spin chiller 11 Hang out poultry and transport them to a cold room 12 Pack and weigh 13 Buy poultry and transport them to the emporium 14

Depicting Activities Network Diagram
After defining the activities of the production process, their interrelationships were identified to construct a network diagram ( Figure 3). A network diagram is a chart in which the activities appear ordered and sequenced on a chronological basis, according to their interconnections. To design the diagram, the "activity on node (AON)" [31] method was used, as shown in Figure 3. According to this method, activities are shown on nodes and vectors indicate their chronological succession. It is obvious that "wash and fumigate poultry house (activity 5)" is independent of incubation house (activities 1-3); therefore, activity 5 is connected from the start to activity 6. Furthermore "Wash and fumigate tray setter and hatcher (activity 4)" must be done after activity 3 and it is independent of activities 6-13.

MOPSO Algorithm: A Proposal to Problem Solving
The algorithms used in this study were performed by MATLAB 2014a software (MathWorks, Natick, MA, USA). In the MOPSO algorithm, each solution for the problem identified is equal to one particle. The proposed algorithm was as follows: for each particle, considering the number of activities, 14 numbers were defined randomly between 0 and 1 ( Figure 4). As the MOPSO works by a continuous range of variables, the solutions must be made randomly in a continuous range, and then they must be changed to an integer range for calculation of the objective function using Equation (1). Therefore, it is important to consider the existing constraints of the defined number of approaches for each activity. In Equation (1), the selected number between 0 and 1 was converted into an integer number in order that a random method could be selected for each activity.
where i is activity number; int is integer number for each activity; R is randomized number between (0,1); and M is number of methods for each activity. After determining the approach of performing each activity, time, cost, and quality for each approach of performing activities were gathered in the fuzzy form. Then each particle has a fuzzy time, fuzzy cost, and fuzzy quality. The defuzzification (non-fuzzy) is made by using the gravity center method for each particle's time, cost, and quality [32]. Defuzzification is the process of obtaining a single number from the output of the aggregated fuzzy set. It is used to transfer fuzzy inference and results in a crisp output. In other words, defuzzification is realized by a decision-making algorithm that selects the best crisp value based on a fuzzy set. There are several forms of defuzzification including center of gravity (COG), mean of maximum (MOM), and center average methods, the most popular of which is the (COG) method. The COG method returns the value of the center of the area under the curve (Equation (2)).
x is a point in the horizontal axis (i = 1, 2, 3, …), and ( ) is the membership value of the resulting conclusion set according to ( Figure 2). In the next stage, considering the 3-dimensionality of objective function space, the solutions must be determined that in all three dimensions are not overcome against the other solutions in Equation (3).
where x1, x2, and x3 are dimensions of time, cost, and quality for x activity, respectively, and y1, y2, and y3 are dimensions of time, cost, and quality for y activity, respectively. So far, the objective function was defined, and it was mentioned how to calculate the objective function for each solution (particle) and compare them with each other to determine non-dominated solutions (particles). The non-dominated particles are reserved in a repository. In this paper, velocity update for a particular particle is done by selecting one of the non-dominated particles as a leader. Thereby, an initially grid must be performed because particles have no preference for each other to being selected as the leader.
In this way, every dimension of the problem, the space between the biggest and smallest member of the repository, is determined and inflated in a certain extent by an inflation rate (in order that the biggest and smallest member of the repository enter into the grid as well). The determined grid number in the algorithm parameters was considered according to the following: this inflated space was divided equally, and the high and low limit for each part of the grid was determined, as seen in Figure 5.
In the next stage, it was determined which member of the repository belongs to which part of the grid and then each part's particle number was determined. More particles in one part's grid means a higher density of particles, and its selection probability is low when they are selected by roulette wheel (known as fitness proportionate selection) to associate a probability of selection.  This probability is defined by Equation (4); then one part of the grid with a roulette wheel, randomly one particle of that part, is selected as a leader and then for each particle, a leader is selected by the previous method. In the following, a velocity is defined for each particle. Velocity is zero at the beginning of the algorithm and updated in every iteration, and the particle is placed in a new position (Equation (5)).

 
where P is selection probability; nG is number of particles in each grid part; and β is selection pressure.
where Vel(i) is velocity of particle(i); W is inertia weight; C1 is personal learning coefficient; BestP(i) is the best particle of each iteration; Pop(i) is every particle of the population; LeaderP is the best particle until the current iteration; and C2 is the global learning coefficient. To increase exploration in the algorithm, a mutation was created in new position of each particle; considering the mutation rate, some activities (variables) were selected randomly in each particle and the approach's number of those activities was changed randomly with respect to its constraint. In the case that the mutation created the more optimal particle, it was accepted. Otherwise, the algorithm was continued by the non-mutated particle. It is worth mentioning that the mutation rate is reduced by an increase of the algorithm iteration according to Equation (6).
where PM is percentage of mutation; i is number of iterations; and Mu is the mutation rate. The new position of each particle was evaluated according to the objective function. If the objective function's value of the new position of the particle is better than the best particle (global memory) of all particles, it is substituted for the best particle (global memory) of all particles. Then, the non-dominated particles are identified and added to existing particles in the repository. The non-dominated particles in the repository are separated and the grid stages are performed again. If the number of particles within the repository is more than the repository capacity, the additional particles must be removed from the repository. Obviously, particles with the least preference for becoming a leader must be removed from the repository. To do this, initially using Equation (7), for each part of grid, a probability is considered such that every part of the grid that has more particles, using a roulette wheel, has a higher selection chance and then one of the particles of that part of grid is eliminated at random. On the contrary, when every part of the grid has fewer particles, it has a lower selection chance by the roulette wheel because the particles should support all the searching space and the chance of selecting particles for removal should be increased in the crowded part of grid. This stage was repeated until the number of existing particles in the repository was equal to the repository capacity.
where P is probability of selection in each grid; nG is the number of particles in each grid part; and  is selection pressure. Each member of the repository (best particles) was a possible position in the implementation of the production process of a broiler. In this case, the best particle was selected from existing particles in the repository based on our considered coefficients of time, cost, and quality (regarding the importance of factors like time, cost, and quality for each manager, the coefficients can be different).
In this study, the relative weight of these criteria such as time, cost, and quality were considered as 0.34, 0.34, and 0.3, respectively (Equation (8)). The objective function coefficients of Equation (8) are determined according to the manager's opinion and considering the existing work conditions. These coefficients determine which factors (cost, time, and quality), taking work conditions into account, must have higher importance.
where OF is objective function; WT is weight of time; T is total time of process; WC is weight of cost; C is total cost of process; WQ is weight of quality; and Q is total quality of process. This procedure was continued until there was no considerable changes in the function utility, according to the factor coefficients. Other parameters were needed for this algorithm, which are shown in Table 3. Parameters were obtained through trial and error for obtaining the best solution along with logical running time.

MOICA Algorithm: A Proposal for Problem Solving
This algorithm was tested using MATLAB 2014a software (MathWorks, Natick, MA, USA). As Figure 6 shows, MOICA begins with creating a country. Each country has an initial population of possible solutions for this problem and contains 14 variables, where each variable has different approaches when considering constraints. Every country is either an imperialist or a colony. Imperialists are assigned to the most powerful countries. Every imperialist makes empires by taking control of colonies. In the following, the imperialists attract all their colonies and create revolution in the imperialist and colonies (like the mutation in the genetic algorithm).
The basis of this algorithm is the competition among empires. In this process, all colonies move from weak empires to the more powerful ones. At the end, MOICA converges to one empire that remains. In this state, the best solution of the optimization problem is calculated when the objective function's value is the same for all colonies and the corresponding imperialists. Like MOPSO, at first the MOICA algorithm selected 14 numbers randomly between 0 and 1 for each country (Equation (1)). The numbers were converted into integers in order that a random approach's number was selected for each activity. After determining the approaches of performing for each activity, time, cost, and quality for each approach of the production process activities were gathered in the fuzzy form. Then each country had a fuzzy time, fuzzy cost, and fuzzy quality. In this stage defuzzification (non-fuzzy), was made by using the gravity center method for each country's time, cost, and quality. In the next stage, to sort countries, considering the 3-dimensionality of the objective function space, the solutions must be ranked, and the crowd distance of each country must be computed. The country set that is in the same rank constitutes the Pareto front surface (used in the multi-objective genetic algorithm). It is worthwhile mentioning that the more the crowding distance of a country is, the emptier the area of the objective function shows, as determined using Equation (9).
where d is crowding distance; i is member of population; D is number of objectives; and f is objective function value. In that case, the countries were sorted based on the ranking (in ascending order) and then they were sorted based on the crowding distance (in descending order).
In the next step, the best countries were separated and selected as imperialists, and the remaining countries were assigned to imperialists by roulette wheel. Each imperialist that had the fewest components in sorting had an increased chance for selection by roulette wheel and the taking of more countries (Equation (4)). After determining and assigning countries to each empire, the countries in each empire moved towards an imperial country (belongs to each empire) with the assimilation coefficient determined using Equation (10) (Figure 7).

( ) ( ) C R (Imp(i) Col(i))
where Col is Colony(i) (country(i)); C is assimilation coefficient; R is random number; and Imp is Imperialist (i). In the next stage, considering the revolution percent, in the countries of each empire a revolution was made by random and in view of the revolution rate, some of the variables were selected randomly and the approach of those variables were changed. Like country as colony, the revolution could also be created in the imperialist countries; if the revolution is not destructive, a more optimal country is produced than the previous imperial country, and the new country is substituted for the previous imperialist country. In the following, new countries were formed by assimilation and revolution, then they competed with their imperialist country and if the imperialist countries were defeated, they were replaced. Then empires had their turned to compete. For competition between empires, Equation (11) must be used for an entire empire, and a value is defined to empires that compare them with each other.
where Val is value of the jth empire; Imp is the imperialist of the jth empire; Z is the colony's (country's) mean objective function's value coefficient; Col is ith colony of the jth empire; and nCol is the number of colonies for the jth empire. Like MOPSO, to sort empires, crowd distance was used. After sorting empires, the first and last empire were identified as the most powerful and weakest empire, respectively. In the following, weakest countries were sorted according to the previous method and the weakest country of the weakest empire was selected, as seen in Figure 8. In the next stage, the weakest country was eliminated from the weakest empire and added to one of the remaining empires selected by the roulette wheel (every empire that has the fewest components in sorting has an increased chance for selection by roulette wheel and taking the weakest country, Equation (4)).
This trend caused some of the empires to be eliminated while performing the algorithm and finally, the imperial countries remained in the empires, which showed us the problem solutions. Parameters used in this algorithm are presented in Table 4. Parameters were obtained through trial and error for obtaining the best solution along with logical running time. A counter in the objective function was applied to measure the reference numbers to objective function in both algorithms in Table 5.

Number of Functional
Evaluation Algorithm 4050 Multi-objective particle swarm optimization 3450 Multi-objective imperialist competitive algorithm

Results and Discussion
For reaching absolute solution in this study, around 270 million possible combinations that contained 14 activities with an average of four subcontracting alternatives have been checked. Project performance such as project time, cost, and quality changed with each possible combination. Searching in the large solution space of a possible solution was improved by the present multiobjective optimization model. The proposed methods can reduce the large solution space by finding a non-dominated solution over successive generations.
In this research, the two algorithms MOPSO and MIOCA were used to trade off between cost, time, and quality of the broiler production process. To compare the two algorithms in solving the problem, the initial population and iteration number in both algorithms were assumed equal; these two algorithms were used in this study by using Equation (1). Tables 6 and 7 show the best solutions for each amount of fuzzy cut based on the lowest objective function. As the tables show, the more (α) increased, the lower the objective function amount was because the uncertainty and risk decreased. The best solutions are given in Tables 6 and 7, which show subcontracting plans and project objectives for each of the solutions. Furthermore, these tables show notable solutions with the optimum values for the three objectives in the Pareto-optimal front and the solution of the best compromise. The results in Tables 6 and 7 indicate that the (α = 1) exhibits the minimum objective function, and the (α = 0) exhibits the maximum objective function. These results make it clear that when (α) is one, the time and cost are the lowest amount and quality is the highest amount due to the fact that time and cost are increased, and quality is decreased by decreasing the (α). This assertion is concluded by center of gravity defuzzification. When fuzzy cut is moved to the base of the fuzzy triangle number (α→0) by the manager, the output number of defuzzification is increased. The objective function is increased by increasing the time and cost and decreasing the quality based on Equation (8). According to the results, the amount of objective function in MOPSO is lower than that in MOICA for all fuzzy cut (α); therefore, MOPSO had better performance than MOICA for this problem.
If the contractor is allowed to increase cost by maximizing the project quality, then he/she should select the solution that has maximum quality value and that minimizes time, but that does not globally minimize the cost. In this case, he/she should give a higher weight to quality and a lower weight to time and cost in the objective function. If a contractor would like to minimize the project cost, then he/she needs to select a solution to achieve the minimum cost value by incurring additional duration with low project quality. In this case, he/she should give a higher weight to cost and a lower weight to time and quality in the objective function. Furthermore, if a contractor has a priority for minimizing project time, he/she should give a higher weight to time and a lower weight to cost and quality in the objective function.  According to Figures 9 and 10, in algorithm MOPSO, lesser amounts were obtained in equal iterations for the objective function, and the algorithm MOPSO had a better search ability to find the optimal solution for the concerned problem than the algorithm MOICA. According to Table 5, the number of functional evaluation (NFE) in MOPSO was lower than that in MOICA, although the numbers of the first population and its iteration in both algorithms were equal. Additionally, considering the convergence of algorithms and the achievement of a constant value, MOPSO showed better performance in contrast to MOICA. Figures 9 and 10 show that under different (α), MOPSO achieved a constant value in fewer iteration numbers. This showed higher efficiency of the MOPSO than MOICA.  As was mentioned before, among all the solutions obtained from the algorithms, one of the solutions that, in view of the utility function according to Equation (8), had the least amount was selected for performing the activities, and so each activity was carried out via a method that led to an optimum general answer.
Mungle et al. [33] identified 25 solutions for the same benchmark problem using fuzzy clustering-based optimization using the MOGA, and the best compromise solution was also identified. In order to demonstrate that the currently proposed methodology can generate more efficient schedules for construction projects, the reported results from Mungle et al. [33] were ranked based on the evidential reasoning (ER) approach. The weights of the attributes were assumed to be equal in order to make the comparison valid. The best solutions identified by them [33] were the three solutions with the highest utility scores using the ER approach, along with their respective attribute values.
Rahimi et al. [34] used multi objective particle swarm optimization for a discrete time, cost, and quality trade-off problem and compared it to the genetic algorithm. They reported that the performance of PSO rather than GA in the same condition was better, and that PSO was more efficient. Their result was the same with this research in terms of better performance of PSO.
Zhang and Xing [35] focused on solving the fuzzy time, cost, and quality trade-off (TCQT) problem through a fuzzy-multi-objective particle swarm optimization (FMOPSO) based on a fuzzy multi-attribute utility methodology, in which constrained fuzzy arithmetic operations are adopted to reflect the inequality constraints among the involved fuzzy numbers. The efficiency or the convergence performance of the FMOPSO was compared with a GA method that also adopts the fuzzy multi-attribute utility and constrained fuzzy arithmetic proposed in their research. Their resulting comparison showed that the proposed FMOPSO can solve the TCQT problem with almost the same efficiency as the GA method based on the same fuzzy multi-attribute utility methodology and the constrained fuzzy arithmetic.

Conclusions
This study showed that MOPSO's performance was better than that of MOICA in the integer domain. Under the uncertainty conditions (α = 0), the time, cost, and quality were computed as 1793.8 h, $260,571.7, and 46.66%, respectively, by MOPSO. The MOICA results were 1792.5 h, $260,585.7, and 51.19%, respectively. The risk in the production process was managed by adjusting the alpha value according to the working conditions. Based on the method that is determined for each activity by algorithms, the least time and cost and highest quality, as much as possible, are obtained.