An improved gray wolf optimization to solve the multi-objective tugboat scheduling problem

With the continuous prosperity of maritime transportation on a global scale and the resulting escalation in port trade volume, tugboats assume a pivotal role as essential auxiliary tools influencing the ingress and egress of vessels into and out of ports. As a result, the optimization of port tug scheduling becomes of paramount importance, as it contributes to the heightened efficiency of ship movements, cost savings in port operations, and the promotion of sustainable development within the realm of maritime transportation. However, a majority of current tugboat scheduling models tend to focus solely on the maximum operational time. Alternatively, the formulated objective functions often deviate from real-world scenarios. Furthermore, prevailing scheduling methods exhibit shortcomings, including inadequate solution accuracy and incompatibility with integer programming. Consequently, this paper introduces a novel multi-objective tugboat scheduling model to align more effectively with practical considerations. We propose a novel optimization algorithm, the Improved Grey Wolf Optimization (IGWO), for solving the tugboat scheduling model. The algorithm enhances convergence performance by optimizing convergence parameters and individual updates, making it particularly suited for solving integer programming problems. The experimental session designs several scale instances according to the reality of the port, carries out simulation experiments comparing several groups of intelligent algorithms, verifies the effectiveness of IGWO, and verifies it in the comprehensive port area of Huanghua Port to get the optimal scheduling scheme of this port area, and finally gives management suggestions to reduce the cost of tugboat operation through sensitivity analysis.


Introduction
Against the backdrop of open world ports and economic development, international trade competition has become increasingly fierce.The expansion of port scale and the continuous growing operational demands have become major trends and currents [1].The demand for tugboats during ship arrivals and departures at ports is increasing.At the same time, sustainable development and low-carbon environmental issues in maritime transportation have become hot topics of concern [2,3].
Tugboats are crucial for components in maintaining port operations, are costly and subject to increasing wear and tear, while port resources are limited [4].As a result, optimizing tugboat scheduling has become a crucial aspect to reduce port operating costs [5].By optimizing the tugboat scheduling process reasonably, a significant amount of tugboat resources can be saved, accelerating the efficiency of port vessel operations and reducing carbon emissions to a certain extent, thereby promoting energy conservation and environmental protection within ports.

Literature review
The tugboat scheduling problem (Tug-SP) originated from the Job-shop Scheduling Problem (JSP).Traditional tugboat scheduling models mainly rely on the experience of scheduling personnel, supplemented by some integer programming methods [6][7][8].This kind of scheduling model has more advantages in solving small-scale scheduling, but with the continuous growth of port operational demands and scale, traditional methods have become inadequate, necessitating more scientific and effective approaches to solve the tugboat scheduling problem.Current research focuses on establishing tugboat scheduling models that are in line with practical situations and employing suitable metaheuristic algorithms for solution.
Existing research related to tugboat scheduling can be broadly categorized into two phases, as shown in Table 1.In the first stage, researchers mainly concentrated in the initial exploration of intelligent tugboat scheduling research, placing greater emphasis on model construction while neglecting the actual needs of ports.For example, researchers such as Liu, Xu, and Kang described Tug-SP as a multiprocessor task scheduling problem, optimized the port tug tugboat scheduling problem based on the objective of maximizing the tugboat's running time using different intelligent algorithms [9][10][11], and Wang et al. investigated the tugboat schedule ng problem considering a multi-service model with multiple waypoints, with the optimization objective of minimizing the cost, and solved it with an adaptive large-domain algorithm and applicable to large-scale problems [12].This significantly promoted the study of Tug-SP models and represented a substantial improvement over traditional scheduling models.However, there are also limitations because, in practical port operations, tugboat efficiency is reflected by more than just maximum operating time.Port managers need to consider scheduling issues from various perspectives, including cost, efficiency, and environmental protection.A single objective function often fails to capture the complexity of real-world situations.
The second stage of the research takes into account other indicators reflecting the efficiency of the model, and the solution algorithms are gradually intelligent, but the applied algorithms have certain shortcomings in solving such problems as Tug-SP.Goli et al. studied the non-permutation flow-shop scheduling problem with the optimization objective of minimizing the completion time and the total energy consumed at the same time, and proposed Multi-objective Ant Lion Optimizer (MOALO), Multi-objective Keshtel Algorithm (MOKA) and Multiobjective Keshtel and Social Engineering Optimizer (MOKSEA) for solving the problem, and the experiments proved that these types of algorithms have good performance on the JSP [13].Tirkolaee et al. jumped out of the JSP and used different heuristic algorithms to solve the Municipal Solid Waste (MSW) and Location Allocation Routing Problem (LARP) respectively, considering multiple metrics such as cost and time, to achieve better results [14,15].Moharam et al. proposed a Discrete Chimpanzee Optimization Algorithm (DChOA) for the delay/loss (TL) penalty scheduling problem in discrete optimization problems, which is more advantageous than many other algorithms for solving the TL problem [16].Zhong et al. establish a dual-objective tugboat scheduling model with maximum completion time and total fuel consumption as the optimization objectives, and choose NSGA-II algorithm for the solution, which achieves better results on the example of Guangzhou port [17].
Due to the current deficiencies in solving high-dimensional problems and inadequate convergence accuracy in tugboat scheduling algorithms, we have noticed a remarkable new intelligent algorithm in the field of scheduling in recent years-the Grey Wolf Optimization (GWO) algorithm proposed by Mirjalili, a scholar from Griffith University, Australia, in 2014 [18].The advantages of this algorithm lie in its simple structure, few parameter settings, fast convergence speed, strong global optimization capability, and applicability to solving high-dimensional problems.As a result, GWO has been widely applied in various fields such as job-shop scheduling, engineering manufacturing, image information, and energy optimization [19][20][21][22].
However, GWO also has its limitations.It encounters conflicts when solving integer programming problems and is prone to getting stuck in local optima in the later stages of the algorithm.Many scholars have proposed improvements from various aspects such as control parameters, individual updates, and incorporating other strategies.Jordehi established an intelligent home energy management scheduling model and used GWO to solve it under the condition of integer encoding.By incorporating a penalty mechanism, the algorithm's convergence accuracy was improved, to some extent, mitigating the conflicts in GWO when solving integer problems [23].Liang et al. targeted the JSP and introduced a reverse learning strategy to increase population diversity.They combined it with simulated annealing to compensate for GWO's tendency to fall into local optima, thereby enhancing the algorithm's performance [24].Zheng et al. proposed an optimal-worst reverse learning strategy for optimizing scheduling algorithms, enabling timely escape from local optima [25].Furthermore, Zuo et al. built upon the optimal-worst reverse learning strategy and proposed a trend-optimizing reverse learning strategy for the Differential Evolution algorithm, further enhancing the algorithm's ability to escape from local optima [26].

