Combining an Extended SMAA-2 Method with Integer Linear Programming for Task Assignment of Multi-UCAV under Multiple Uncertainties

: Uncertainty should be taken into account when establishing multiobjective task assignment models for multiple unmanned combat aerial vehicles (UCAVs) due to errors in the target information acquired by sensors, implicit preferences of the commander for operational objectives, and partially known weights of sensors. In this paper, we extend the stochastic multicriteria acceptability analysis-2 (SMAA-2) method and combine it with integer linear programming to achieve multiobjective task assignment for multi-UCAV under multiple uncertainties. We ﬁrst represent the uncertain target information as normal distribution interval numbers so that the values of criteria (operational objectives) concerned can be computed based on the weighted arithmetic averaging operator. Thus, we obtain multiple criteria value matrices for each UCAV. Then, we propose a novel aggregation method to generate the ﬁnal criteria value matrix based on which the holistic acceptability indices are computed by the extended SMAA-2 method. On this basis, we convert the task assignment model with uncertain parameters into an integer linear programming model without uncertainty so as to implement task assignment using the integer linear programming method. Finally, we conduct a case study and demonstrate the feasibility of the proposed method in solving the multiobjective task assignment problem multi-UCAV under multiple uncertainties.


Introduction
Unmanned combat aerial vehicles (UCAVs), compared with airborne weapons and airborne platforms (e.g., cruise missiles, manned fighters, and bombers), have demonstrated advantages in tactical flexibility and multitasking ability due to having zero casualties, excellent stealth performance, high flight height, long endurance, multiple hard points for missiles, and a relatively low life-cycle cost [1].Benefiting from the rapid development of control theory and intelligent technology, UCAVs are increasingly likely to become the frontline fighter aircrafts in modern or future wars [2].As one of the cores of mission planning, cooperative task assignment for UCAVs has been an emerging issue in terms of its application in attacking ground targets.The purpose of task assignment is to match tasks to available UCAVs without conflicts under the constraints of the battlefield environment, operational capabilities of the UCAVs and weapons on board, and target types in order to obtain the maximum operational benefit, such as maximizing damage to hostile targets and minimizing the user's own consumption.With regard to research on task assignment, most of the existing methods assume that the parameters related to mission scenarios are deterministic.These methods are mainly divided into two categories: methods based on exact algorithms and methods based on heuristic (metaheuristic) algorithms [3].Exact algorithms are mostly based on the mathematical properties of task assignment models and attempt to obtain the optimal solution for given objectives.They include the mixed integer linear programming (MILP) formulation [4][5][6], the branch and bound tree search algorithm [7,8], and the dynamic programming algorithm [9] among others.Most research to date on the application of heuristic (metaheuristic) algorithms have either extracted some specific search rules based on the properties of the problem to obtain optimal or suboptimal solutions rapidly or introduced some local search mechanism to the basic algorithm framework to improve the solution quality, such as tabu search algorithms [10], auction algorithms [11], genetic algorithms [12][13][14][15], ant colony algorithms [16], and particle swarm optimization [17,18].The algorithms mentioned above-whether belonging to exact algorithms or heuristic (metaheuristic) algorithms-have demonstrated the ability to provide optimal or suboptimal solutions for task assignment problems of UAVs or UCAVs in various mission scenarios but may become difficult to apply when there are uncertain parameters related to mission scenarios.
From what we know, uncertainties in the task assignment model for multi-UCAV exist in at least three aspects: (1) the target information obtained by various sensors, (2) the weights of operational objectives, and (3) the weights of these sensors.Considering that sensors may be affected by the dynamic battlefield environment or the antireconnaissance measures taken by hostile forces in the reconnaissance process, it is reasonable to believe that there are errors in the target information obtained by various sensors, indicating that the target information is uncertain.When there are multiple operational objectives to consider, the commander needs to weigh the importance of these objectives but can usually only provide the relative importance of operational objectives instead of the exact weights of these objectives.For example, the commander may, at most, express a clear preference for operational benefit over operational cost but would not provide the exact weights of the two operational objectives.Therefore, the weights of operational objectives are also thought to be uncertain.In addition, the uncertainty contained in the weights of sensors also needs to be considered.In this paper, the weight of a sensor is obtained by integrating its subjective weight and objective weight in a certain way.The subjective weight is usually determined in advance by experts according to the performance parameters of the sensor, while the objective weight is heavily dependent on the accuracy and credibility of the target information provided by the sensor.As the target information is uncertain, the objective weight of the sensor is unknown, which indicates that the integrated weight of the sensor is partially known.
In order to reflect the combat scenarios as accurately as possible, the above three aspects of uncertainties need to be taken into account when formulating the task assignment problem.However, in this case, the problem would become extremely complicated and can neither be solved by exact algorithms nor by heuristic (metaheuristic) algorithms.Therefore, in order to address the challenge, it is necessary and worthwhile to propose new methods.In this way, the optimal or suboptimal task assignment scheme under severe uncertainties can be found.

