Variable Neighborhood Strategy Adaptive Search to Solve Parallel-Machine Scheduling to Minimize Energy Consumption While Considering Job Priority and Control Makespan

: Environmental concerns and rising energy prices put great pressure on the manufacturing industry to reduce pollution and save energy. Electricity is one of the main machinery energy sources in a plant; thus, reducing energy consumption both saves energy costs and protects our planet. This paper proposes the novel method called variable neighborhood strategy adaptive search (VaNSAS) in order to minimize energy consumption while also considering job priority and makespan control for parallel-machine scheduling problems. The newly presented neighborhood strategies of (1) solution destroy and repair (SDR), (2) track-transition method (TTM), and (3) multiplier factor (MF) were proposed and tested against the original differential evaluation (DE), current practice procedure (CU), SDR, TTM, and MF for three groups of test instances, namely small, medium, and large. Experimental results revealed that VaNSAS outperformed DE, CU, SDR, TTM, and MF, as it could ﬁnd the optimal solution and the mathematical model in the small test instance, while the DE could only ﬁnd 25%, and the others could not. In the remaining test instances, VaNSAS performed 16.35–19.55% better than the best solution obtained from Lingo, followed by DE, CU, SDR, TTM, and MF, which performed 7.89–14.59% better. Unfortunately, the CU failed to improve the solution and had worse performance than that of Lingo, including all proposed methods.