Objectives and contributions
This paper investigates the problem of harbor tugboat scheduling.Effective tugboat scheduling can save harbor cost and improve harbor navigation efficiency.A multi-objective tugboat scheduling model is constructed under taking into account tugboat operation time, port operation costs and port navigation efficiency.The improved Gray Wolf Optimization Algorithm is used to solve the problem, and as a result, the optimal scheduling plan is obtained, which provides effective decision-making reference and optimization direction for port managers.
The main contributions of this paper are as follows: • This paper addresses the tugboat scheduling problem, which has been widely used in reality but has not received much attention.
• A new multi-objective tugboat scheduling model with optimization objectives that are more in line with realistic scenarios, including tugboat running time, tugboat fuel cost and tugboat power overflow, is proposed to better simulate the real situation.
• In order to solve the problem of integer coding conflicts and the tendency to fall into local optimization in GWO, the GWO algorithm is improved in this paper.This includes the introduction of a convergence parameter based on cosine mode, a dynamic weight adjustment mechanism with an inverse learning strategy for trend optimization to complement the current shortcomings of GWO.
• In the experimental part, the Taguchi method was tested to determine the optimal parameter combination of the algorithm.Simulation comparisons of the solution effects of IGWO are performed using multiple arithmetic cases of different scales by combining IGWO with CPLEX, Improved Particle Swarm Optimization (IPSO) [27], Discrete Chimp Optimization Algorithm (DChOA), Adaptive Large Neighborhood Search Algorithm (ALNS), GWO and IGWO are evaluated by comparing them and the results show that IGWO has better solution performance on TSP.
The rest of this paper is organized as follows: Section 2 presents the tugboat scheduling model, including a detailed model description, parameter settings, and value ranges.
Section 3 introduces the improved Grey Wolf Optimization algorithm (IGWO), including the principles of GWO, the rationale and effects of the proposed improvements, and the process of applying IGWO to solve the Tug-SP problem.
Section 4 designs numerous simulation experiments and applies multiple algorithms to solve the problem, providing a comparative analysis to validate the performance of IGWO.
Section 5 is designed to further validate the effect and improvement direction of IGWO, to obtain the optimal scheduling plan to provide some reference value for port scheduling, and to design a sensitivity analysis scheme to consider the effect of cost on the tugboat scheduling model.
Section 6 is to optimize the tugboat scheduling plan for the integrated port area of Huanghua port, and to test the sensitivity of the unit usage cost of two types of tugboats in this port area, so as to provide further adjustment experience for the managers.
Section 7 sets out the conclusions and outlook.