Related Work
To our knowledge, Mataric et al. first studied task assignment for multirobots in uncertain environments in 2003 [19].The uncertainty in their research mainly included sensor noise and the randomness of the moving direction of robots.To simplify, they used a random number from a normal distribution to represent the sensor noise and a finite but small probability that robots would move in a random direction instead of the intended direction at each time step to model the randomness of the movement.Since then, various methods have been proposed to solve task assignment problems under different types of uncertainty.
In general, uncertainty is mathematically described by one or a combination of the following theories: probability theory, fuzzy set theory, rough set theory, and interval set theory [20].In the existing literature on task assignment under uncertainty, as far as we know, probability theory has been the most widely used.The simplest way to model uncertainty based on probability theory is to assume that an uncertain parameter follows a specific distribution.Based on this assumption, other methods can be applied to task assignment.Mataric et al. assumed that sensor noise followed a normal distribution and that each robot had a small probability of moving in other directions than the intended direction.The authors then adopted the auction algorithm to assign tasks to robots [19].Kim et al. used two certain probability values to measure the likelihood that a UAV would sense targets within its sensing range and the likelihood that a UAV would fail during the mission.The authors then proposed a distributed task assignment method based on resource welfare [21].Jia et al. assumed that the flight velocity of a UAV was a stochastic variable that corresponded to a normal distribution with alterable variance according to the actual flight distance between each pair of targets.The authors then solved the task assignment problem in the framework of a two-stage stochastic programming model using a genetic algorithm [3].A similar way of dealing with uncertain parameters can be found in Reference [22], where the task duration time was considered to be subject to a uniform distribution.
In multistage task assignment problems, the Markov decision process (MDP) or partially observable Markov decision process (POMDP) is applied based on the assumption that the uncertain parameters that are closely related to the state transition are subject to either predefined distributions or predefined probability values, as can be seen in References [23][24][25][26][27].For example, in Reference [24], it was assumed that each target could move to any of its neighboring cells at each time step with the same probability and that a UAV could detect a target with a predefined probability.Based on the above assumptions, the possible states of UAVs and targets, the actions that UAVs could take, and the corresponding rewards of UAVs could be obtained in turn.Then, the expected cumulative reward earned during finite time steps could be computed.Thus, the task assignment scheme could be acquired by maximizing the expected cumulative reward.
In addition to probability theory, interval set theory has also been applied to task assignment problems.Liu et al. proposed an interval Hungarian algorithm for task assignment of UAVs when the utility of a UAV performing a task is an interval number.The proposed algorithm was an extension of the classic Hungarian algorithm and was able to compute the maximum interval of deviation for each entry in the assignment matrix, which would retain the same optimal assignment [28,29].Hooshangi et al. used the interval VIKOR method to sort tasks, the characteristics of which were specified as interval numbers.They then applied an auction algorithm to assign tasks to agents [20].The task assignment problem in this research was, to some extent, regarded as a decision problem with multiple uncertain criteria, and the simulation results demonstrated that the proposed method was able to solve task assignment problems with more uncertainty than the other methods mentioned above.Despite this, when the criteria values and the corresponding weights are both uncertain, the interval VIKOR method would become unavailable.
The methods introduced above have been proven to be effective in solving the task assignment problem where uncertainty exists either in the action of agents (UAVs) or in task rewards.However, when multiple uncertainties are simultaneously contained in the task assignment problem, as described in Section 1, none of these methods can be applied.Therefore, we need to put forward new methods to solve the problem.As a multicriteria decision support method, the stochastic multicriteria acceptability analysis (SMAA) is able to solve problems with unavailable criteria weight information and uncertain criteria values [30].Hu et al. investigated the application of the SMAA-2 method, a variant of the SMAA method, in task assignment of multi-UAV when both the target information and the weights of operational objectives were uncertain.The proposed method demonstrated the availability of the SMAA-2 method in solving the task assignment problem that contained uncertainties from two aspects.However, the method was limited in its application [31] because the premise of applying the SMAA-2 method was that there could only be one criteria value matrix: the evaluation of all targets in terms of all concerned operational objectives.Note that in this paper, the criteria in the SMAA-2 method are equivalent to the operational objectives.In the task assignment model described in Section 1, there are various sensors and each sensor has the ability to acquire battlefield information, so there are multiple criteria value matrices.In this situation, if the SMAA-2 method is to be applied, these matrices would be first aggregated to generate a final criteria value matrix.Considering that the weights of sensors are partially known, the aggregation process cannot be directly performed.In summary, the original SMAA-2 method is difficult to use in such a situation, so it is in urgent need of improvement.
Inspired by References [30,31], we extend the original SMAA-2 method in this research by integrating a novel aggregation method so that multiple criteria value matrices can be aggregated to generate the final criteria value matrix based on which the SMAA-2 method can be applied.Then we combine the extended SMAA-2 method with integer linear programming (ILP) to solve task assignment problems for multi-UCAV where there are multiple sensors (or information sources) to provide uncertain information, describing multiple attributes of multiple targets.Furthermore, the weights of the criteria concerned are uncertain and the weight information of information sources is partially available.In subsequent sections, we will refer to the above combination method as ESMAA-2-ILP for convenience.The flow of the method is illustrated in Figure 1.
As can be seen from Figure 1, to achieve task assignment for multi-UCAV, we first formulate the task assignment problem with a variety of uncertain parameters to obtain the criteria formulas and then, based on the target information from K information sources, we construct K criteria value matrices for each UCAV through the above criteria formulas.For each UCAV, we apply the extended SMAA-2 method to compute the holistic acceptability index of attacking each target.Thus, we can obtain a holistic acceptability index matrix, which describes the holistic acceptability that each UCAV attacks each target.On this basis, according to the meaning of holistic acceptability indices, we convert the task assignment model established above into an integer linear programming model with an objective function to maximize the holistic acceptability index of the task assignment scheme.Finally, based on integer linear programming, we obtain the task assignment scheme.The main contributions of this research are summarized as follows:

•
We propose a novel method to aggregate multiple criteria value matrices from different information sources (ISs) under the condition that the criteria weight information is unavailable and the weight information about ISs is partially available.

•
We integrate the aggregation method into the original SMAA-2 method.Thus, we propose an extended SMAA-2 method.

•
We propose to convert the task assignment model for multi-UCAV with a variety of uncertain parameters into an integer linear programming model to maximize the holistic acceptability index of the task assignment scheme according to the meaning of holistic acceptability indices, and then apply integer linear programming to achieve task assignment.
The rest of this paper is organized as follows.We formulate the task assignment problem with a variety of uncertain parameters in Section 3. In Section 4, we briefly introduce the original SMAA-2 method.Section 5 gives details of the extended SMAA-2 method.In Section 6, we convert the task assignment model established in Section 3 into an integer linear programming model.We examine the proposed method through a case study in Section 7. Section 8 concludes the paper.

Problem Formulation
Here, we consider a combat scenario where M UCAVs under the command of one command center are to attack N geographically dispersed ground targets with different values.Each UCAV is loaded with a certain number of the same kind of air-to-ground missiles to destroy targets with different lethality probabilities.Prior to the strike, K sensors, including reconnaissance satellites, reconnaissance aircrafts, and reconnaissance UAVs, conduct multiple reconnaissances on targets and then transmit acquired target information to the command center.Based on the received target information, the command center computes K criteria value matrices for each UCAV and then aggregates them to generate a final criteria value matrix according to which commanders are able to achieve task assignment for multi-UCAV.The operational process described above is briefly shown in Figure 2. The command center receives the target information transmitted by three sensors (one reconnaissance satellite, one reconnaissance aircraft, and one reconnaissance UAV) and then assigns four targets to four UCAVs.

Problem Formulation
Here, we consider a combat scenario where M UCAVs under the command of one command center are to attack N geographically dispersed ground targets with different values.Each UCAV is loaded with a certain number of the same kind of air-to-ground missiles to destroy targets with different lethality probabilities.Prior to the strike, K sensors, including reconnaissance satellites, reconnaissance aircrafts, and reconnaissance UAVs, conduct multiple reconnaissances on targets and then transmit acquired target information to the command center.Based on the received target information, the command center computes K criteria value matrices for each UCAV and then aggregates them to generate a final criteria value matrix according to which commanders are able to achieve task assignment for multi-UCAV.The operational process described above is briefly shown in Figure 2. The command center receives the target information transmitted by three sensors (one reconnaissance satellite, one reconnaissance aircraft, and one reconnaissance UAV) and then assigns four targets to four UCAVs.We assume that each UCAV can only attack one target in one flight and that each target can only be attacked by one UCAV (for convenience, we assume M = N).Furthermore, we consider maximizing the damage to targets (i.e., operational benefit) while minimizing the flight consumption of UCAVs (i.e., operational cost) and then formulate the task assignment problem as follows: ( ) subject to We assume that each UCAV can only attack one target in one flight and that each target can only be attacked by one UCAV (for convenience, we assume M = N).Furthermore, we consider maximizing the damage to targets (i.e., operational benefit) while minimizing the flight consumption of UCAVs (i.e., operational cost) and then formulate the task assignment problem as follows: subject to In the above model, the objective function, Equation ( 1), aims to maximize the combat reward, a weighted sum of operational benefit O B and operational cost O C .The corresponding weights ω B and ω C are determined by commanders and should satisfy constraint in Equation (7).The constraints expressed in Equations ( 4)-( 6) guarantee a one-to-one relationship between the UCAVs and the targets.The constraints in Equations ( 8)-( 12) define the method to aggregate criteria value matrices.
Solutions for the above model are relatively easy to find when all parameters are deterministic.However, due to the complexity of actual combat, many parameters in the model are uncertain.
Considering errors in the target information obtained by sensors, the following parameters v ) are assumed to be uncertain.Moreover, in our research, criteria weight information is unavailable and sensor weight information is partially available, which means ω B , ω C and As can be seen, to solve the task assignment problem, we need to first aggregate K criteria value matrices under serve uncertainty, and then conduct the stochastic multicriteria acceptability analysis for each UCAV in terms of all targets.Therefore, we suggest extending the original SMAA-2 method by combining a novel aggregation method.To facilitate a better understanding of the extended SMAA-2 method, we give a brief introduction to the original SMAA-2 method in the next section.