Introduction
The manufacturing industry is facing a great deal of pressure with regard to saving energy and reducing emissions, since it is an energy-intensive industry. In 2020, its energy consumption reached 1443.1 trillion British thermal units (Btu), of which electricity consumption was 200.7 trillion Btu [1]. In the same year, energy-related carbon dioxide emissions by the manufacturing industry reached 1057 million metric tons [2]. Studying energy-efficiency scheduling under electricity cost is of great significance for manufacturing firms to improve their energy efficiency. In a factory, machinery is the main energy-consuming unit [3][4][5][6]. Reducing machine energy consumption economically and environmentally improves sustainable manufacturing. There are many potential energyreduction approaches in a manufacturing plant, such as developing more energy-efficient machines and processes. In the production line, production planning and scheduling generally impact energy efficiency. Energy efficiency can also be increased by the suitable utilization of machines in the shop floor [6][7][8][9][10][11][12]. In the production process, the electricity consumption of each job is sometimes different, indicating that an operation with high electricity consumption is arranged for low processing time jobs, and an operation with low electricity consumption is arranged for high processing time jobs, so that the electricity cost can be decreased [5,6,13]. Paying too much attention to energy cost as an outstanding optimization objective may cause a long makespan and imbalanced workload on machines, as well as lead to bottlenecks, late deliveries, and even untimely machine failure [14][15][16][17].
Production scheduling is explicitly important in a modern manufactory, consisting of planning and sequencing jobs into machines, since mass production heavily relies on a large number of machines working in parallel. Parallel machines can process several jobs simultaneously without affecting each other. Parallel machine scheduling is defined as sequencing and assigning jobs into machines when similar types of machines are available and jobs can be scheduled in these machines. A variety of sequencing and/or processing restrictions often exist when decision makers try to minimize some related objective functions [14][15][16][17][18][19][20][21][22][23][24]; however, most enterprises still use advanced machines running alongside outdated ones. In contrast to the old ones, modern machines are usually adjusted to work at a high speed and save energy. Speeding up old machines to operate as fast as modern ones results in them consuming more electric power and releasing more pollutants [25]. Sequencing and assigning jobs to a particular machine that can process different jobs at different production rates is considered to be the parallel-machine-scheduling problem. Solutions or algorithms constructed to schedule parallel machines are important because, when implemented in all sorts of parallel-machine problems, they can obtain good solutions for the overall production-planning process in a factory [18,19].
Although most scheduling problems consider production efficiency, cost, and quality as optimization problems, considering the costs of energy and gas emissions-known as the green scheduling problem-has attracted the attention of researchers [9][10][11]14,25]. In addition, various constraints have been addressed in parallel-machine-scheduling problems, such as setup time, machine-available time, ready time, release date, due date, and delivery time [14,26]. In this study, the due-date constraint was considered with makespan control, since late-delivery cost and production-overhead cost per working hour may be a thorny issue for entrepreneurs. Therefore, job priority and makespan control should be seriously considered for production planning as well. On the basis of the above discussion, many scholars and manufacturing firms should be aware of the importance of energyefficiency and job-priority concerns, on top of minimizing the makespan. There are a few studies on the parallel-machine-scheduling problem with a concept of energy consumption that consider job priority and makespan control. Therefore, this addressed scheduling is both an objective optimization problem and an NP-hard problem. The differential evolutionary (DE) algorithm is suitable to solve this kind of problem because it can obtain non-dominated solutions in a single run and has been successfully applied to optimization problems [17,20,24,[27][28][29]. In addition, metaheuristics based on local search methods, such as the variable neighborhood strategies adaptive search (VaNSAS), have been successfully applied to solve many combinatorial optimizations problems [30][31][32][33][34][35][36][37][38][39][40], which inspired us to develop a parallel-machine-scheduling model, and to propose VaNSAS and new neighborhood strategies: (1) solution destroy and repair (SDR); (2) track-transition method (TTM); and (3) multiplier factor (MF).
Nowadays, the cost of energy consumption is key in terms of production efficiency and environmental sustainability for industrial firms. A significant amount of research has been done on parallel-machine-scheduling problems that are strongly NP-hard. We refer to Chen [41], Pinedo [42] and Behera [43] for the discussion on the complexity of parallelmachine scheduling. Therefore, this study investigates the parallel-machine-scheduling problem for minimizing energy consumption while considering job priority and makespan control. In order to obtain the near Pareto front, we developed VaNSAS: (1) solution destroy and repair (SDR); (2) the track-transition method (TTM); and (3) multiplier factor (MF). On the basis of the problem characteristic and evolutionary algorithm, the proposed VaNSAS uses a new encode scheme that can convert discrete optimization problems into continuous optimization problems, and applies VaNSAS with three different techniques to further improve quality.
The remainder of the study is organized as follows: Section 2 reviews the literature related to parallel scheduling, considering energy consumption, job priority, and makespan control, including DE and VaNSAS applications; Section 3 describes a formal definition of the addressed problem and a mathematical model for the proposed problem; Section 4 describes the proposed VaNSAS, SDR, TTM, and MF with multiple operators in order to handle the scheduling problem; Section 5 presents and analyzes the experiment results of our proposed VaNSAS, SDR, TTM, and MF; Section 6 provides the conclusion and suggestions for future studies.

Literature Review
Environmental and sustainable energy is a critical topic. CO 2 emissions into the atmosphere are strongly related to energy consumption in human activities; they are generated by traditional fuels from natural resources which are becoming depleted [44]. Liu [39] studied scheduling problems that were related to an environment where the objective of minimizing CO 2 emissions was investigated with total weighted tardiness (TWT). The non-dominated sorting-based genetic algorithm II (NSGA-II) and the ε-archived genetic algorithm (ε-AGA) were presented to solve two batch-scheduling problems. In their subsequent article, Liu and Huang [10] demonstrated both the effectiveness of an adaptive multi-objective genetic algorithm (AMGA) in finding the Pareto optimal set, and the efficiency of NSGA-II in bicriterion scheduling on a batch-processing machine with dynamic job arrivals. Regarding low carbon emissions in the parallel-machine-scheduling problem, Pan et al. [40] introduced the advantages of a novel imperialist competitive algorithm (ICA) to minimize total tardiness, total energy consumption, and CO 2 emissions.
Various energy concerns (e.g., energy consumption, energy efficiency, carbon footprint, and energy cost) are investigated in production-scheduling problems. There are also several studies on energy-efficient production systems. For example, Mouzon et al. [45] improved product operation methods, such as several dispatching rules and mathematical programming, in order to minimize the total completion time and energy consumption of manufacturing equipment in production-scheduling problems. Angel et al. [46] introduced a randomized approximation algorithm for the problem of energy-consumption scheduling for unrelated parallel machines with the average weighted completion time. Sobottka et al. [47] demonstrated the development and evaluation of a digital method for multiobjective optimization problems considering traditional business aims and energy efficiency in a metal-casting manufactory. The improved method, including hybrid simulation-based multi-criterion optimization as a multi-stage hybrid heuristic and metaheuristic method, was employed in a heat-treatment process, which requires order batching and sequencing on parallel machines under complex restrictions.
In the context of energy cost, a number of works attempted to reduce it during production when energy and electricity prices vary with time of use. For example, Fang et al. [7] proposed a mixed integer programming model for optimizing the operating schedule of a flow shop considering both productivity-and energy-related criteria. Zeng [23] studied the uniform parallel-machine-scheduling problem with electricity cost and time-dependent or time-of-use electricity tariffs, where electricity price changes by the working hour within a day. Firstly, a bi-objective mixed-integer linear-programming model was constructed for this problem. Then, the proposed method (an insertion algorithm) was applied to minimize total electricity cost and the number of operated machines. Zhou et al. [19] presented a multi-objective differential evolution algorithm to solve the parallel-batchprocessing machine-scheduling problem (BPM) in the presence of dynamic job arrivals and a time-of-use pricing scheme. The objective was to simultaneously minimize makespan and minimize total electricity cost (TEC). Recently, Nanthapodej et al. [48] introduced the hybrid differential devolution algorithm and adaptive large neighborhood search, and demonstrated a superior performance in finding high quality solutions within a short computation time of the proposed algorithm for solving the parallel-machine-scheduling problem, with minimizing total energy cost (TEC) as key for environmental sustainability, and controlling machine load balance as an indicator of production efficiency.
Although more attention is now paid to energy cost in parallel-machine-scheduling problems, some constraints for job priority and makespan control are also important, such as starting time, completion time, due date, and delivery time, since late-delivery cost and production-overhead cost per working hour impact manufacturers. Some studies attempted to avoid late-or express-delivery charges, such as lateness and tardiness issues. For instance, Chaudhry et al. [49] considered the minimization of total tardiness in identical parallel-machine-scheduling problems by using a genetic algorithm (GA), and compared it with branch-and-bound and particle-swarm-optimization methods. Then, Pei et al. [50] demonstrated minimizing maximal earliness and number of tardy jobs in parallel-machinescheduling problems. That article proposed the hybrid variable neighborhood search (VNS)-gravitational search algorithm (GSA), which is a combination of VNS and GSA to find a solution. Solution-search efficiency is also an attractive goal for many studies. Maecker and Shen [51] also applied the VNS algorithm with the fast evaluation technique (FET) to improve the computational efficiency of solution finding in the identical parallelmachine problem with machine-dependent delivery times in order to minimize total weighted tardiness.
In a continuous search space for function optimization, differential evolution (DE) is a type of evolutionary algorithm that is widely implemented [27]. The DE algorithm has gained much attention from scholars since it was first presented by Storn and Price [28]. There are a variety of DE algorithms in industrial applications studied for optimization problems. For instance, the DE algorithm was implemented for solving productionscheduling problems. Wang et al. [29] demonstrated hybrid differential evolution (HDE) in scheduling problems with splitting jobs. The study proposed a global search method with block mutation and block crossover. Experiment results revealed that the proposed HDE performed better than the traditional DE did. Zhou et al. [52] showed the performance of the hybrid DE algorithm for the uniform parallel-machine-scheduling problem with arbitrary job sizes, non-identical capacities, and different speeds. The objective was to minimize makespan. The HDE algorithm was proposed for solving large-scale problems. In this algorithm, individuals were represented as a discrete job sequence. The proposed algorithm and novel mutation were designed on the basis of this representative. Wu and Che [53] proposed a memetic differential evolution algorithm for an energy-efficient biobjective-unrelated parallel-machine-scheduling problem in order to minimize makespan and total energy consumption. Efficient speed adjusting and job machine swap heuristic were introduced and integrated into the algorithm as a local search method with an adaptive meta-Lamarckian learning strategy. Zhou et al. [19] presented a multi-objective DE algorithm for effectively solving the parallel-batch-processing machine-scheduling problem while considering energy cost on a large scale. In this algorithm, individuals were encoded into job permutations and discrete mutations; crossover operators were designed on the basis of the encoding structure, and a Pareto selection operator was proposed to select what the individuals for the next population are. Defining job permutation, a heuristic was employed to group jobs into batches and schedule them on BPMs.
In addition, a procedure to minimize the total energy cost of a schedule without compromising the makespan was proposed by Li et al. [54]. They investigated the parallelmachine-scheduling problem with different-colored families, sequence-dependent setup times, and machine-eligibility restriction of the dyeing process in textile manufacturing. Generally, the dyeing optimization problem is NP-hard, and the HDE algorithm was proposed to solve real-world data problems. The proposed HDE algorithm, a special encoding and decoding scheme, was constructed to deal with the machine-eligibility constraint, and chaos theory was applied to determine the parameter settings of the DE algorithm. Kusoncum et al. [17] presented the traditional DE with embedded heuristics to obtain nearoptimal solutions in realistically sized machine scheduling, including capacitated machine restrictions and sequencing-independent setup-time considerations. Results revealed that the proposed method's performance tended to find a new optimal solution during the simulation, while the local-search-based heuristics were trapped at some local optima, and the DE was insufficient for search intensification.
Variable neighborhood strategies adaptive search (VaNSAS) is a new type of metaheuristics that was first introduced by Theeraviriya et al. [36], with the concept of solving combinatorial optimization problems. The primary idea of the VaNSAS is to allow for algorithms to search in many different areas to obtain the best possible solution by using several searching methods. VaNSAS is very flexible to use in many optimization problems, such as the location and routing problem (LRP). Theeraviriya et al. [36] studied the LRP in the Thai rubber industry, and proposed VaNSAS to solve the problem with the objective function of fuel-cost minimization, including a realistic constraint that allowed for a vehicle to collect products by visiting rubber farms more than once in cases where they had more rubber than vehicle capacity. Computation results showed that VaNSAS could find solutions for all problem sizes in much less processing time than that needed by the exact method. After that, Pitakaso et al. [37] presented VaNSAS with another LRP, the green 2-echelon location-routing problem (G2ELRP), which is a variant of the capacitated location-routing problem (CLRP) and the 2-echelon location routing problem (2ELRP). The G2ELRP aims to reduce overall fuel consumption on the basis of distance and road conditions in both echelons. A new constraint considered in the G2ELRP is that a customer can be served more than once. The finding demonstrated that VaNSAS could solve this case and be applied to other industries.
Kusoncum et al. [17] introduced another application of VaNSAS in a sugar mill as a computational tool for scheduling sugarcane-vehicle-unloading systems. The objective was to minimize the makespan of parallel-machine-capacity scheduling with a cyclic sequence, and machine restriction included sequencing-independent setup time. Numerical results showed that, during the simulation, VaNSAS could find optimal solutions. In a manufacturing case study regarding the garment industry, the VaNSAS proposed by Jirasirilerd et al. [22] presented a better solution and less computation time in order to minimize cycle time for a simple assembly line, balancing the Type 2 problem while considering the number and types of machines operated in each workstation. Recently, Pitakaso et al. [38] applied VaNSAS to minimize the cycle time while considering the limited number of machine types in a particular workstation for the special case of the simple assembly line balancing Type 2 problems, where multi-skilled workers have a set of competencies that allow them to work on more than one machine in a workstation. Results showed that VaNSAS was able to reduce the cycle time and increase assembly line effectiveness.
The studies discussed above show that global search metaheuristics (e.g., VNS, GA, DE and VaNSAS) are effective in solving optimization problems such as parallel-machine scheduling. Due to the importance of the environmental and economic impact, studies on energy-cost concerns, job priority, and makespan control are obviously attractive and important in terms of theoretical and application value; however, no previous articles had investigated this kind of problem and employed the VaNSAS algorithm to find the solution. As a result, this paper's objective is to minimize total energy consumption while considering job priority and makespan control for the parallel-machine-scheduling problem. In order to solve this problem, we developed a mathematical model, and introduced VaNSAS, SDR, TTM, and MF algorithms to further improve solution-search efficiency.

Problem Description
This paper considers the problem of energy consumption concerning scheduling, while considering job priority and makespan control for parallel machines, to improve energy consumption and production efficiency for environmental sustainability [9]. Production costs such as overhead and late delivery are also significant and should be mentioned. Therefore, the objective of this study is to minimize total energy cost, including considering due-date and lateness constraints. In this context, a set of I different jobs {j 1 , j 2 , ..., j I } were scheduled on M parallel machines {m 1 , m 2 , ..., m M }. We assumed that each job j had deterministic processing time p jm on each machine. The jobs could be assigned to any machine. Due to the job characteristics, each machine required different levels of energy consumption e jm to process each job. All jobs were available at time zero. The completion time of job j is denoted by C jm , which is the time when the processing of the job was completed on machine m. Only after the machine starting the process could no idle time be inserted into the schedule with no preemption. Additionally, in each time period, only one job could be processed on machine m.
On the basis of these assumptions and the following notations, a mathematical model is proposed to formulate this problem in order to minimize energy consumption while considering job priority and makespan control. This is defined as the minimization of energy-consumption, late-delivery, and production-overhead costs.

Mathematical Formulation
On the basis of the characteristics of minimizing energy consumption with job priority and makespan control for a parallel-machine-scheduling problem, the mathematical model is formulated. The details of index, parameters, decision variables, objective function, and constraints are as follows: Objective function The objective function in Equation (1) attempts to minimize energy used to produce all jobs, the priority cost for jobs that deliver late, and the production-overhead cost per working hour. Equations (2) and (3) show that all jobs must be assigned to at most one machine, and the start of all assigned jobs must be Job 0 (dummy job). Equation (4) is the relationship constraint of X ijm and Y im . Equation (5) confirms that job i is processed before job j in machine m only when job j is assigned to machine m. Equation (6) ensures that a job can be assigned to at most one machine. Equation (7) calculates the starting time of job j in machine m, and Equation (8) is the calculation of the completion time of job j. Equations (9) and (10) define the starting and completion time of the dummy job. Equation (11) calculates the lateness of job i, while Equation (12) is used to ensure that the lateness of job i does not exceed the predefined value (L MAX ). Equation (13) decides whether job i is late, and Equation (14) ensures that the number of tardy jobs does not exceed the predefined number (T MAX ).
The mathematical model was formulated and tested on Lingo software. The optimal solution was obtained for a small problem, albeit with intolerable computation time. Solving medium or large problems desired more computation time, even for an incomplete solution. Therefore, in this study, we developed metaheuristics to solve a medium or large problem as a realistically sized problem.

Proposed Method
When a problem becomes larger and more complicated, solving it may not be possible by mathematical methods. Therefore, we introduced a variable neighborhood strategy adaptive search (VaNSAS), a new type of metaheuristics, which successfully improved solution-search efficiency in previous studies by adapting some VaNSAS mechanisms. This study aims to solve the parallel-machine-scheduling problem in order to minimize energy consumption with job priority and makespan control, and proposes VaNSAS for improving solution-search efficiency. The goal of VaNSAS is to allow algorithms to search for the best possible solution in different areas by using several searching methods. The solution-searching steps are to find more diversification and intensification at all times, depending on the black-box methods that were designed. The search method in VaNSAS can be a basic local search, well-known heuristics, modified metaheuristics, and a newly designed local search; therefore, VaNSAS is very flexible for applying new ideas to the algorithm and can be implemented to many optimization problems.
VaNSAS is composed of two main steps: (1) generating the initial solution (the set of tracks); and (2)

Generating an Initial Solution
Indirect encoding was used to solve the proposed problem while we randomly generated the set of tracks. Table 1 shows an example of five tracks to operate in VaNSAS. The tracks shown in Table 1 were composed of 13 positions. The first 10 positions were the track elements that represented the jobs assigned to Machines A, B, and C. We used indirect encoding, so the decoding method is essential to let the tracks reflect the real solution. The decoding method is explained as follows: 4.1.1. Decoding Methods

1.
Separately apply the rank order value (ROV) between job and machine vectors. The ROV of jobs is called the job sequence at position i (SJi), and machine sequence in position m is called SMm.

2.
Assign the first job in Position 1 of SJi to Machine 1 of SMm, assign the second job in Position 2 of Sji to the machine in Position 2 of SMm, and continuously assign the remaining jobs in position i + 1 to machine m + 1 until m = M, where M is the total number of machines. In the context of job or machine restrictions, job j is assigned to machine m only when job j is allowed to be produced on machine m. If job j is not allowed to be produced on machine m, job j is assigned to the next order of machine.

3.
Continue assigning jobs in position i + 1 to the machine that has minimal total processing time among all M machines until all jobs are assigned to a machine. While assigning a job to a machine, the maximal number of tardy jobs and maximal lateness must be controlled to be less than the maximal predefined number.

4.
Calculate completion time and energy used in Step 3.

Decoding-Method Example
An example value of a scheduling problem with 10 jobs and 3 machines is provided in Table 2. Remarks: N/A is a job-machine restriction in which the job could not be produced on that machine. P im is the processing time of job i on machine m, and E im is the energy consumed to produce job i on machine m. D i is the due date of job i and a j is penalty cost of job j.
Step 1: Separately apply ROV to the job and machine tracks from the smallest to the largest value; results are shown in Table 3.
Step 2: Apply Jobs 7, 3, and 10 to Machines B, C, and A, respectively.
Step 3: Assign Job 8 to B, but it is not allowed to produce Job 8 on Machine B; thus, assign Job 8 to Machine C instead of Machine B. Since it has the lowest total processing time, Jobs 9, 8, 5, 6, 3, and 7 are assigned to Machines A, B, A, C, B, and B, respectively. The results of this assignment are shown in Table 4.  In this assignment, the job-production sequence of Machine B produces jobs {7, 4, 9, 6}, and Machines C and A produce jobs {3, 8, 2} and {10, 1, 5}, respectively. Then, the total energy consumption of this plan is THB 247. Table 5 shows the sequencing and scheduling of jobs and machines.  Table 5 illustrates that the makespan equals 98 m. The minimum labor cost was assumed to be THB 100 per m, production-overhead cost per working unit was THB 9800, and the number of tardy jobs was 5, so the penalty cost was THB 2094. Therefore, the total cost for this assignment was THB 247 + THB 9800 + THB 2094 = THB 12,141.

Performing the Track-Touring Process
The track selects one out of a certain number of neighborhood strategies with the probability function shown in Equations (15) and (17), proposed by Pitakaso et al. [39].
where G bt is the probability of black box b selecting in iteration t before it is adjusted by the edge boundary. C is equal to the total number of black box c, and b is the index of black box b or c. S bt is the weight to select black box b in iteration t. N bt−1 is the number of tracks that select black box b in the previous iteration. A bt−1 is the average objective function of all tracks that select black box b in the previous iteration. I bt−1 is a binary decision variable. It is equal to 1 if the box contains the iterative best solution of the last iteration; otherwise, it is equal to 0. F is a predefined random variable that lies between 0 and 1. K is predefined factors that are located between 1 and 5. P bt is the probability to select black box b in iteration t after the edge boundary. G Max and G min are the maximal and minimal probabilities that are allowed to select a black box, respectively. The black boxes (neighborhood strategies) applied in this research were: (1) solution decompose and repair (SDR); (2) track-transition method (TTM); and (3) multiplier factor (MF).

Solution Decompose and Repair (SDR) Method
This neighborhood strategy comprises three steps: (1) destroy the current solution by using N-job-string removal algorithm; (2) select the repair methods; (3) perform the repair method; and (4) redo Steps 1-3 until it meets the termination condition. a.
Destroy Method In this section, the destroy method was employed to disassemble the initial solution so it would become an incomplete solution that made the solution move to other search areas, and a new solution was thereby obtained. This paper applied N-job-string removal as the destroy method, the value of the considered solution on the basis of a randomly generated job sequence before removing jobs from list I. The N-job-string removal algorithm is shown in Algorithm 2. Remove job in position e + N-1 from list B 6.
Insert removed job into list L b. Repair Method After the destroy procedure deconstructed the initial solution, the repair procedure was performed to reconstruct the solution by randomly using one of two repair methods: (1) best insertion and (2) random insertion.
b.1 Best Insertion Best insertion was used to repair a solution by determining the processing time to move the job from list L into an empty machine for operating that job. The best-insertion algorithm is shown in Algorithm 3. For instance, there was a list of jobs {2,10,8,5}. Producing Job 2 consumed 14, 41, and 39 energy units for Machines A, B, and C, respectively. Since Job 2 was removed from machine C, only two choices remained, namely machines A and B. Therefore, Job 2 was placed into Machine A due to it needing the least energy to produce Job 2. After that, Jobs 10, 8, and 5 were continuously executed with the same mechanism until all jobs in B had been reassigned. b. 2 Random Insertion Random insertion is a method used to repair an incomplete solution by finding a random machine to operate the job under conditions to reconstruct the solution. Algorithm 4 shows the random insertion algorithm.

Track-Transition Method (TTM)
The track-transition method (TTM) is composed of three steps: (1) randomly select the original track from the pool of the tracks that were not selected by the TTM as the black box; (2) randomly generate a new track; (3) find the average value of the track for each track, denote this number as lower boundary (LB), and denote the track to which this number belongs as TLB; (4) find the minimal number of values in the tracks of the maximal value of the track that is not a member of track TLB, denote this number as the upper boundary (UB), and let the track to which the UB belongs be TUB; (5) generate transition rates (TRs) for every element of the track; (6) transit the original track to a new track by using Equation (18), while the value in track i in position h of the new track (VT N ih ) is constructed. A track that is neither TLB nor TUB is called an unplaced track (UT).
where VT * h is the type of track, and * can be a TLB, UP, or TUB track. For example, if we have Tracks 1 and 2, and a random track as shown in Figure 1a, then the average value in the positions of Track 1, Track 2, and TR is 0.34, 0.46, and 0.51, respectively, as shown in Figure 1b. Therefore, LB = 0.34, and TLB is Track 1. The maximal numbers of VT of Track 2 and the random track are 0.83 and 0.97; therefore, UP = 0.83 and TUB = Track 2. As a result, Equation (18) was modified as shown in Equation (19), and the result of the transition is shown in Figure 1b. The result of the new track is shown in Figure 1c.
After the track was generated, the decoding method shown in Section 4.1.1 was performed to find the answer for the proposed problem.

Multiplier Factor (MF)
MF is the neighborhood strategy of which the basic idea is to pull out the current solution from the local optimum by multiplying the current value in that position of the track by the learning multiplier factor. The new value in the track is obtained when multiplied by the multiplier using Equation (20): where R ih is the random number that corresponds to position h of track i, VT MF ih is the track after applying the MF strategy, and VT N ih is the track before applying the MF. An example of the MF method is shown in Table 6. After the MF operation, VT MF ih uses the decoding method to obtain the solution for the proposed problem. The MF is iteratively applied to the track that selects this strategy. The predefined number of iterations was previously determined. In our method, 500 iterations were set as the stopping criterion. The MF algorithm is shown in Algorithm 6. Let S be the set of all feasible solutions and consider a solution Z ijt ∈ S. A neighborhood strategy associates each Z ijt ∈ S with a neighborhood N k (Z ijt ) ⊆ F of the solution Z ijt . For this paper, the three neighborhood structures are SDR, TTM and MF. The time complexity of neighborhoods (N 1 (Z ijt ), N 2 (Z ijt ) and N 3 (Z ijt ), respectively) are determined both by its respective structure and by the solution it is being applied to. The size of neighborhood

Updating Track and Other Information
The track is updated by using Equation (21): where Z ijt+1 is the value of track i, element j iteration t + 1. α and β are predefined parameters. In this research, we defined α and β as equal to 0.3. Z 2jt is the first randomly selected track, and Z 3jt is the second randomly selected track. We iteratively performed steps in Section 4.2; the number of iterations needed to simulate depends on the predefined parameter of number of iterations (IT). The proposed methods were tested with the case study and randomly generated data. Results were compared with the current practice procedure. Details of the current practice procedure are below.

Current Practice Procedure
The current practice procedure (CU) in vegetable-farm case-study data assigns vegetables (jobs) to grow in all farms (machines), and it is composed of four steps, shown below.
Step 1. Sort jobs according to energy used from least to most used, and name this list job list (JL).
Step 2. Calculate average number of jobs per machine and call this number AM.

AM =
Total number o f jobs Total number o f machines (22) Step 3. Assign jobs to machines according to JL. Job j is assigned to the machine that uses the least energy. If that machine has more jobs than AM, the next machine that uses the least energy is selected and continuously performs until all jobs are assigned.
Step 4. Calculate the objective function of the assignment from Step 3.

Differential Evolution Algorithm
The differential evolution algorithm (DE) was used to compare it with the proposed method. DE was constructed by the following steps: (1) randomly select two other tracks that are not the track that selected DE as the black box; (2) randomly generate a new track; (3) perform the mutation process using Equation (23); (4) perform the recombination process using Equation (24); (5) perform the selection process using Equation (25) and repeat steps 1 to 5 until the termination condition is met. The DE algorithm is shown in Algorithm 7. Algorithm 7. Differential evolution (DE) pseudocode.
Set NP, CR, F, NP (size of track) Generate initial solution Begin For G = 1 to G max when G = iterations and G max = Maximum iteration Randomly generates the set of initial solution (tracks) For N = 1 to NP Perform mutation process using v i,j,G = x r 4 ,j,G + F x r 2 ,j,G − x r 3 ,j,G (23) Perform the recombination process using Perform selection process using formula

Computational Framework and Result
The proposed metaheuristics were coded in C++ and simulated on a computer with Intel (R) Core i7-3520M CPU @ 2.90 GHz Ram 8.00 GB. We tested our algorithms on four test instances: small, medium, large, and case study. The simulation was executed five times, and the best solution among all was selected. Details of the test instances are shown in Table 7.  Table 7 shows that we tested 37 sample data (12 small, 12 medium, and 12 large sample data, and 1 case study). For small test instances, the proposed methods were tested and compared with the optimal solution obtained from Lingo V.11 (mathematical method). For the medium and large samples, the proposed method was compared with the lower bound that was obtained by Lingo v.11 within 72 h of computation. P jm is the processing time of job j of machine m (in minutes). E jm is the energy used to produce job j on machine m and then converted into money units (THB). B is a constant number (THB), which is the production-overhead cost per hour of production in the factory. The proposed method was compared with the result of a traditional DE algorithm and the current practice method (CU). Details of the investigated algorithms are shown in Table 8.  Table 8 shows that six algorithms were tested in the provided sample datasets. The performance of VaNSAS and other proposed methods was compared with that of a traditional DE.
The first experiment was executed with the small and medium test instances. The stopping criterion for Lingo was the time period when it found the optimal solution; then, it collected the best solution and computation time. The stopping criteria for all proposed methods were set as 10% of the lowest computation time of Lingo. In this case, Lingo's lowest computation time to find the optimal solution was 651 s; thus, the stopping criteria of all proposed methods were set to 65 s. The percentage difference of all proposed methods with the obtained result from Lingo was found using Equation (26). Five replications were executed for each proposed method, and the best solution among the five runs is shown in Table 9. The statistical test is shown in Table 10.
where F opt is the solution obtained by Lingo during the predefined computation time, and F H is the solution obtained by the proposed methods. The solutions are shown in Tables 9-14.

Max. Error
Best Sol.

Max. Error
Best Sol.

Max. Error
Best Sol.

Max. Error
Best Sol.

Max. Error
Best Sol.

Max. Error
Best Sol.

Max. Error
Best Sol.

Max. Error
Best Sol.

Max. Error
Best Sol.

Max. Error
Best Sol. The computation results in Tables 11 and 12 illustrate that VaNSAS performed better than other proposed and traditional methods did. VaNSAS performed as well as the exact method did, and also showed the highest performance of solution stability over other algorithms as shown in Table 11. While others were significantly different from the exact method as seen in the statistical test result (Table 12), DE, CU, SDR, TTM, and MF were an improvement from the exact method by 14.59%, 4.37%, 11.19%, 11.04%, and 11.24%, respectively. However, VaNSAS was at 19.55%, which means that VaNSAS obtained the solution at 19.55% lower cost than that of the solution from Lingo. Lingo consumed 2880 min computation time, and VaNSAS used only 30 min.

Max. Error
The next experiment was tested with a large random dataset. This group of test instances was composed of 12 test instances, the number of jobs was randomly selected to be from 80 to 134, and the number of machines was set from 5 to 8 (L-1 to L-12). In this set of test instances, we included the case study, which had 201 jobs and 11 machines (C-1). In this case, Lingo was executed for 72 h to find the lower-bound solution; then, the termination condition as the computation time of all proposed heuristics was 45 min. The simulation was executed five times; the best solution and the maximum error of large instances are shown in Table 13 and the statistical results are shown in Table 14.
The computation results in Tables 13 and 14 show that VaNSAS performed better than other proposed and traditional methods did. The performance of VaNSAS is similar to the exact method, while others were significantly different from the exact methods. DE, CU, SDR, TTM, and MF were 9.88%, 1.42%, 7.89%, 8.01%, and 9.33%, respectively, different from the exact method. In addition, VaNSAS was 16.35% different from the exact method, which means it could find 16.35% lower cost than how much Lingo could. The maximum error of each method indicates that VaNSAS outperformed other methods in terms of solution stability. In addition, Lingo required 2880 min of computation time while VaNSAS only required 45 min.

Conclusions and Outlook
Green machine-scheduling problems, such as optimization problem related to energy concerns, are paid more attention by a wide array of industries, since environmental awareness is part of industrial manufacturing sustainability. This paper presents a novel method called variable neighborhood strategy adaptive search (VaNSAS) to solve the parallel-machine-scheduling problem in order to minimize energy consumption while considering job priority and makespan control. Although VaNSAS successfully improved solution-search performance in previous studies [17,22,[36][37][38], none had accounted for energy consumption, late delivery charge, and production overhead. The advantage of applying VaNSAS in this study was that its algorithms search for the best possible solution in many different areas by using several searching approaches, thereby moving to find more diversification and intensification at all times depending on the designed blackbox methods. In addition, we improved the solution-search performance of VaNSAS by presenting new neighborhood search strategies: (1) solution destroy and repair (SDR); (2) track-transition method (TTM); and (3) multiplier factor (MF). The proposed methods were performed in three sets of test instances and one case study, and compared with the exact method.
Computation results show that the proposed VaNSAS outperformed all traditional and other proposed methods. It performed as well as the exact method did, illustrated in the small problem, while the traditional DE could find only 25%; CU, SDR, TTM, and MF could not find the optimal solution at all, while VaNSAS could find 100% of the optimal solution. Even by increasing the problem size, VaNSAS still gave promising results, as it could improve the solution quality by 16.35-19.55% of the solution obtained from the exact method with 30-45 min of computation time, while the exact method required 48-72 h. In the medium-sized and large samples, DE, SDR, TTM, and MF could also improve the solution quality by 7.89-14.59% more than the exact method. In addition, the current practice gave a worse solution than that of the exact method, including all proposed methods, by 1.42-4.37%.
The excellent solution-search efficiency of VaNSAS in this study and its advantage of a flexible neighborhood search scheme allow researchers to further develop and design new mechanisms for improving solution quality. Additionally, it would be interesting to implement the proposed VaNSAS in various problem areas, such as production planning in a real manufacturing environment, while considering other additional factors or constraints, e.g., the study of the capability of each type of machine, job restriction, and customer-order scheduling.