Problem description
Tugboat berthing and unberthing tasks can be divided into three stages: the approach stage, the assist stage, and the return stage.The approach and return stages are considered as the process in which a tugboat operates without any load, while the assist stage is considered as the operational process.Taking the example of two tugboats assisting in berthing, the specific process is illustrated in Fig 1.
The assumptions for constructing the model are as follows: 1. Task scheduling: Prior to the start of tugboat operations, the port receives the daily vessel operational requirements and specific berthing/unberthing details.Tasks are assigned to each vessel in chronological order, taking into account the available number of tugboats, their power, and other relevant information.Tasks that are close in time (with similar start and end times) can be combined, utilizing shared tugboats.
2. Tugboat operations: After completing the assist stage of a task, the tugboat randomly returns to a base and waits for the next task (which may not be immediately following).
Returning to the nearest base may not be advantageous for achieving the globally optimal solution.
3. Tugboat bases: Initially, the tugboats are randomly distributed among the tugboat bases.
The time required for tugboats to enter or exit the bases is negligible.
4. Tugboat power: Each tugboat has a fixed power for propulsion.Tugboats with higher power have higher empty cruising speeds and fuel consumption.The cost of empty cruising per unit distance and the cost of assistance per unit time are also higher.

Navigation rules:
The port has a one-way channel, and vessels cannot sail in parallel.
6. Task time windows: Based on the task scheduling information, the time interval for each vessel to complete the berthing/unberthing can be determined.

Parameter settings
The parameter settings for the tugboat scheduling problem, along with their descriptions, are presented in Table 2, while Table 1 is a two-tier table, where the first and third columns represent the parameter in the model, and the remaining two columns provide corresponding explanations.In addition, the individual tugboats in the model are meant to be similar tugboats and have different performance in terms of power, speed, and fuel consumption.In addition, the individual tugs in the model are meant to be similar tugs and have different performance in terms of power, speed, and fuel consumption.L1 k ij is related to the current position of the tugboat and the position of the target ship entering the harbor, and L2 k ij is related to the current position of the tugboat and the berth where the ship is docked.

Model establishment
Firstly, in order to improve the efficiency of port operations, it is necessary to minimize the operational time of tugboats during the berthing and unberthing process of vessels.The operational time of tugboats in the port refers to the total time spent on both idle operations and assisting operations for completing all tasks.The first objective function is shown in Eq (1).
Additionally, the control of fuel costs for tugboats is considered.When tugboats assist vessels in berthing and unberthing operations, the main costs are based on power of tugboats and the duration of task execution.Therefore, in this study, the cost during the assistance phase is measured based on the assistance cost of different tugboats.The remaining phases involve unloaded travel and do not involve assistance phase, which can be measured by the cost per unit distance traveled by the tugboats.The cost of tugboat operations (during the assistance phase) is influenced by two factors.On the one hand, it depends on the power of the operating tugboat, with higher power resulting in higher costs.On the other hand, it is determined by the duration of the operation, which is governed by the berthing and unberthing schedule of the target vessel.The second objective function is shown in Eq (2).
Furthermore, the overflow power of the tugboats is calculated.For tugboat-assisted berthing and unberthing tasks, tugboats should arrive within the time window for the target vessel to avoid delays in cargo handling and occupation of the designated berth, resulting in demurrage fees.Early arrival has less impact on the port, therefore, the waiting cost is not considered in this study.The third objective function is shown in Eq (3).
In summary, considering the optimization objectives of minimizing the maximum  (4).
The constraints of the model are represented by Eqs ( 5) to (14), Eq (5) represents the number of tugs allocated to each task.Eq (6) indicates the initial positions of the tugs before the first task.Eq (7) states that the tugs completing a task must meet the minimum power requirement for the current task.Eq (8) means that the tugs assigned to the first task should arrive at the starting point of the task before the start of the first task.Eq (9) denotes that after completing the previous task, the tugs should return to the docking base and reach the starting point of the current task in time.Eq (10) represents the requirement for the tugs to return to a specific base before proceeding to the next task.Eq (11) indicates that the number of tugs at a base should meet the required number of tugs to complete the tasks.Eq (12) states that a tug can only dock at one base at a given time.Eq (13) implies that in the absence of any tasks at the current stage, the tug remains at the previous base.Eq (14) signifies that the variables are integer values ranging from 0 to 1.