The Stochastic Multicriteria Acceptability Analysis-2 (SMAA-2) Method
Stochastic multicriteria acceptability analysis (SMAA), first proposed by Lahdelma et al. in 1998 [30], is a multicriteria decision support method to solve discrete problems with unavailable criteria weight information and uncertain criteria values.Based on the exploration of the weight space, the method is able to describe the valuations that would make each alternative the preferred one.During this process, uncertain criteria values are always represented by probability distributions, thus the method can compute confidence factors to describe the reliability of the analysis.The SMAA-2 method, also proposed by Lahdelma et al., extends the original SMAA by analyzing the sets of weight vectors for any rank from best to worst for each decision alternative [32].It is introduced as follows: Here, we consider a multicriteria decision-making problem, including a set of N alternatives Set_alt = {A 1 , A 2 , • • • , A N }, from which one is to be chosen.All alternatives are evaluated in terms of n criteria Set_cri = {C 1 , C 2 , • • • , C n }.The method maps the ith alternative A i to a utility value U i by U i = U(ξ i , ω), where ] is a weight vector for each decision-maker to quantify their subjective preferences for different criteria, and U(•) is a utility aggregation function.In SMAA-2, any type of utility aggregation function is available [33].In this paper, we use the utility additive form and then obtain the following: where u ij (•) and u ij ξ ij are the utility value function and the utility value of alternative A i in terms of criterion C j , respectively.In SMAA-2, the criteria value ) is assumed to be a stochastic variable with an assumed or estimated joint probability density function f (ξ) in the criteria value space Ω.Similarly, the preferences of decision-makers are assumed to be unknown or partially known and can be represented by a weight distribution with density function g(ω) in the feasible weight space W, defined as follows: When the weights are totally unknown, g(ω) is thought to be a uniform distribution density function and is expressed as follows: where vol(W) = √ n/(n − 1)! is the volume of the feasible weight space, an (n − 1)-dimensional simplex.
Based on the above assumptions, stochastic criteria and weight distributions can be mapped into utility distributions U(ξ i , ω) through the utility function.The utility distributions are employed to compute the rank of each alternative, which is defined as an integer from the best rank (=1) to the worst rank (=N) by means of the rank function [34]: where ρ(ture) = 1 and ρ( f alse) = 0.Then, the favorable rank weight space W r i (ξ) is defined by restricting rank(ξ i , ω) = r, which means any weight ω ∈ W r i (ξ) can make alternative A i rank r.Thus, the favorable rank weight space W r i (ξ) can be expressed as follows: By analyzing W r i (ξ), the SMAA-2 method obtains four descriptive measures to assist decision-makers in choosing the best alternative.These measures are the rank acceptability index, the central weight vector, the confidence factor, and the holistic acceptability index.
The rank acceptability index b r i is defined as the expected volume of the favorable rank weight space W r i (ξ) and measures the variety of different weights which guarantee a ranking r for alternative A i .This index is a multidimensional integral over criteria distributions and the favorable rank weight space, which can be computed numerically as follows: Obviously, b r i is between 0 and 1; 0 indicates that alternative A i will never acquire a given rank, and 1 indicates that it will always acquire the given rank r with arbitrary weights.The larger the index, the more the decision-makers are willing to accept the given rank r for alternative A i .The best alternatives are those with high acceptability indices for the best ranks.
The central weight vector ω c i is defined as the expected center of gravity of the favorable weight space W 1 i (ξ) = {ω ∈ W : rank(ξ i , ω) = 1}, and it is computed numerically as an integral of the weight vector over the criteria distributions and the weight distributions: The central weight vector ω c i is the best single vector representing the preference of a typical decision-maker who favors alternative A i , and can help decision-makers understand the correspondence between different weight vectors and different decision-maker choices.
The confidence factor p c i is the probability for alternative A i to be the preferred one when the central weight vector ω c i is adopted.It is an integral over criteria distributions: Symmetry 2018, 10, 587 The confidence factor p c i is used to measure whether the criteria data are accurate enough to discern alternative A i based on the central weight vector ω c i .For any given weight vector and alternative, the corresponding confidence factor can be computed in a similar manner.The confidence factor can be described as the proportion of the stochastic criteria space that makes the alternative the best one with the given weight vector.
The holistic acceptability index e i is a comprehensive index to examine the overall acceptability of alternative A i by combining all rank acceptability indices as follows: where α r is a metaweight.Lahdelma et al. [32] suggested choosing the alternative with the highest holistic acceptability index as the best one.In addition, considering the requirements for metaweights, i.e., they should be nonnegative, normalized, and nonincreasing when the rank increases ( Lahdelma et al. proposed three kinds of metaweights, namely linear weights α r = (N − r)/(N − 1), inverse weights α r = 1/r, and centroid weights α r = ∑ N i=r 1/i/∑ N i=1 1/i [32].Figure 3 illustrates the differences in linear, inverse, and centroid weights for 10 ranks.It can be observed that linear weights allocate more weight to the mediocre ranks, while inverse and centroid weights emphasize the best ranks.Compared to centroid weights, inverse weights are more evenly distributed over the worst ranks, which indicates that the aggregation by inverse weights is more insensitive to the order among the worst ranks.Therefore, in this research, we use centroid weights.
the central weight vector c i ω is adopted.It is an integral over criteria distributions: The confidence factor c i p is used to measure whether the criteria data are accurate enough to discern alternative i A based on the central weight vector c i ω .For any given weight vector and alternative, the corresponding confidence factor can be computed in a similar manner.The confidence factor can be described as the proportion of the stochastic criteria space that makes the alternative the best one with the given weight vector.
The holistic acceptability index i e is a comprehensive index to examine the overall acceptability of alternative i A by combining all rank acceptability indices as follows: where r  is a metaweight.Lahdelma et al. [32] suggested choosing the alternative with the highest holistic acceptability index as the best one.In addition, considering the requirements for metaweights, i.e., they should be nonnegative, normalized, and nonincreasing when the rank increases ( 120  The multidimensional integrals in Equations ( 19)-( 21) are almost impossible to compute analytically because the distributions ( ) f Ω and ( ) g W vary according to their application and can be arbitrarily complex.For this reason, Tervonen et al. proposed to employ the Monte Carlo simulation to compute approximate values for high-dimensional integrals, where the required number of iterations is set as 10,000 to achieve accuracy with 95% confidence.More details about Monte Carlo simulation-based SMAA-2 computations can be found in Reference [34].The multidimensional integrals in Equations ( 19)-( 21) are almost impossible to compute analytically because the distributions f (Ω) and g(W) vary according to their application and can be arbitrarily complex.For this reason, Tervonen et al. proposed to employ the Monte Carlo simulation to compute approximate values for high-dimensional integrals, where the required number of iterations is set as 10,000 to achieve accuracy with 95% confidence.More details about Monte Carlo simulation-based SMAA-2 computations can be found in Reference [34].
Since it was first proposed, the SMAA-2 method has been applied in various fields, such as the evaluation of combined heat and power units [35], task assignment for multi-UAV [31], strategic forest planning [36], reservoir flood control operations [33], project portfolio optimization [37], and more.However, this method still suffers from a notable shortcoming: The application of the method is based on one criteria value matrix.When there are multiple criteria value matrices from different ISs, the SMAA-2 method is difficult to apply directly.In view of the above shortcoming, we extend the SMAA-2 method in this research.The details are introduced in the next section.

The Extended SMAA-2 Method
As can be seen from Figure 1, the extended SMAA-2 method first aggregates multiple criteria value matrices from different information sources (ISs) under the condition that the criteria weight information is unavailable and the weight information for the ISs is partially available in order to acquire the final criteria value matrix.Then, the holistic acceptability indices are computed through the Monte Carlo simulation based on the matrix.To achieve matrices aggregation under serve uncertainty, we propose a novel aggregation method.To compute the holistic acceptability indices, we integrate the aggregation method into the original SMAA-2 method.
This extension is based on the assumption that when one certain criterion of an alternative has been evaluated multiple times, the criteria values obtained are stable and tend to be some of the most possible values.According to the normal distribution interval number theory proposed by Wang et al. [38,39], we assume that the target value v (k) j and the distance d provided by the sensors (or information sources) are represented by normal distribution interval numbers.

Normal Distribution Interval Number and Weighted Arithmetic Averaging Operator
Definition 1. Assuming a = [a − , a + ] is an interval number, if any x ∈ [a − , a + ] follows the same normal distribution N µ, σ 2 , then a is called a normal distribution interval number [38], denoted by β = (µ, σ).If we let I be the set of all normal distribution interval numbers, according to the '3σ rule' of normal distribution (if x ∼ N µ, σ 2 , then P(µ − 3σ ≤ x ≤ µ + 3σ) = 0.9974) [40] and µ and σ are determined as follows: The results are as follows: Based on Definition 1 and the principle that the linear combination of a finite number of independent normal random variables still follows a normal distribution, we can define the computation rules of normal distribution interval numbers as follows: Definition 2. Assuming β 1 = (µ 1 , σ 1 ) and β 2 = (µ 2 , σ 2 ) are arbitrarily two normal distribution interval numbers, then (1) Furthermore, according to the addition commutative law and the multiplication allocating law of rational numbers, it is easy to prove that (1) Based on the distance measure of interval numbers in References [41,42], we define the distance between two normal distribution interval numbers as follows: Definition 3. Assuming β 1 = (µ 1 , σ 1 ) and β 2 = (µ 2 , σ 2 ) are arbitrarily two normal distribution interval numbers, then the distance between them can be computed as follows: Symmetry 2018, 10, 587 where minup = min(µ 1 + 3σ 1 , ) is a group of normal distribution interval numbers, and assuming WAA: is always satisfied for arbitrary , then WAA is called the weighted arithmetic averaging operator for normal distribution interval numbers, and ω is called the weighted vector of is the corresponding weighted vector, then Proof.Here, we use the mathematical induction to prove Theorem 1. ( and ω 2 β 2 = (ω 2 µ 2 , ω 2 σ 2 ) are two normal distribution interval numbers.Therefore, we have: (2) Assuming that when m = k, Theorem 1 holds: (3) When m = k + 1, according to Definition 2: That means Theorem 1 holds when m = k + 1.Therefore, it can be concluded from Equations ( 1) and (2) that Theorem 1 holds.
) are a group of normal distribution interval numbers, then the matrix S = β ij N×n is called a normal distribution interval number matrix.
Here, we consider a multicriteria decision-making problem, including a set of N alternatives Set_alt = {A 1 , A 2 , • • • , A N } from which one is to be chosen.All alternatives are evaluated in terms of n criteria Set_cri = {C 1 , C 2 , • • • , C n } by a set of K experts (or information sources) Set_is = {IS 1 , IS 2 , • • • , IS K }.Thus, there are K criteria value matrices S (1) , S (2) , • • • , S (K) , each of which has N rows and n columns, corresponding to N alternatives and n criteria, respectively.According to the above assumption, i.e., the criteria values are represented by normal distribution interval numbers, the kth criteria value matrix is a normal distribution interval number.In order to facilitate criteria value matrices aggregation, all matrices need to be normalized by Reference [43]: In the following, unless otherwise stated, all criteria value matrices are assumed to be normalized.

Definition 6. S
where λ = [λ 1 , λ 2 , • • • , λ K ] is the integrated weight vector of K information sources (ISs), satisfying λ k ∈ [0, 1] and ∑ K k=1 λ k = 1.The integrated weight of the kth IS λ k is computed as follows: where τ ∈ [0, 1] is an adjustment factor and η k ∈ [0, 1] and ρ k ∈ [0, 1] are the subjective weight and the objective weight of the kth IS, respectively, satisfying ∑ K k=1 η k = 1 and ∑ K k=1 ρ k = 1.The subjective weights are usually predefined and would not be adjusted once determined.The objective weight of each IS is determined by the accuracy and credibility of the criteria values it provides [44].The iterative algorithm for computing objective weights is introduced in the following subsection.

Iterative Algorithm for Computing Objective Weights
As objective weights are related to the accuracy and credibility of criteria values, measuring their accuracy and credibility is critical.In this research, we use the deviation between the criteria value Symmetry 2018, 10, 587 13 of 25 matrix of one IS and the final criteria value matrix as the measurement [44] and then proposed an iterative algorithm to compute the objective weights of the ISs.

Definition 7.
If we let and ) be the criteria value matrix of the kth IS and the final criteria value matrix, respectively, then the deviation between these two matrices can be computed: where is the weight of the jth criterion, satisfying ω j ∈ [0, 1] and ∑ n j=1 ω j = 1, and dis β ij , as defined by Definition 3.
As can be seen from Definitions 7 and 3, a smaller deviation Dev S (k) , S (0) represents a higher similarity between S (k) and S (0) , which means that the kth IS should be given a higher objective weight.Based on this perception, the objective weight ρ k is then computed: The iteration end condition is set as ∆ρ ≤ ζ, where ζ is the error margin, a very small number (i.e., 10 −8 ), and ∆ρ is computed as follows: In Equation (39), ρ are the objective weight of the kth IS at the lth and the (l − 1)th iteration, respectively.For convenience, the initial objective weight vector is set as Once the objective weight vector is determined, then the integrated weight vector can be computed by Equation (36).On this basis, all criteria value matrices can be aggregated to obtain the final criteria value matrix by Equations ( 34) and (35).The pseudo-code of the aggregation method based on the iterative algorithm is presented in Algorithm 1.On the basis of the output of Algorithm 1 (i.e., the final criteria value matrix S (0) ), the SMAA-2 method can be subsequently applied.

Input:
number of alternatives, criteria, and information sources: N, n, K; criteria value matrices: S (1) , S (2) error margin used in computing the objective weight vector: ζ.Algorithm 1: Cont.

Integrating the Aggregation Method into the Original SMAA-2 Method
By integrating the aggregation method into the original SMAA-2 method, the stochastic multicriteria acceptability analysis can be performed to solve the task assignment problem formulated in Section 3. The implementation of the extended SMAA-2 method is also mainly based on the Monte Carlo simulation with 10,000 iterations to achieve accuracy with 95% confidence and consists of two phases of which the first is the computation of the rank acceptability index b r i , the central weight vector ω c i , and the holistic acceptability index e i , and the second is the computation of the confidence factor p c i .The pseudo-code is specified in Algorithm 2. RAND Ω (S), function returning a random criteria value matrix from criteria distribution f (Ω); RAND W , function returning a random weight vector from weight distribution g(W); RANK(U), function returning a vector of ranks corresponding to U.

Input:
number of alternatives, criteria, and decision-makers: N, n, K; criteria value matrices: S (1) , S (2) , • • • , S (K) ; subjective weight vector and initial objective weight end for end for // Phase 2, compute p c i // Initialize p c i , and generate the final criteria value matrix S (0) As can be seen from Algorithm 2, for each UCAV, we need to compute the holistic acceptability indices in terms of each target (here, a target is equivalent to an alternative) by the extended SMAA-2 method.In the next section, we convert the task assignment model established in Section 2 into an integer linear programming model.Thus, the holistic acceptability indices can be applied to task assignment.

Converting the Task Assignment Model Based on Holistic Acceptability Indices
The task assignment model established in Section 3 involves M UCAVs and N targets.By applying the extended SMAA-2 method M times, we obtain a holistic acceptability index matrix E = e ij M×N , where e ij (i = 1, 2, • • • , M; j = 1, 2, • • • , N) is the holistic acceptability index for assigning the jth target T j to the ith UCAV UCAV i .Then, we can convert the established task assignment model into an integer linear programming model in which the objective function aims to maximize the holistic acceptability index of the task assignment scheme.The converted model is expressed as follows: x ij e ij (40) subject to the following: The reasons the previous model can be converted in this way are explained as follows: The value of the objective function ) in the previous model is an absolute value of combat reward.However, due to uncertain parameters, the value is not able to be computed directly.According to the SMAA-2 method, Lahdelma et al. suggested choosing the alternative with the highest holistic acceptability index [32].In other words, e ij in the task assignment model, which describes the overall acceptability for assigning T j to UCAV i , can be interpreted as a score of the assignment.For UCAV i , the larger the score, the more UCAV i should attack T j to achieve a larger combat reward given the criteria distributions f (Ω) and the criteria weight distributions g(W).Thus, the holistic acceptability index e ij can be treated as a relative value of combat reward.Therefore, it is reasonable to convert the objective function F(x) in the previous model to a linear combination of e ij , expressed as Equation (40).
In terms of constraints, as the constraints expressed in Equations ( 7)-( 12) on the criteria weights and the aggregation method of the criterion value matrices in the previous model have been satisfied in computing the holistic acceptability indices, the converted model only contains the three constraints represented by Equations ( 41)-( 43) that specify the range of the decision variable x ij .Thus, there are no uncertain parameters in the converted model, which indicates an easier process of searching for solutions.

