https://orcid.org/0000-0001-8102-8846
Figures
Abstract
Recently, there has been considerable research on combining multi-agent simulation and particle swarm optimization in practice. However, most existing studies are limited to specific engineering fields or problems without summarizing a general and universal combination framework. Moreover, particle swarm optimization can be less effective in complex problems due to its weakness in balancing exploration and exploitation. Yet, it is not common to combine multi-agent simulation with improved versions of the algorithm. Therefore, this paper proposes an improved particle swarm optimization algorithm, introducing a multi-level structure and a competition mechanism to enhance exploration while balancing exploitation. The performance of the algorithm is tested by a set of comparison experiments. The results have verified its capability of converging to high-quality solutions at a fast rate while holding the swarm diversity. Further, a problem-independent simulation-optimization approach is proposed, which integrates the improved algorithm into multi-agent systems, aiming to simulate realistic scenarios dynamically and solve related optimization problems simultaneously. The approach is implemented in a response planning system to find optimal arrangements for response operations after the Sanchi oil spill accident. Results of the case study suggest that compared with the commonly-used shortest distance selection method, the proposed approach significantly shortens the overall response time, improves response efficiency, and mitigates environmental pollution.
Citation: Du S, Fan W, Liu Y (2022) A novel multi-agent simulation based particle swarm optimization algorithm. PLoS ONE 17(10): e0275849. https://doi.org/10.1371/journal.pone.0275849
Editor: Ali Safaa Sadiq, Nottingham Trent University School of Science and Technology, UNITED KINGDOM
Received: May 31, 2022; Accepted: September 25, 2022; Published: October 13, 2022
Copyright: © 2022 Du et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: All relevant data are within the paper.
Funding: The author(s) received no specific funding for this work.
Competing interests: The authors have declared that no competing interests exist.
Introduction
The merging of optimization and simulation technologies has seen rapid growth recently. The optimization of simulation models refers to finding the best values of some decision variables for a system where the performance is evaluated based on the output of a simulation model of this system [1]. However, due to the complexity, analytical simulation methods may be unsuitable for large-scale systems in manufacturing, supply chain management, financial management, etc. Agent-based modeling (ABM) is a bottom-up modeling method with high autonomy and interactivity, which is particularly suitable for complex system research [2]. Given the advantages of ABM, researchers have been solving optimization problems by integrating multi-agent simulation (MAS) with a wide range of methods, including game theory, reinforcement learning, and swarm intelligence (SI). Applications can be found in various fields such as intelligent transport systems, efficient electric grids, and smart buildings. Among the methods mentioned above, SI has been proved to be robust and efficient for optimization problems. As a classic category of SI, particle swarm optimization (PSO) has the advantage of simple parameter setting, fast convergence rate, high stability, and strong scalability, making it an excellent option for complex systems.
However, despite a high convergence rate, the standard formation of PSO (SPSO) [3] is likely to stagnate at a local optimum, which reduces the solution accuracy. Thus, it is necessary to modify the algorithm for the application in complex problems where stable optimal solutions are indispensable. Different approaches have been considered to improve SPSO, among which setting parameters is especially representative [4]. Changing coefficients is a typical way of parameter setting. The performance of SPSO is affected by the values of its coefficients, i.e., the acceleration coefficients and the inertia weight. Thus, several attempts were made to tune the values of these parameters, such as the random weight method (RPSO) [5], the constriction factor approach (CPSO) [6], and the fuzzy adaptive approach (FAIPSO) [7]. Apart from the coefficients, the population size of the swarm is also a vital parameter. The basic idea for adjusting population size is to concentrate the individuals in the most promising area during the exploitation phase of the algorithm. De Oca et al. [8] conducted this idea by introducing incremental social learning (IPSO). On the other hand, as SPSO came from a simplified simulation of bird clustering behaviors, some researchers tried to improve it by mimicking other bird activities simultaneously. Neshat et al. [9] modified SPSO based on the predation phenomenon in bird swarms (PPSO). Besides these approaches, Li et al. [10] decoupled exploration and exploitation to make the algorithm more effective in large-scale optimization.
In MAS, agents can be either homogeneous or heterogeneous. Homogeneous agents are usually treated as particles when PSO is applied [11]. Yang et al. [12] presented a multi-agent-based model which simulated human behaviors in a multi-exit evacuation environment, where individuals were homogeneous agents and PSO was introduced to simulate their movement. Ali et al. [13] proposed a collective motion and self-organization control of a swarm of 10 unmanned aerial vehicles, which were divided into two clusters of five agents each. PSO was adopted to provide the best agents of the cluster. While homogeneous systems often simulate certain collective activities of agents with the same properties, heterogeneous systems are more practical for realistic scenarios due to the diversity of agents. Kanaga and Valarmathi [14] proposed a multi-agent model using PSO to solve patient scheduling problems, using a PSO Agent to perform the algorithm. Yu et al. [15] designed a Web services selection model, using PSO to select the optimal Web services combination in E-business. The model was realized by MAS. Thiel et al. [16] simulated Hanoi dwellers’ choices between traditional markets and a hypothetical new one, aiming to find the optimal location for the ready-to-build local supermarket by using PSO to maximize sales volume. Mahad et al. [17] proposed an intelligent MAS approach to optimize power supply in community-based multi-microgrids systems, where various layers of autonomous and intelligent agents took decisions based on the PSO method.
Although there are many studies on improving SPSO, there are not so many attempts to combine proposed algorithms with MAS and apply them to practical problems. In fact, utilizing improved algorithms in simulation-optimization systems is necessary because complex tasks require high-performance optimizers, and SPSO may not meet this requirement. Besides, existing studies on the combination of MAS and PSO may vary in engineering fields, yet they all designed the system based on specific problems of their own research field without summarizing a general and universal framework, which is field-independent and can be implemented in various settings. Therefore, this paper focuses on a more in-depth exploration of improving PSO and integrating it into MAS to solve practical problems. The main contributions of this paper are as follows:
- Proposing an improved particle swarm optimization algorithm (CoPSO) according to the biological nature of birds. It has a multi-level structure and a competition mechanism, aiming to improve performance by enhancing exploration while balancing exploitation.
- Proposing a general simulation-optimization approach (MAS-CoPSO) to integrate CoPSO into MAS. It is supposed to achieve a very close combination between simulation and optimization, and the application should not be limited to certain engineering fields and specific problems.
- Conducting a series of experiments to substantiate the validity of CoPSO.
- Implementing MAS-CoPSO in a case study of the Sanchi oil spill accident to solve the response planning problem.
Methods and materials
This section elaborates on the methodologies used in this study. Firstly, a brief introduction to MAS technology is given. Detailed descriptions of the proposed CoPSO algorithm and the MAS-CoPSO approach follow afterward.
Multi-agent simulation
The fundamental element in MAS is the agents, which can be software or physical entities. Although different agents exhibit distinct behaviors, they share some common properties. For instance, they all have a certain degree of autonomy that enables them to work even without human intervention, they can communicate and interact with each other, and they are capable of perceiving and reacting to the changes in the environment as well as determining the proper behaviors to achieve the final goal [18]. In a simulation system, agents can be either homogeneous or heterogeneous depending on the specific problem background. Heterogeneous agents can deal with diverse scenarios and complex tasks, making it more suitable for complicated realistic systems [11]. Based on the aim of proposing a general framework for the combination of MAS and CoPSO, this paper mainly focuses on heterogeneous systems due to their universality.
Improved particle swarm optimization algorithm
PSO mimics the behavior of bird swarms. The fundamental idea is to assume a potential solution to the optimization problem as a bird without quality and volume, i.e., a particle. Particles are described by their positions and velocities. Generally, the position of a particle represents a specific solution while the velocity represents the searching direction and scope, which leads the particle to fly through the solution space. Every particle will get a fitness value by calculating the objective function. It will also adjust the value according to the experience of itself and the neighbors.
Supposing the search space is D-dimensional, the ith particle is represented as , where
, d ∈ [1, D],
and
are the lower and upper bounds of the dth dimension, respectively. The velocity of the ith particle is represented as
, where
, d ∈ [1, D], Vmin and Vmax are the minimum velocity and maximum velocity specified by the user, respectively. In SPSO [3], during each iteration, particle i updates its position and velocity as follows:
(1)
(2)
where w is the inertia weight, r1 and r2 are random values between 0 and 1, c1 and c2 are acceleration constants which control how the particle moves,
is the best previous position of the ith particle, and gbest = (gbest1, gbest2, …, gbestD) is the best position among all the particles in the swarm.
PSO encourages exploration during the early stage of the iterations while encouraging exploitation in the latter stage. However, Angeline [19] has pointed out that the balance between exploration and exploitation is subtle and difficult to control in PSO, leading to stagnation or premature convergence. Inspired by the features of bird flocks, this paper proposes CoPSO to improve the performance of PSO, introducing a multi-level structure and a competition mechanism.
Multi-level structure.
According to the biological nature of birds, their ability to sing and create a complete song is learned from their parents. Experiments have shown that if a bird is reared in silence, it can only scream. In addition, birds’ singing skills are not set in stone. Basically, as they grow up, their expertise in singing also gets improved [20]. Based on this phenomenon, CoPSO has a multi-level structure, dividing the particles in a swarm into different singing levels. It is stipulated that the closer a particle is to a possible global optimum, the higher its level. The level of the ith particle is denoted by li, which will increase if the particle doesn’t update its personal best position pbesti in μ continuous iterations. Since birds of different ages learn to sing at different rates, particles of different levels have different acceleration constants.
Fig 1 illustrates how the particles’ levels change during the iterations in two-dimensional unimodal problems. At the early stage of the algorithm, particles are like newborn birds with little knowledge about singing. They will explore the search space and try to find promising areas. Once found, they might temporarily stop updating their personal best positions, leading to an upgrade. At the latter stage of the algorithm, some particles can be very close to the global optimum. Consequently, they will find it difficult to locate better positions, and their levels are also relatively higher.
The points represent the particles. (a) The early stage. (b) The latter stage.
Competition mechanism.
While a particle’s high level is likely to come from finding the global optimum, there are exceptions. Note that falling into local optima can also stop particles from updating their personal best positions. Since particles’ levels are linked to their ages (the higher level and the more skillful in singing, the older), the idea of the competition mechanism inside bird swarms is utilized to improve the algorithm. The survival resources of swarms and the lifespan of birds are both limited. Therefore, a bird swarm is constantly renewed by the death of old birds and the born of new ones. If a particle’s level is too high, it becomes too old to possess the resources of the swarm. The population size is set as a constant, once a particle’s level reaches the upper bound Lmax, it will be replaced by a randomly generated new one. The new particle has the lowest level since it’s a newborn.
Fig 2 depicts how the competition mechanism functions in two-dimensional multimodal problems. If particles prematurely converge to a local optimum, they will upgrade around it. It can be seen from Fig 2a that few particle gets close to the global optimum. However, once the particles in the local area become too high-leveled, they will be replaced by new birds. As these new particles are randomly generated, they have the potential to explore unknown search space and finally locate the global optimum, which is shown in Fig 2b.
The points represent the particles. (a) Particles are trapped in the local optimum. (b) Particles jump out of it.
Since particles tend to move closer to the global optimum and oscillate around it as the algorithm runs, in practice, the value of μ will be increased as the level rises, which allows an elaborate control over the balance between exploration and exploitation. The pseudo-code of CoPSO is shown in Algorithm 1.
Algorithm 1: Pseudo-code of CoPSO
Input: The value of Lmax, the value of μ, and the maximum number of iterations Imax
Output: The final solution gbest and its fitness value f(gbest)
1 Initialize gbest;
2 for each particle i do
3 Generate xi and vi randomly;
4 pbesti = xi;
5 if f(xi) is better than f(gbest) then
6 gbest = xi;
7 end
8 end
9 while Imax is not met do
10 for each particle i do
11 Update vi according to Eq (1);
12 Update xi according to Eq (2);
13 if f(xi) is better than f(pbesti) then
14 pbesti = xi;
15 end
16 if f(xi) is better than f(gbest) then
17 gbest = xi;
18 end
19 if pbesti hasn’t changed in μ continuous iterations then
20 li+ = 1;
21 end
22 if li == Lmax then
23 Regenerate xi and vi randomly;
24 li = 0;
25 end
26 end
27 end
Problem-independent simulation-optimization approach
The architecture of the proposed MAS-CoPSO approach is shown in Fig 3, where the MAS module and the CoPSO module are the main components. For the optimization problem brought up by the user (e.g., decision making, scheduling, route planning, etc.), a suitable simulation time step should be chosen first according to specific backgrounds. In every time step, MAS is supposed to simulate the problem scene, during which the agents’ behaviors might change the values of related variables and parameters. Simulation statistics are then put into the optimization module, whose objective can be maximizing the efficiency or minimizing the cost of certain activities. CoPSO algorithm is supposed to optimize the objective function while satisfying all the constraints. The optimization result will then play as the input of MAS, instructing the simulation process in the next time step. This iteration will go on until meeting the stop criteria set by the user. The outcome of MAS-CoPSO is the final solution to the problem. For complex systems, the solution may not be strictly optimal (actually it is intractable to find the strict global optimum for these problems). However, MAS-CoPSO can compromise between the solution quality and computational cost, making it possible to get a satisfying result in a reasonable period of time.
Successful simulation of the problem scene relies on a proper setup of the environment (which reflects the background of the problem) and an appropriate design of all the agents. In MAS, agents have autonomous behaviors and complex interactions, and the relationship between them can be various (e.g., cooperation, competition, control, etc.). Value changes of related variables and parameters link simulation and optimization together. For CoPSO, simulation results will influence the calculation of the objective function’s fitness value and the expression of constraints. Optimization results of CoPSO will function as the instructions for simulation in the next time step, i.e., how to adjust agents’ behaviors and tune the values of variables, parameters, and environmental settings. These two modules compose a cohesive framework through consistent information exchange and close integration.
Experiments and analysis of the optimization algorithm
To substantiate the validity of the proposed CoPSO algorithm, this study conducts numerical comparison experiments and computational complexity analysis. Exhaustive descriptions and discussions are given in this section.
Comparison experiments
The performance of CoPSO is compared with a set of classic algorithms:
- SPSO [3]: the standard formation of PSO, introduced in detail in the previous section.
- RPSO [5]: PSO with a random inertia weight factor designed for tracking dynamic systems.
- CPSO [6]: PSO with a new constriction factor related to the original acceleration coefficients.
- IPSO [8]: PSO with an increasing population size, i.e., whenever the algorithm can not find a satisfactory solution, add a new particle to the population.
Benchmarks and parameter setting.
Five well-known benchmark functions commonly used in literature [21] are applied to evaluate the performance of CoPSO, both in terms of solution quality and convergence rate. They are non-linear minimization problems that present different difficulties to the optimizers. Table 1 lists the functions, the problem dimension n, the global minimum fitness value fmin, and the search space ranges (which are also the initial ranges in this research).
The main parameters are set based on the recommendation in [22]. In CPSO, the maximal velocity Vmax and the minimal velocity Vmin are 100,000 and -100,000, respectively (since it is believed that Vmax and Vmin are not even needed in CPSO [6]). For other algorithms, Vmax is 10% of the search space’s upper bound while Vmin is 10% of the lower one. The acceleration constants c1 and c2 are both 2.05 for CPSO. For SPSO, RPSO, and IPSO, the figure is 2.0. For CoPSO, Lmax is 3, the initial value of μ is 5 (added by 1 for each level up), and the gap between the acceleration constants of each level is 0.2 (c1 and c2 are the same with a maximum value of 2.0). A fixed inertia weight value of 0.6 is adopted for SPSO, IPSO, and CoPSO. The maximum iteration number is 1000 for all algorithms, meaning that the maximum particle number is 1000 in IPSO. In other cases, the population size is 80. A total of 50 repeating runs are conducted for each experiment. In addition, as recommended by Yang et al. [23], Wilcoxon rank sum tests are conducted for the statistical analysis between CoPSO and the other algorithms with a significance level of 0.05.
Results and discussion.
Table 2 summarizes all the experimental results of 50 independent runs. Wilcoxon rank sum tests further verify the significance of these numerical results. The numbers in bold represent the comparatively best values. In general, CoPSO greatly outperforms SPSO, RPSO, CPSO, and IPSO on most metrics. For the exceptions, it still shows equivalent performances to the comparison algorithms. These results indicate that the multi-level structure and the competition mechanism can improve the solution quality of CoPSO.
It is believed that CoPSO can preserve excellent exploration and exploitation abilities while keeping a good balance between these two. In other words, CoPSO is able to jump out of local optima and fully exploit the promising areas at a fast rate. To substantiate this hypothesis, swarm diversity is used to measure an algorithm’s exploration ability. Supposing the size of swarm S is N, the problem dimension is D, and xi represents the position of particle i. As recommended by Yang et al. [23], the diversity of swarm S is computed as follows:
where D(S) is the swarm diversity of S and
is the average position of the swarm. Further, the converging speed can reflect whether an algorithm can compromise between exploration and exploitation properly. During the iterations, average swarm diversities and average global best fitness values (of 50 runs, in logarithmic form) for each algorithm on five benchmarks are recorded in Fig 4.
Note that the vertical axes are in logarithmic form. (a) The results on f1. (b) The results on f2. (c) The results on f3. (d) The results on f4. (e) The results on f5.
Generally, CoPSO converges faster with better solutions than all other algorithms in most cases. As for swarm diversity, CoPSO holds a similar tendency on all benchmarks, i.e., decreasing rapidly in the early stage of the evolution and staying at a relatively high value (or even the highest) in the latter stage. This can be explained by the structure of CoPSO. Particles of different levels are in different regions of the search space. Those of lower levels are more likely to be far away from the global best position while those of higher levels can be very close to it. CoPSO treats different levels differently by enhancing exploration for lower levels and promoting exploitation for higher ones. Thus, the algorithm has more potential to locate promising areas in the early stage and then converges to them very quickly, causing the swarm diversity to decline. However, particles may not upgrade because they find a global optimum, but because they are trapped in local optima. In this case, the competition mechanism will force the swarm to escape from local areas, which increases the diversity simultaneously.
Despite the overall similarity, the behaviors of CoPSO differ slightly on different kinds of test functions. f1 is a simple unimodal function with only one global minimum, making it possible to achieve fast convergence. Although all algorithms converge exponentially to the optima, CoPSO greatly surpasses the others due to better exploitation of the promising areas at the beginning and high-intensity exploration in the latter stage. For f2, CoPSO has a competitive performance compared to RPSO and even shows better potential in the end as a result of higher swarm diversity. f3 is a noisy quartic function. f4 is a classic optimization problem with a global minimum inside a long, narrow, parabolic-shaped valley. f5 is a multimodal function where the number of local optima increases exponentially as the problem dimension increases. For all these complicated benchmarks, CoPSO always keeps an excellent balance between exploration and exploitation, which helps it to achieve the fastest converging rate or the shortest stagnation at local optima. The other algorithms lack the adjusting ability of CoPSO. For instance, IPSO puts too much emphasis on exploration (holding the highest swarm diversity on most benchmarks) while RPSO and CPSO are too focused on exploitation, which degenerates their overall performances.
Computational complexity analysis
According to Yang et al. [23], given a fixed number of fitness evaluations, the computational complexity of an evolutionary algorithm is generally calculated by analyzing the extra cost in each generation without considering the cost of function evaluations, which is problem-dependent. As CoPSO inherits the simple structure of SPSO, its computational complexity will be analyzed by comparison with SPSO.
Assuming the population size is N and the problem dimension is D, updating particles’ states takes O(N × D) time in each iteration for SPSO. By comparison, adding line 19 to line 25 in Algorithm 1 costs extra time for CoPSO. Specifically, particles tend to update pbesti frequently at the early stage of the iterations, so it only costs O(N) extra time to check li and update related parameters. As li increases, the algorithm may have to regenerate some new particles, which costs O(N × D) extra time in the worst case. However, this situation is only possible in certain generations. Besides, the value increase of μ also controls the time cost of this part. As for the space complexity, since CoPSO needs to store the information of li and the time particles stay at the same pbesti, it requires O(N) extra memory than SPSO (which takes O(N × D) space to save particles’ positions and velocities). To conclude, with proper parameter setting, the computational complexity of CoPSO can be controlled within an acceptable range.
Table 3 lists the average computing time for a single run (in 50 repeating experiments) with the same parameters as the numerical experiments. SPSO, RPSO, and CPSO take nearly the same time on different functions due to their similar algorithmic structures. However, IPSO has a growing population, dramatically increasing fitness evaluations at the latter stage of the iterations. Consequently, the computational cost of IPSO is the highest. Based on the complexity analysis above, CoPSO takes acceptable extra time in each iteration compared to SPSO, verified by the experimental result that CoPSO costs slightly higher than SPSO on the benchmarks. In conclusion, CoPSO surpasses the other algorithms in solution quality while remaining computationally efficacious.
Case study of the simulation-optimization approach
To demonstrate that MAS-CoPSO can be used in practical applications, a case study of a real marine oil spill accident is carried out. This section presents the accident background, the implementation of MAS-CoPSO in a response planning system, and the mathematical formulation of the optimization problem. Further, an in-depth analysis of the system’s performance is given.
Background
Frequent marine oil spills have posed severe tests to the environment. Reducing the loss caused by oil spills has become an important issue. The capability of the MAS-CoPSO approach is tested through a case study of the Sanchi oil spill accident, which is caused by a collision between Panamanian oil tanker Sanchi and Chinese bulk carrier Changfeng Crystal. The accident happened on January 6, 2018, and the collision site is about 160 nautical miles east of the Yangtze River Estuary. Fig 5 illustrates the approximate location where the accident happened and the areas around it. Tons of condensate oil carried by Sanchi leaked into the ocean, causing significant economic loss and environmental pollution. Researchers have pointed out that oil tankers are transporting around 90% of all the oil around the world [24]. Furthermore, ship collision is the primary cause of many catastrophic marine oil spills, which suggests the representativeness of the case study.
Map services and data available from the U.S. Geological Survey.
Zhong and You [25] have concluded that after spill accidents, oil leaked into the ocean forms scattered slicks around the spill site and undergoes various physical and chemical changes due to environmental forces. The most dominant processes are spreading, evaporation, emulsification, and dispersion. These processes will dramatically change the volume, area, thickness, and viscosity of those slicks, which will affect the efficiency of oil spill response operations as a consequence.
There are some equipped response teams berthed at docks near the spill site. They will conduct a wide range of actions to clean up the ocean, among which booms, skimmers, chemical dispersants, and in situ burning are the most frequently used methods. Booms are generally the first equipment deployed after an oil spill and are often used to protect shorelines, divert oil to certain areas, or concentrate oil for further recovery and burning. Skimmers are adopted to recover oil or oil-water mixtures from the water surface. The characters of the slick will affect the working efficiency of skimmers. Chemical dispersants can reduce the oil-water interfacial tension, which should be initialized as soon as possible. In situ burning indicates controlled burning of the oil at or near the site. In this case, it is assumed that the spraying of dispersants and the concentration of leaked oil have been completed quickly after the accident. Besides, in situ burning will emit toxic substances and waste the leaked oil. Thus, the primary response operation is recovery.
This study focuses on supporting response planning after the accident. A simulation-optimization system is built based on the proposed MAS-CoPSO method. The overall objective is to clean up the polluted sea area in the most efficient way, i.e., as soon as possible. The system is supposed to realize interactive spill simulation and response optimization while considering the influences of the oceanic environment. It should also give the most reasonable allocations for available resources and the schedule of response teams’ actions.
Implementation
The implementation of MAS-CoPSO in the system is described in Fig 6, which is an instantiation of Fig 3. The MAS module and the CoPSO module work collectively, composing a dynamic system for marine oil spill response planning. The CoPSO algorithm is encapsulated as a function to be called before each time step. Meanwhile, agents of the simulation module consistently update their behaviors and interact with each other, realizing the oil spill fate modeling and response operation modeling. The agents must follow specific rules in reality. For instance, the property changes of oil slicks are calculated by weathering models considering spreading, evaporation, emulsification, and dispersion. Skimmers and ships of response teams should obey the rules for oil recovery. The performance of a response team can be affected by the actions of other teams and the characteristics of oil slicks. For example, when a skimmer works on a slick to collect oil, the volume and area of the slick, the evaporated and dispersed oil rates, the viscosity and the water content of the oil, will all be affected. If more than one team are instructed to head for the same slick, they will collaborate.
To complete the global objective of minimizing response time while cleaning up the spill site to a large extent (the stop criteria, only a small portion of leaked oil remained in the ocean), CoPSO is utilized to maximize the decreased volume of oil in each simulation time step. Meanwhile, MAS should simulate equipment location, response process, oil spill fate, and transportation. The system can be further updated to suit different purposes and requirements.
Mathematical formulation of the optimization problem
Since the weathering process dramatically affects oil properties, it is necessary to simulate it as accurately as possible. Furthermore, the objective function and constraints of the optimization problem can be derived from the model.
Simulation of the weathering process.
First of all, the initial area of a slick, A0 (m2), can be calculated by the equation shown as follows [26]:
(3)
where ρω (g ⋅ cm−3) is the density of seawater, ρ0 (g ⋅ cm−3) is the density of oil, νω (0.801 × 10−6m2 ⋅ s−1 under 30°C) is the kinematic viscosity of seawater, V0 (m3) is the inital volume of the slick, g (m ⋅ s−1) is the acceleration of gravity, k2 and k3 are constants with values of 1.21 and 1.53, respectively.
The area changing rate of the slick due to spreading can be modeled as follows [27]:
(4)
where A (m2) is the current surface area of the slick, t (s) refers to the time since the accident happened, V (m3) is the current volume of the slick, K1 (s−1) is a dominant physicochemical parameter of the oil with a default value of 150.
Evaporation is assumed as the most influential environmental effect, causing a loss of up to 20-50% of the spill. The rate that oil evaporates from the sea surface is modeled by the following equation [28]:
(5)
where FE (%) is the current volume fraction of the oil that has been evaporated with an initial value of FE(t=0) = 0, TK (K) is the oil temperature, Aev and Bev are empirical constants with fixed values of 6.3 and 10.3, respectively. Kev (m ⋅ s−1) is the mass transfer coefficient for evaporation which can be calculated as follows [29]:
where vwind (m ⋅ s−1) is the wind speed. TO and TG are the initial boiling point and the gradient of the oil distillation curve, respectively. Their values can be calculated through functions of the oil API (American Petroleum Institute) degree as follows:
Emulsification can also influence the slick characteristics dramatically, especially for the oil viscosity. The initial oil viscosity μ0 (%) can be calculated using the following equation [29]:
(6)
where AC (%) is the asphaltene content of the parent oil, which is a constant. As the viscosity changes over time, its changing rate is given by [30]:
(7)
where μ (%) is the current viscosity, C3 is a constant for the final water content, C4 is an oil-dependent constant, YW (%) is the fractional water content in the emulsion, it changes over time, which can be computed with the following equation [27]:
where Kem is an empirical constant between 1 × 10−6 and 2 × 10−6, YW has an initial value of YW(t = 0) = 0.
Finally, dispersion is also a vital factor in the loss of oil. The volume of oil naturally dispersed changes over time, which can be modeled by [27]:
(8)
where VD (m3) is the total dispersed oil volume until now with an initial value of VD(t=0) = 0 and ζt (cP) is the oil-water interfacial tension.
Objective function and constraints.
In the response process, the optimization problem is formulated on the basis of each simulation time step. Suppose that the ith response team is recovering oil on the kth slick during time step Tstep, define ORRi (m3 ⋅ hr−1) as the oil recovery rate of team i (i.e., the amount of oil that the team can recover per hour) and STk (m) as the current oil thickness of slick k. ORRi is determined by the properties of the team’s skimmer and STk:
(9)
where α and β are empirical coefficients obtained from experimental tests, reflecting the features of the skimmer. STk can be calculated by:
(10)
where Vk (m3) and Ak (m2) are the current oil volume and the area of the kth slick, respectively. Ak can be calculated according to Eqs (3) and (4), Vk can be computed as follows:
(11)
where Vk0 (m3) is the remaining oil volume of this slick at the end of the last time step (i.e., the initial oil volume of the current time step), EVk (m3) is the volume loss caused by evaporation, DVk (m3) is the volume loss caused by dispersion, and SVi (m3) is the volume recovered by team i. Notice that SVi may not exist, which indicates that no response team is working on slick k right now. The system can get the value of EVk, DVk, and SVi according to Eqs (5), (8) and (9), respectively. During the simulation, for simplification, a short Tstep should be chosen and the system will use the value of Vk0 for the calculations regarding Vk.
The aim of the system is to clean up the polluted sea area as soon as possible. Therefore, based on the proposed method, given a fixed simulation time step, the oil volume loss will be maximized during each step, which is composed of the evaporation loss, the dispersion loss and the recovered volume. Suppose there are m oil slicks and n response teams during Tstep, Vstep represents the total oil volume loss in this period. FEk and μk refer to the evaporated oil fraction and the oil viscosity of slick k, respectively. The objective function and constraints are given by:
(12)
s.t.
(13)
(14)
(15)
(16)
(17)
(18)
Performance of the system
To examine the efficiency of MAS-CoPSO, its performance is compared with the shortest distance selection approach (SDS). Ye et al. [31] have suggested that SDS is a simple strategy commonly used in marine oil spill emergency response, which indicates a process that allows a response team to choose the nearest oil slick as the target for oil recovery. After meeting the cleaning requirements, they will choose the second nearest oil slicks to continue.
Therefore, when conducting SDS, the distances between ships and slicks are the judgment criteria for response planning. When initializing the simulation process, each response team (corresponding to a certain response agent) is assigned with the closest oil slick (corresponding to a certain spill agent) to clean. Before each time step, instead of calling the CoPSO module, SDS will check whether the slick is cleaned up, if the oil volume is below the set threshold, the ship will head to the closest one in the remaining slicks, otherwise, it will stay at the current location until finishing the cleaning task. Compared to MAS-CoPSO, SDS is more straightforward, which explains why it is commonly used in real accidents. However, SDS makes plans without considering the oceanic forces and the interactions between agents, which according to previous discussions, can dramatically influence the overall efficiency. Thus, comparative experiments are conducted for further analysis.
Simulation settings.
Chen et al. [32] introduced the process and consequences of the accident in detail. Based on their description, it is assumed that the spilled oil was split into ten slicks within the East Sea with a total volume of 113,000 tons. Table 4 lists the oil volumes and locations (described by the distance and direction from the approximate collision site) of these slicks. The ship-mounted skimmers belonging to five different response teams are the only available nearby cleanup means to be applied. Assume that all the teams have enough storage space for recovered oil. As the allocation process needs specific transportation time, the speed of ships is set at 30 km ⋅ h−1. Skimmers of different teams have different recovering efficiencies, which are affected by α and β in Eq (9). Table 5 lists the values of the empirical coefficients for different teams. Fig 7 illustrates the simulation process. The red cross represents the collision site, the black spills represent the oil slicks (the size reflects the oil volume of the corresponding slick), and the ships represent the response teams carrying the recovery task. The position and size of a slick will change over time. A ship can either be heading to a slick or working on it.
As mentioned above, spreading, evaporation, emulsification, and dispersion are the main oceanic effects. Table 6 lists the inputs for the modeling of these processes, including environmental parameters and the characteristics of the leaked oil (condensate oil). As [26] suggested, the drifting speed of oil slicks is about 2.5-4.5% of the wind speed, so it is set at 0.03 m ⋅ s−1.
The model is built in AnyLogic®. 50 repeating runs are carried out with 50 iterations per time step (Tstep is one hour). Once the volume of the remaining oil is reduced to 10% of the original figure, the simulation will stop.
Results and discussion.
Table 7 presents the simulation results of MAS-CoPSO and SDS with the same stop criteria and time step. The operation time for achieving an oil recovery rate of 90% is 129.82 hours based on the optimal outcomes given by the MAS-CoPSO system, which is only 74.2% of the figure for SDS. Shorter response time means higher operation efficiency, less volume of wasted oil, and slighter damage to the environment. The results also indicate that total recovered oil holds a proportion of 65.54% in MAS-CoPSO, which is approximately 1.2 times higher than the number for SDS.
Fig 8 illustrates the variations of remaining oil volume for each slick in MAS-CoPSO and SDS scenarios during the entire simulation process. The optimization outcomes of MAS-CoPSO are highly related to the thickness of the oil slicks, and the system intends to keep a balanced volume level for each slick by optimizing the time for allocation and the response efficiency (see Fig 8a). However, oil thickness does not affect the decision-making of the SDS model, so it will not try to reduce the volume steadily. As a result of different strategies, the curves of MAS-CoPSO are much smoother than those of SDS. In addition, the SDS method might ignore slicks too far from the response teams (e.g., Slick 1), leading to long-time exposure of the oil, emitting more harmful substances (see Fig 8b).
(a) The results of MAS-CoPSO. (b) The results of SDS.
Fig 9 depicts the variations of accumulated recovered volume for each team in two approaches. MAS-CoPSO tries to find optimal schedules for all teams in every time step, so the curves climbs stably (see Fig 9a). Meanwhile, SDS always cleans up one slick before moving to another, leading to discontinuous curves, which also hinders the increase of the overall efficiency (see Fig 9b).
(a) The results of MAS-CoPSO. (b) The results of SDS.
The above results show that the response planning system based on MAS-CoPSO outperforms the commonly used SDS strategy in terms of time consumption, oil recovery rate, and environmental impact. Therefore, oil spill response teams can schedule their operations according to the simulation-optimization outcomes of this system rather than relying on the inefficient SDS approach. Moreover, complex problems and high-intensity interactions can enhance the advantages of MAS-CoPSO, indicating the system’s ability to function in scenarios with higher oil volumes and more response teams. Even though the case study is simplified and focuses on the oil recovery process, the system has the potential to comprehensively support multiple cleanup techniques concerning booms, chemical dispersants, and in situ burning. Besides, modeling of the weathering process can be easily enriched by considering more complicated oceanic forces such as dissolution, photo-oxidation, sedimentation, and biodegradation [25]. Hydrodynamic simulation of oil spill trajectories can also be considered. In addition, the application range of the MAS-CoPSO-based system can be further enlarged by handling uncertainties and risk assessments. With these extensions, oil spill response teams can get a practical simulation-optimization tool to support their planning process.
Conclusion
An improved particle swarm optimization algorithm called CoPSO is proposed in this paper. Inspired by the features of bird flocks, a multi-level design and a competition mechanism are introduced into CoPSO to balance the exploration and exploitation more effectively, increase solution accuracy, achieve fast convergence, and avoid stagnation at local optima. The performance of CoPSO is tested on a series of well-known benchmark functions. Comparisons are made between CoPSO, SPSO, and three classic variants of it. Experimental results have verified its capability.
This paper also proposes a problem-independent simulation-optimization approach called MAS-CoPSO to combine CoPSO with MAS. MAS-CoPSO can achieve a dynamic simulation of various complex systems and give satisfying solutions to related optimization problems. In MAS-CoPSO, the relationship between the simulation module and the optimization module is highly cohesive. The effectiveness of MAS-CoPSO is substantiated through a response planning system for the Sanchi oil spill accident. The system simulates the accident scene with consideration of the oil slicks’ physicochemical evolution. Case study results show that the system can optimize response device allocation, operation scheduling, and time consumption, causing slighter damage to the environment.
MAS-CoPSO can also be applied in other scenarios to solve different optimization problems depending on the users’ requirements. For future studies, we will consider the combination with other methods, such as fuzzy control and neural network.
References
- 1.
Olafsson S, Kim J. Simulation optimization. In: Proceedings of the winter simulation conference. vol. 1. IEEE; 2002. p. 79–84.
- 2.
Macal CM, North MJ. Tutorial on agent-based modeling and simulation. In: Proceedings of the Winter Simulation Conference, 2005. IEEE; 2005. p. 14–pp.
- 3.
Kennedy J, Eberhart R. Particle swarm optimization. In: Proceedings of ICNN’95-international conference on neural networks. vol. 4. IEEE; 1995. p. 1942–1948.
- 4. Bonyadi MR, Michalewicz Z. Particle swarm optimization for single objective continuous space problems: a review. Evolutionary computation. 2017;25(1):1–54. pmid:26953883
- 5.
Eberhart RC, Shi Y. Tracking and optimizing dynamic systems with particle swarms. In: Proceedings of the 2001 congress on evolutionary computation (IEEE Cat. No. 01TH8546). vol. 1. IEEE; 2001. p. 94–100.
- 6. Clerc M, Kennedy J. The particle swarm-explosion, stability, and convergence in a multidimensional complex space. IEEE transactions on Evolutionary Computation. 2002;6(1):58–73.
- 7. Neshat M. FAIPSO: fuzzy adaptive informed particle swarm optimization. Neural Computing and Applications. 2013;23(1):95–116.
- 8. De Oca MAM, Stutzle T, Van den Enden K, Dorigo M. Incremental social learning in particle swarms. IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics). 2010;41(2):368–384. pmid:20875976
- 9. Neshat M, Sargolzaei M, Masoumi A, Najaran A. A new kind of PSO: predator particle swarm optimization. International Journal on Smart Sensing and Intelligent Systems. 2012;5(2).
- 10. Li D, Guo W, Lerch A, Li Y, Wang L, Wu Q. An adaptive particle swarm optimizer with decoupled exploration and exploitation for large scale optimization. Swarm and Evolutionary Computation. 2021;60:100789.
- 11. Rizk Y, Awad M, Tunstel EW. Decision making in multiagent systems: A survey. IEEE Transactions on Cognitive and Developmental Systems. 2018;10(3):514–529.
- 12.
Bo Y, Cheng W, Hua H, Lijun L. A multi-agent and PSO based simulation for human behavior in emergency evacuation. In: 2007 International Conference on Computational Intelligence and Security (CIS 2007). IEEE; 2007. p. 296–300.
- 13. Ali ZA, Han Z, Masood RJ. Collective Motion and Self-Organization of a Swarm of UAVs: A Cluster-Based Architecture. Sensors. 2021;21(11):3820. pmid:34073061
- 14. Kanaga EGM, Valarmathi M. Multi-agent based patient scheduling using particle swarm optimization. Procedia Engineering. 2012;30:386–393.
- 15.
Yu C, Wang G, Liu Y. A multi-agent based architecture for web service selection in E-business. In: 2011 IEEE 8th International Conference on e-Business Engineering. IEEE; 2011. p. 245–250.
- 16.
Thiel D, Hovelaque V, Pham DN. A multi-agent model for optimizing supermarkets location in emerging countries. In: 2012 IEEE 13th International Symposium on Computational Intelligence and Informatics (CINTI). IEEE; 2012. p. 395–399.
- 17. Malik MM, Kazmi SAA, Asim HW, Ahmed AB, Shin DR. An intelligent multi-stage optimization approach for community based micro-grid within multi-microgrid paradigm. IEEE Access. 2020;8:177228–177244.
- 18. Wang L, Wang Z, Yang R. Intelligent multiagent control system for energy and comfort management in smart and sustainable buildings. IEEE transactions on smart grid. 2012;3(2):605–617.
- 19.
Angeline PJ. Evolutionary optimization versus particle swarm optimization: Philosophy and performance differences. In: International Conference on Evolutionary Programming. Springer; 1998. p. 601–610.
- 20. Poesel A, Fries AC, Miller L, Gibbs HL, Soha JA, Nelson DA. High levels of gene flow among song dialect populations of the Puget Sound white-crowned sparrow. Ethology. 2017;123(9):581–592.
- 21.
Vesterstrom J, Thomsen R. A comparative study of differential evolution, particle swarm optimization, and evolutionary algorithms on numerical benchmark problems. In: Proceedings of the 2004 congress on evolutionary computation (IEEE Cat. No. 04TH8753). vol. 2. IEEE; 2004. p. 1980–1987.
- 22.
Eberhart RC, Shi Y. Comparing inertia weights and constriction factors in particle swarm optimization. In: Proceedings of the 2000 congress on evolutionary computation. CEC00 (Cat. No. 00TH8512). vol. 1. IEEE; 2000. p. 84–88.
- 23. Yang Q, Chen WN, Da Deng J, Li Y, Gu T, Zhang J. A level-based learning swarm optimizer for large-scale optimization. IEEE Transactions on Evolutionary Computation. 2017;22(4):578–594.
- 24. Zhang HY, Ji Q, Fan Y. What drives the formation of global oil trade patterns? Energy Economics. 2015;49:639–648.
- 25. Zhong Z, You F. Oil spill response planning with consideration of physicochemical evolution of the oil slick: A multiobjective optimization approach. Computers & Chemical Engineering. 2011;35(8):1614–1630.
- 26.
Fay JA. The spread of oil slicks on a calm sea. In: Oil on the Sea. Springer; 1969. p. 53–63.
- 27.
Mackay D, Buist I, Mascarenhas R, Paterson S. Oil spill processes and models: Environment Canada Manuscript Report No 8. EE-8, Ottawa, Ontario. 1980;.
- 28. Mackay D, Matsugu RS. Evaporation rates of liquid hydrocarbon spills on land and water. The Canadian Journal of Chemical Engineering. 1973;51(4):434–439.
- 29. Buchanan I, Hurford N. Methods for predicting the physical changes in oil spilt at sea. Oil and Chemical Pollution. 1988;4(4):311–328.
- 30. Mooney M. The viscosity of a concentrated suspension of spherical particles. Journal of colloid science. 1951;6(2):162–170.
- 31. Ye X, Chen B, Li P, Jing L, Zeng G. A simulation-based multi-agent particle swarm optimization approach for supporting dynamic decision making in marine oil spill responses. Ocean & Coastal Management. 2019;172:128–136.
- 32. Chen J, Di Z, Shi J, Shu Y, Wan Z, Song L, et al. Marine oil spill pollution causes and governance: A case study of Sanchi tanker collision and explosion. Journal of Cleaner Production. 2020;273:122978.