Principles of gray wolf algorithm
The GWO algorithm is characterized by simple structure, few parameters to be adjusted, easy to implement, etc., in which there are convergence factors that can be adjusted adaptively as well as the information feedback mechanism, which can achieve a balance between local optimization and global search, and GWO has a good performance in terms of solution accuracy and convergence speed.
The solution of multi-objective problems is often limited by the complexity of the solution space, while the global optimization ability and high-dimensional solution performance of GWO can compensate for this limitation to a certain extent, so that GWO has better robustness and applicability in multi-objective problems after parameter adjustment.
(1) Social hierarchy GWO is optimized by simulating the predation behavior of grey Wolf groups.Grey Wolf groups have a strict social hierarchy, as shown in Fig 2, which is divided into four classes.The specific meanings of each social hierarchy in the figure are as follows: Level 1: Alpha wolf (α).The leader of the pack is responsible for leading the pack to hunt prey, which is the optimal solution in the optimization algorithm.
Level 2: Beta wolf (β).Assisting the alpha pack, the sub-optimal solution in the optimization algorithm.
Layer 3: Delta wolf (δ).Obey the orders of alpha and beta, responsible for reconnaissance, lookout, etc.The α and β of fitness differences will decrease to δ.
(2) Hunting process GWO mimics the predatory characteristics of a gray wolf pack, with the goal of tracking, surrounding, chasing, and attacking prey to achieve optimal search.The hunting process of gray wolves is implemented based on Eqs ( 15) and ( 16): While t represents the current iteration count, X !P represents the position vector of the prey, X !ðtÞ represents the position vector of the current gray wolf; A and C are cooperative coefficient vectors that vary according to Eqs (17) and (18): During the hunting process of the wolf pack, a linearly approaches 0 from 2; r1 and r2 are random vectors in the range [0, 1].Clearly, as the iteration progresses, a gradually decreases to 0, causing the gray wolves to approach the prey.Simultaneously, r1 and r2 provide the opportunity for the wolves to escape from local optimal solutions.

Convergence parameters based on cosine variation
The magnitude of the cooperative coefficient vector A directly affects the overall and local optimization capabilities of the algorithm.Eq (17) indicates that A is controlled by the convergence parameter a, which linearly decreases from 2 to 0 with the number of iterations.However, during the operation and search of the algorithm, this linear variation does not fully reflect its optimization process.Hou et al. improved the gray wolf optimization algorithm by introducing a non-linear convergence parameter [28].Based on this theory, this study introduces a convergence parameter expression based on cosine variation: While a initial and a final respectively denote the initial and final values of the convergence parameter a, while n represents the decreasing index, 0�n�1.The convergence parameter varies over time as shown in Fig 3 .From the change in the absolute value of the slope in the graph, it can be observed that the descent pattern of this convergence parameter is slow-fastslow.This indicates that the improved convergence parameter decreases slowly in the initial stage, enabling a global search over a wide range.In the middle stage, it decreases rapidly, allowing faster encirclement of the prey and positioning in the vicinity of the optimal solution.In the final stage, it decreases slowly, enhancing the algorithm's local search capability and improving algorithm accuracy.