Case Study
In this section, we examine the proposed method through a case study.Here, we consider a combat scenario where four UCAVs under the command of one command center are to attack four geographically dispersed ground targets.The initial battlefield situation is shown in Figure 4. Before the mission begins, the command center has received target information that is independently acquired by four different sensors through multiple reconnaissances, as shown in Table 1.Other necessary parameters are shown in Tables 2 and 3.

Sensors
Target Value ( 7.0, 9.0 6.0,8.0 3.0, 7.0 2.0, 5.0  Sensor 3 5.0, 9.0 5.0,8.5 4.0, 9.0 5.0,8.0 6.0,8.0 4.0, 6.5 3.0,8.5 4.5, 7.5    The holistic acceptability index matrices in the above two case are illustrated in Figure 5a,b, respectively.Table 6 gives the specific values of two matrices.It may be noted that the holistic acceptability index e ij in Table 6 does not indicate the actual combat reward of UCAV i attacking T j , but the relative value of the combat reward under the given distributions of operational objectives and the given distributions of operational objective weights.For UCAV i , the larger e ij , the more UCAV i should attack T j to achieve a larger combat reward.Based on the holistic acceptability index matrices in Table 6, task assignment schemes are obtained in the two cases, as shown in Table 7. Figure 6 shows the task assignment diagrams corresponding to the schemes in Table 7.As can be seen from Figure 5a and Table 6, when commanders are more concerned with the operational benefit, ω B ≥ ω C , attacking T 1 can bring the largest holistic acceptability index for each UCAV.The reason for this can be found in Tables 4 and 5.It can be seen from the second column of Table 4 that the average value of T 1 is the largest according to the target value information obtained by each sensor.Moreover, in Table 5, for each UCAV, the utilities (or normalized operational benefits) obtained by attacking T 1 are the largest, all being equal to 1.However, due to the constraints expressed by Equations ( 41) and ( 42) in Section 6, T 1 is assigned to UCAV 2 in order to maximize the holistic acceptability index of the task assignment scheme shown in the left part of Table 7.
When commanders are more concerned with the operational cost, ω B < ω C , the holistic acceptability indices of UCAV 2 attacking T 1 and UCAV 3 attacking T 1 are almost unchanged, both being equal to 1, which indicates that for both UCAV 2 and UCAV 3 , attacking T 1 is the best choice.On the other hand, for UCAV 1 and UCAV 4 , the best choice is attacking T 3 instead of T 1 .The reason for the change is that commanders become more inclined to consider the operational cost in this case.Due to the same reason, the task assignment scheme varies, as shown in the right part of Table 7.The holistic acceptability index matrices in the above two case are illustrated in Figure 5a,b, respectively.Table 6 gives the specific values of two matrices.It may be noted that the holistic acceptability index ij e in Table 6 does not indicate the actual combat reward of i UCAV attacking j T , but the relative value of the combat reward under the given distributions of operational objectives and the given distributions of operational objective weights.For i UCAV , the larger ij e , the more i UCAV should attack j T to achieve a larger combat reward.Based on the holistic acceptability index matrices in Table 6, task assignment schemes are obtained in the two cases, as shown in Table 7. Figure 6 shows the task assignment diagrams corresponding to the schemes in Table 7.As can be seen from Figure 5a and Table 6, when commanders are more concerned with the operational benefit, BC   , attacking 1 T can bring the largest holistic acceptability index for each UCAV.The reason for this can be found in Tables 4  and 5.It can be seen from the second column of Table 4 that the average value of 1 T is the largest according to the target value information obtained by each sensor.Moreover, in Table 5, for each UCAV, the utilities (or normalized operational benefits) obtained by attacking 1 T are the largest, all being equal to 1.However, due to the constraints expressed by Equations ( 41) and ( 42) in Section 6, 1 T is assigned to 2 UCAV in order to maximize the holistic acceptability index of the task assignment scheme shown in the left part of Table 7.
When commanders are more concerned with the operational cost, T instead of 1 T .The reason for the change is that commanders become more inclined to consider the operational cost in this case.Due to the same reason, the task assignment scheme varies, as shown in the right part of Table 7.It can be seen from the above analysis that although the actual operational reward cannot be computed due to multiple uncertainties in the task assignment problem, we can still achieve multiobjective task assignment for multi-UCAV.Specifically, we have transformed the actual operational reward of assigning  6.By maximizing the sum of the holistic acceptability indices of all assignments, we have found the optimal task assignment scheme.In addition, it can also be seen that as the preferences of commanders varies, the holistic acceptability index of the same assignment changes.For example, in both cases, 1 UCAV attacks 1 T , indicating that the actual operational benefit and the actual operational cost have remained unchanged.However, due to the variation of the preferences of commanders, the utility values of the actual operational benefit and the actual operational cost have changed, leading to the holistic acceptability index of assigning 1 T to 1 UCAV falling from 0.9664 to 0.5689.A similar situation could also be seen in other assignments.Ultimately, the changes in the holistic acceptability indices of all assignments result in the adjustment of the task assignment scheme, which indicates that that the preferences of commanders have an important impact on operational decision-making.
In summary, the proposed method in this research is able to solve the multiobjective task assignment problem for multi-UCAV under multiple uncertainties and reveal, to some extent, the influence of the variation in the preferences of commanders on operational decision-making.In our research, we only considered two objectives: the operational benefit and the operational cost.However, more objectives can be considered, such as maximizing the operational benefit, minimizing the operational cost, and minimizing the operational loss of UCAVs.In this case, we would need to establish more complex models and acquire more target information.It can be seen from the above analysis that although the actual operational reward cannot be computed due to multiple uncertainties in the task assignment problem, we can still achieve multiobjective task assignment for multi-UCAV.Specifically, we have transformed the actual operational reward of assigning T j to UCAV i into the holistic acceptability index e ij of the assignment, as shown in Table 6.By maximizing the sum of the holistic acceptability indices of all assignments, we have found the optimal task assignment scheme.In addition, it can also be seen that as the preferences of commanders varies, the holistic acceptability index of the same assignment changes.For example, in both cases, UCAV 1 attacks T 1 , indicating that the actual operational benefit and the actual operational cost have remained unchanged.However, due to the variation of the preferences of commanders, the utility values of the actual operational benefit and the actual operational cost have changed, leading to the holistic acceptability index of assigning T 1 to UCAV 1 falling from 0.9664 to 0.5689.A similar situation could also be seen in other assignments.Ultimately, the changes in the holistic acceptability indices of all assignments result in the adjustment of the task assignment scheme, which indicates that that the preferences of commanders have an important impact on operational decision-making.
In summary, the proposed method in this research is able to solve the multiobjective task assignment problem for multi-UCAV under multiple uncertainties and reveal, to some extent, the influence of the variation in the preferences of commanders on operational decision-making.In our research, we only considered two objectives: the operational benefit and the operational cost.However, more objectives can be considered, such as maximizing the operational benefit, minimizing the operational cost, and minimizing the operational loss of UCAVs.In this case, we would need to establish more complex models and acquire more target information.

Conclusions
In this paper, we investigated the task assignment problem for multi-UCAV under multiple uncertainties, including uncertain target information provided by multiple sensors, uncertain weights of the criteria or objectives concerned, and partially known weights of sensors.To find solutions for the problem, we extended the SMAA-2 method and then combined it with the integer linear programming method.Through a case study, the feasibility of the proposed method was demonstrated.
The most significant contribution of our research is the extension of the original SMAA-2 method.The new method possesses the following strengths and weaknesses:

•
Due to the integration of the WAA operator-based aggregation method, the extended method is able to perform stochastic multicriteria acceptability analysis under circumstances of multiple criteria value matrices.This is also the most prominent advantage of the extended method over the original SMAA-2 method.

•
The extended method employs the iterative algorithm to compute the objective weight of a sensor based on the accuracy and credibility of the target information it provides so as to minimize the subjective judgment errors of experts on the sensor performance.By this means, the actual importance of a sensor can be reflected to the greatest extent so that the final criteria value matrix obtained could be as accurate as possible.Through these measures, the credibility and accuracy of task assignment schemes can be further increased.

•
The foundation of the extended method is the stochastic multicriteria acceptability analysis, so it inherits the advantage of the SMAA technique in that it does not require decision-makers to express their preferences explicitly [30].This property is particularly useful in operational planning because the commander could only express implicit preferences for multiple operational objectives most of the time.

•
The extended method is based on the assumption that the uncertain target information can be represented by normal distribution interval numbers.This assumption may limit the application of the method.When the uncertain target information is represented by fuzzy sets or rough sets, the extended method cannot be applied.

•
The extended method employs a linear utility value function to map criteria values to utility values.Although Lahdelma et al. have proved it may often be safe to apply a linear utility value function, they also mentioned that when the convexity or concavity of the utility value function is severe, the results would be slightly different [45].This means that in some special cases, the linear utility value function may not be the best choice, and we would need to construct a more appropriate one.
In future research, we will focus on improving the weaknesses listed above in the extended SMAA-2 method so that the combination method can be used to solve similar multiobjective resource allocation problems under serve uncertainties in other fields.
Lahdelma et al.  proposed three kinds of metaweights, .Figure3illustrates the differences in linear, inverse, and centroid weights for 10 ranks.It can be observed that linear weights allocate more weight to the mediocre ranks, while inverse and centroid weights emphasize the best ranks.Compared to centroid weights, inverse weights are more evenly distributed over the worst ranks, which indicates that the aggregation by inverse weights is more insensitive to the order among the worst ranks.Therefore, in this research, we use centroid weights.
through the WAA operator as follows:

Algorithm 2 :
Pseudo-code of the extended SMAAof times alternative A i is evaluated into rank j in the Monte Carlo simulation; N ite , number of iterations, set as 10000; r = [r 1 , r 2 , • • • , r N ], vector of ranks of alternatives; U = [U 1 , U 2 , • • • , U N ], utility value vector of alternatives; Agg(ω), function returning the final criteria value matrix S (0) based on the criteria weight vector ω (details are presented in Algorithm 1)

3 UCAV attacking 1 T 2 UCAV and 3 UCAV , attacking 1 T 1 UCAV and 4 UCAV
are almost unchanged, both being equal to 1, which indicates that for both is the best choice.On the other hand, for , the best choice is attacking 3

Figure 5 .
Figure 5.The holistic acceptability index matrices when (a)

Figure 6 .
Figure 6.Task assignment diagrams when (a) holistic acceptability index ij e of the assignment, as shown in Table
error margin used in computing the objective weight vector: ζ; number of iterations of the Monte Carlo simulation: N ite .

Table 6 .
The holistic acceptability indices when

Table 6 .
The holistic acceptability indices when ω B ≥ ω C and ω B < ω C .

Table 7 .
Task assignment schemes when ω B ≥ ω C and ω B

Table 7 .
Task assignment schemes when B ijx .