Dynamic weight correction and update mechanism
In the original Grey Wolf Optimization (GWO) algorithm, each individual has the same learning degree towards the alpha, beta, and delta wolves, which hinders the ability of the new generation individuals to become more excellent through learning from the alpha wolf.To address this limitation, Meng et al. proposed a nonlinear dynamic learning weight mechanism [29].however, it still couldn't fully exploit the learning capability of individuals.In practical situations, the excellence of alpha, beta, and delta wolves is ranked as α>β>δ.Therefore, a dynamic weight mechanism based on step length Euclidean distance is introduced initially to dynamically adjust the learning weights according to the differences among the alpha, beta, and delta wolves.The specific formulas are shown in Eqs (18) and (19).
The tugboat scheduling problem is typically addressed using real-number encoding.However, the updating process in the Grey Wolf Optimization (GWO) algorithm often involves non-integer solutions.Moreover, the solution individuals in tugboat scheduling problem are subject to constraints, such as the non-repetition of tugboats for the same task and the minimum power requirement for each task.After the application of the GWO algorithm, the individuals may easily produce infeasible solutions.To address this issue, a cross-repair mechanism is proposed in combination with the aforementioned dynamic weight adjustment mechanism to optimize non-integer solutions.The specific steps of the dynamic weight correction and updating mechanism are as follows: 1. Based on the positions of the current wolf and the contemporary alpha, beta, and delta wolves in the population, as well as their respective cooperative coefficient vectors A and C.
2. Calculate the distance vectors between the current wolf and each of the alpha, beta, and delta wolves, and perform preliminary individual updates based on the dynamic weight mechanism using Formulas (18) and (19).
3. Correct the preliminary updated individuals to integer solutions (which are likely to still be infeasible solutions).4. Obtain the set of eligible tugboats based on the minimum power requirement for each task.
5. Take the intersection between the infeasible integer solutions obtained in step 3 and the set of eligible tugboats from step 4, and select the individual that is closest to the infeasible integer solution to complete the update of the current individual.If there are missing parts in the individual, they are randomly selected from the remaining tugboats.
6. Perform the final update on all individuals in the Grey Wolf population following the above steps, thereby completing the overall update of the Grey Wolf population.
Taking the example of 5 tasks and 8 tugboats, Fig 4 demonstrates the process of cross-repair updating.The individual is composed of the tugboats required for all tasks.In this process, the current individual is cross-fused with the parameters A, C, and the contemporary alpha, beta, and delta wolves of the current generation.After rounding, the next generation's Grey Wolf position vectors X1, X2, and X3 are obtained.The vectors are then centered and adjusted, and their intersection with the set of eligible tugboats is taken (with the remaining missing parts of the individual randomly selected).After multiple individual updates, the final result is a new generation of Grey Wolf individuals.
The application of the dynamic weight correction updating mechanism can make the new generation of individuals self-adaptive learning, and the learning rate of better wolves is

Towards excellent inverse learning strategy
The inverse learning strategy, introduced by Tizhoosh in 2005 [30], is a novel algorithm that enhances the global optimization capability of various algorithms, thereby improving their search efficiency.Building upon the theory of elite inverse learning, Li et al. incorporated a reverse learning mechanism for the best and worst individuals [31], as calculated in Eq (20).
where x ij represents the current individual, lb represents the upper bound of the current individual, ub represents the lower bound of the current individual, and x r represents the reverse solution.
Ordinary elite inverse learning strategies only apply reverse learning to the best individuals in the population.However, the worst individuals in the population also contribute to diversity accumulation.Considering the characteristics of the Grey Wolf Optimization (GWO) Their respective inverse solutions are calculated using Eq (23).These reverse solutions are then sorted based on their fitness values, and the lower half of the sorted individuals with the lowest fitness values are selected as new solutions.The algorithm randomly replaces individuals in the wolf pack with these new solutions, thereby enhancing the algorithm's ability to escape local optima.

IGWO solving process
Based on the above theoretical methods, this paper adopts the new IGWO algorithm to solve the tugboat scheduling model.The specific solution process is shown in Fig 7 .The solution process is divided into four stages, which are as follows: 1.The first stage: First obtain the relevant information of the port task, and enter the task data, including the number of tugboats required for each task, the required working time, and the required tugboat power.The second is the data of the existing tugboats in the port, including the number of tugboats, power, idle speed, idle cost and navigational aid cost, two of which are jointly set by local regulations, the Maritime Authority and the port company.Finally, the data of the tugboat base, including the number and capacity of the tugboat base, as well as the location of the base.
2. The second stage: tasks will match a certain number and power of tugboats, which restricts the use of tugboats, and according to the ordering of tasks and the constraints of model conditions, the set of tugs currently available is obtained.Finally, according to the current available tugboat set, the initial grey Wolf population is randomly generated, which ensures that there will be no illegal solution in the initial solution.
3. The third stage: In the process of solving the gray Wolf algorithm, it is necessary to calculate the fitness of individual gray wolves in the gray Wolf population.The closer the individual is to the objective function, the greater the fitness is.All individuals are sorted according to the fitness, and the three gray Wolf individuals with the largest fitness are screened out, which are the first wolves of this generation α,β,δ.Next, all individuals need to learn from the modern leading Wolf to update themselves, and the updating method adopts the 4. The fourth stage: the update in the third stage satisfies the convergence process of the gray Wolf algorithm, but the algorithm is still prone to local optimization.Therefore, the optimal reverse learning strategy is adopted in the fourth stage to update the gray Wolf population twice.Although the optimal reverse learning strategy can help the algorithm jump out of the local optimal, it is prone to illegal solutions after updating, so the illegal solutions are corrected at last.The specific method is to take the intersection of the available tugboat set and the current gray wolf population to generate a legal new solution.If the length of the new solution is insufficient, it is randomly selected from the set of available tugboats to make up for it.

Tugboat scheduling case study simulation
To validate the effectiveness of the IGWO algorithm in solving the tugboat scheduling problem, simulation experiments are conducted to compare different instances and multiple intelligent algorithms.The code for each experiment is implemented using Python 3.7 (with the CPLEX library called through Python functions).The simulation runs on a laptop with an 11th Gen Intel(R) Core(TM) i5-11400H @ 2.70GHz processor, 16GB RAM, and Windows (64-bit operating system).

Experimental parameter settings
To assess the performance of IGWO in harbor tugboat scheduling, we generate tugboat scheduling cases based on different scales.The scale of each instance is controlled by the number of tasks and tugboats.In addition to the parameters controlling the scale mentioned in Table 2, there are several common parameter ranges for all cases.These ranges are used to define constraint conditions and calculate various cost expenditures during the scheduling process.The parameter ranges are determined by referring to relevant data on China's port charging standards [32], taking into account the specific model proposed in this study and the actual conditions of the one-way channel in the harbor.The parameter ranges are shown in Table 3. Table 3 shows the ranges of some parameters.Columns 1 and 3 are the parameters, columns 2 and 4 indicate the parameter ranges, and columns 3 and 6 denote the units of the parameters.The explanations for the parameters have been given earlier in the text and will not be repeated here.

Feasible solution setup
Given the numerous constraints in the model, it is necessary to set up the generation of feasible solutions.The constraints to be considered can be divided into three aspects: the demand of each task for tugboats, the operational limitations of the tugboats themselves, and the capacity limitations of the tugboat base.The setting of feasible solutions for the multi-objective tugboat scheduling model is shown in Algorithm 1, which first obtains the task ship and tugboat data through the port, and generates the set of available tugboats through this data to clarify the number of tugboats required for each task.On this basis, a single feasible solution is generated, and the population of feasible solutions is generated according to the population size parameter N P .

Algorithm parameter optimization
The parameters that can be adjusted in the IGWO algorithm include: the population size N P , the number of iterations N t , and the convergence parameter a, in which the convergence parameter a has been adjusted based on the cosine law, so the next two parameters can be adjusted.
Different combinations of algorithm parameters will affect the solution performance of the algorithm, so the Taguchi method is used to test to determine the best combination of parameters.One of the two parameters, three levels of parameter combination experiments are shown in Table 4.
Next, the results of each parameter combination were tested at three different scales, respectively, and the scale of the test cases is shown in Table 5.The test results were measured using the signal-to-noise ratio (SNR), which is calculated as shown in Eq (24): Where y ij represents the test results of each group.Orthogonal experiments were carried out according to three levels of each parameter, and each group was tested 10 times, and the test results are shown in Table 6, and it can be seen by observing the table that the SNR value reaches the minimum in the 5th group of experiments, so this paper sets the algorithm parameters as the number of populations N P = 200, and the number of iterations N t = 1000.

Small-scale simulation experiments
To verify the applicability of the IGWO algorithm in small-scale tugboat scheduling problems, 15 sets of different small-scale cases were designed for solving using the IGWO algorithm.
A comparison experiment was conducted with CPLEX, a mathematical solver, as shown in Table 7.The first column of the table represents the 15 sets of generated instances.The second column indicates the number of tasks in each instance, which is slightly larger than the number of target vessels since each task requires 1-5 tugboats.On average, the number of target vessels is 2.5 times the number of tasks.The third column represents the total number of available tugboats for each instance.The fourth and fifth columns show the results of CPLEX and IGWO, respectively, for each instance.The last column presents the deviation between the two solutions.
Among all the small-scale cases, IGWO achieved the same optimal value as CPLEX for the first 6 instances, with a deviation of 0.00%.For instances 7-12, the solutions of both methods were very close, with a deviation controlled within 0.9%.In the last 3 instances, CPLEX encountered a memory overflow issue and failed to obtain a solution, while IGWO's solving process remained unaffected.At present, most of the ports in reality are small in scale, and it can be seen from the above examples that IGWO has good solving performance in the problem of small-scale tugboat scheduling.Therefore, the IGWO proposed in this paper has important practical significance.

Medium to large-scale simulation experiments
Similarly, we need to go to observe the performance of IGWO in medium and large-scale instances, because CPLEX cannot solve larger-scale problems, so the other three improved intelligent algorithms as well as the original Gray Wolf optimization algorithm were selected for experiments, designed 10 groups of medium and large-scale instances, and analyzed the strengths and weaknesses of the algorithms based on a number of indicators, as shown in Table 8.Table 8 consists of the following columns: Instance Number, Task Quantity, Tugboat Total, and Solution Metrics including best value (Bes.), average value (Avg), worst value (Wor.) and solution time (Time.).Columns 5 to 8 represent the solution results of each intelligent algorithm for the given model, and column 9 indicates the gap in the optimal value (Gap.) between IGWO and other algorithms.Gap. is calculated using Formula (25), where a larger Gap.indicates the superiority of IGWO over other algorithms.
Comparing the simulation results of different algorithms at different scales, it can be seen that the IPSO algorithm enables the algorithm to complete the fine search at a later stage through the nonlinear inertia weighting strategy, but it is limited by the algorithm's own problem that the global optimal object of individual learning is only a single individual, which is not conducive to searching in the high-dimensional solution space.For the DChOA algorithm, although the algorithm improves the parameters and operators, the problem solving effect is not ideal for high-dimensional solution space, it can be seen that as the size of the solution increases, the dimensionality of the solution space is also increases, and the solving accuracy of DChOA highlights the lack of precision, and it is more likely to produce the phenomenon of premature maturity.The ALNS algorithm, due to the improved domain search operator, has a better performance compared to the previous algorithms, but the convergence speed is slower, and the GWO algorithm can easily fall into the local optimum without adjusting the convergence parameter alpha and the learning strategy, resulting in a lower solution accuracy.For the IGWO algorithm, observation of the experimental results shows that in examples 2-5, due to the relatively small size of the instances, the solution accuracy of IGWO is 10%-14% ahead of the other algorithms, and the remaining 7 instances still maintain a 2%-10% lead.In addition, the solution time of IGWO is significantly shorter than other algorithms, and its solution speed dominates at all scales.This indicates that IGWO is more stable and outperforms all other compared algorithms, including the original GWO algorithm, in terms of overall performance.This is due to the tuning of the convergence parameters and the dynamic learning strategy, the better balance between the algorithm's global and local convergence capabilities, and the introduction of the inverse learning strategy, which enables the algorithm to jump out of the local optimum in time.

Optimal solution tracking
Section 3 verified IGWO's solution accuracy and stability was accomplished.In order to further enhance the algorithm's performance and apply it to practical port operations, this Section presents a redesigned set of four instances with a large span of scales.The algorithm parameters remain the same as in Section 3. The aim is to track the optimal values and corresponding scheduling solutions of each intelligent algorithm, obtain their iteration curves, analyze the directions for improving IGWO, and ultimately provide a reference for practical port operations based on the optimal scheduling solution obtained through IGWO.

Validation of optimal solutions
As shown in The IGWO algorithm outperforms other comparative algorithms in terms of solution accuracy but has some limitations.Further improvements could focus on enhancing its ability to escape local optima in the later stages of the algorithm.

Optimal scheduling solution
To apply this algorithm to practical port operations, case (a) was selected as an example.The Gantt chart in Fig 9 represents the optimal scheduling plan generated by IGWO, providing a clear view of the assigned tugboat numbers and the sequence in which they perform tasks.It serves as a reference for tugboat scheduling in port operations.For example, tugboats assigned to serve Task 1 are Tugboat 1, Tugboat 14, Tugboat 15, and Tugboat 17.The sequence of tasks performed by Tugboat 12 is Task 2-Task 7-Task 18-Task 19.Additionally, it can be observed that Tugboat 1 is frequently assigned tasks, while Tugboats 3-8 remain idle.This is related to the existing tugboat capacities and task requirements.To address this, further optimization can be achieved by introducing constraints on continuous tugboat operations in the model.

Sensitivity analysis
This section focuses on the sensitivity of the objective function to the tugboat's unloaded travel cost and pilotage cost.The optimal value (Bes.) is used as the evaluation criterion.The simulation experiments in this chapter continue to utilize the four instances introduced in Section 4.
Firstly, a sensitivity coefficient is defined with a range of [0, 2] and a step size of 50%.It is used to test the impact of tugboat's unloaded travel cost (Ct i ) and navigation assistance cost (Co i ) on each set of cases.The sensitivity performance of Ct i and Co i in the four sets of cases is shown in Fig 10,while (a) shows the curve of the optimal value with respect to the sensitivity coefficient of idle cost, while panel (b) depicts the curve of the optimal value with respect to the sensitivity coefficient of navigation assistance cost.
From Fig 10 it can be seen that sensitivity is directly proportional to the slope of the curve.Analyzing Fig 10-(A), it can be observed that in small-scale instances, changes in Ct i do not significantly affect the variation of the optimal value, only exerting a relatively larger impact on the optimal value in larger-scale cases.Conversely, in Fig 10-(B), Co i has caused significant disturbances to the optimal value across instances of various scales.In conclusion, the tugboat scheduling model established in this paper is significantly influenced by the tugboat assistance cost, and it is particularly affected by the tugboat unloaded travel cost in large-scale scenarios.This provides some recommendations for port managers.They should strive to minimize the tugboat assistance cost while also considering the unloaded travel cost, especially in large-scale ports, in order to achieve cost control and operational efficiency.

Huanghua port example study
In order to validate the effectiveness of the algorithm and the model, an example study is conducted with the integrated harbor area of Huanghua Port in China, which is about 9.5 nautical miles in length and is now in operation with 11 berths, 1 tugboat base, and 2 mission start points.
The data of May 1, 2022 is selected for the study of this port area, and the information of the tugboat configuration of this port area and the incoming ship data of that day are shown in S2 Table, the tugboat scheduling of that day in this port area is optimized by using the IGWO algorithm, and the parameters of the algorithm are the same as in the previous section, which results in the optimal operation scheme as shown in Table 9.The 21 vessels on that day are divided into 8 tasks, and the optimal order of tugboats to perform each task can be clearly seen through this table.In addition, the same method was used to test the sensitivity of tugboat idling cost Ct i and aids to navigation cost Co i to the optimal value Bes. as shown in Table 10, where the first column of the sensitivity coefficients indicates the proportion of the change in the two factors, and the second and third columns indicate the extent to which the change in Ct i and Co i , respectively, affects the optimal value (Bes.).In this example, the effect of the change in Ct i is small, while the effect of Co i is large.Therefore, if conditions permit, port managers should try to reduce the cost of tugboat navigation assistance Co i , which can more efficiently reduce the cost of port operation as well as improve the efficiency of port navigation.

Conclusion and future work
This paper proposes a novel tugboat scheduling model based on different power levels of tugboats and vessels.Multiple objective functions are established, including in-port operation time, various tugboat costs, and overflow power.The study explores the application of the Improved Grey Wolf Optimization (IGWO) algorithm to the tugboat scheduling problem.The algorithm's performance is verified through simulation experiments at various scales.The results show that IGWO optimizes several other comparative algorithms, and finally the research in this paper is applied to the integrated port area of Huanghua Port to obtain the optimal scheduling scheme with effective management insights.
The proposed tugboat scheduling model in this study simulates the real-world scenario of tugboats assisting ships in entering and leaving ports.Through IGWO optimization, it generates efficient tugboat scheduling schemes that lead to fuel and manpower savings for ports while also reducing carbon emissions to some extent.This offers new insights and methods for promoting the sustainable development of maritime transportation.
This paper establishes a new tugboat scheduling model based on the task scheduling theory, but the model and algorithm can be further optimized in the future, model-wise the objective function can be replaced according to the actual needs of the port and the constraints of continuous operation of tugboats and inter-port scheduling can be added.In terms of algorithms, deep reinforcement learning is noteworthy, but such algorithms have high environmental requirements, and whether they are applicable to the tugboat scheduling problem needs to be studied.

Fig 8 ,
while (a), (b), (c), and (d) correspond to instances with task-tractor ratios of 20-32, 40-48, 60-72, and 80-96 respectively.According to the iteration curve, it can be observed that both HNGA and GWO algorithms suffer from the disadvantage of easily converging to local optima.Although IPSO improves the inertia weight, its performance is poor for moderate-scale problems.

S1
Appendix.Port billing schemes.The S1 Appendix describes the rules and regulations related to port billing in China.(DOCX) S1 Table.Tugboat scheduling simulation dataset.The S1

Table 2 . Parameters of the tugboat scheduling model. Parameter type Parameter Description Set
tugboats, fuel costs associated with idle and assisting operations, and the overflow power of tugboats, a new multi-objective tugboat scheduling model is established.The overall objective function of the model is shown in Eq I Set of tugboats, i2I = {1,2,3,. ..,n}J Set of tasks, j2J = {1,2,3,. ..,m}K Set of bases, k2K = {1,2,3,. ..,l} i Location of tugboat i at base k after completing the previous task st i Time of arrival of vessel for task j at the starting point et j Time of arrival of vessel for task j at the ending point v i Empty cruising speed of tugboat i Decision variables x ij 0-1, 1 if tugboat i is assigned to task j, otherwise 0 y k ji 0-1, 1 if tugboat i is at base k after completing task j, otherwise 0 https://doi.org/10.1371/journal.pone.0296966.t002operational time of

Table. Data of Huanghua Port Comprehensive Port Area
. The S2 Table contains three sub-tables representing the basic data of the port area, the ships on a particular day, and the pre-processed ship data.(XLSX)