**5. Conclusion**

In this chapter, we have presented the mathematical investigation of static pattern formation across a resistively coupled RD-MCNN structure *via* the application of the theory of local activity. The considered networks have identical cells in a compact form, each of which is composed of a DC bias voltage source, a bias resistor, a locally active memristor, and a parallel capacitor. The memristor was represented through a generic model, which helped to reduce the numerical complexity and to shorten the simulation time of large-scale networks. A useful AC equivalent circuit could be derived for the locally active device, which facilitated further calculations related to the small signal model of the network cell. We have adopted a systematic circuit theoretical approach that we applied for the stability analysis of the isolated cell and for the extraction of its parameter values for its operation on the edge-of-chaos and sharp-edge-of-chaos domains. In this way, we have performed a simple, fast, and robust analysis, which at the same time allowed an efficient interpretation of the results. All the calculations were performed in parametric form, which allowed a deep investigation of the results, making it possible to extract the related sets of parameters in terms of the cell characteristics, namely the DC operating point and the capacitor value, as well as the value of the

**Figure 9.**

*Simulation results of various m* � *n RD-MCNN structures for Rb* ¼ 2*k*Ω*, Rc* ¼ 0*:*2*k*Ω *and Vb* ¼ 1*:*2*V, where m* ¼ *n* ¼ 31 *in (a), m* ¼ *n* ¼ 41 *in (b), m* ¼ *n* ¼ 51 *in (c), and m* ¼ *n* ¼ 61 *in (d). The corner cells have the same initial conditions as the center cell for all the structures. Although the four structures share a square geometry, they differ in size, which leads to clearly distinguishable patterns.*

coupling resistance. The proposed RD-MCNN was shown to generate diverse patterns while we have extensively investigated the effect of design parameters, initial conditions, and the size of the network on the patterns generated therein. Furthermore, the mathematical and circuit theoretical approach employed during the investigation of pattern formation dynamics in the RD-MCNN under focus in this chapter can be easily adopted to RD-MCNNs with different cell designs.

Finally, the content of this work can be extended to the investigation of the role of the network geometry and variety in cells' initial conditions on pattern formation. Furthermore, the impact of memristor non-idealities such as variability in memristor behavior from cycle to cycle, as well as from device to device, and endurance degradation and short- and long-term reliability issues shall be further examined to achieve a robust design. Such non-idealities, essentially, can affect the NDR characteristics of the memristors and narrow the width of the NDR region, or the devices can even get stuck in one of the high-resistance or low-resistance locally passive regimes where they act as dead cells without any dynamics. Preliminary simulation results show that endowing memristors with narrower NDR width negatively affect the patterns' color contrasts, while the presence of the high- or lowresistance locally passive devices results in the formation of color clusters within the patterns. Interestingly, similar color clusters emerge also when some of the array memristors' initial conditions are randomly selected, which suggests a possible strategy to compensate for the disturbing action of locally passive cells on the formation of predefined patterns through a dynamic conditioning of their initial conditions. Similarly, the low color contrast quality in the patterns, which occur due *Pattern Formation in a RD-MCNN with Locally Active Memristors DOI: http://dx.doi.org/10.5772/intechopen.100463*

to the existence of memristors featuring a narrower NDR region, can be improved by reprogramming dynamically the bias voltage sources in these "defect" cells. In conclusion, proper control strategies including cells' dynamic biasing and initial condition conditioning, as well as reconfiguration of the network array, may be considered as a solution to overcome the harmful effects which memristor nonidealities may induce on the emergence of predefined patterns in the proposed RD-MCNNs. Future studies will also be devoted to envision an application, where the capability of our cellular medium to generate a variety of steady-state static patterns may be useful to our modern society.


**5. Conclusion**

482 Lipoproteins – Role in Health and Diseases

PON, nmol/ mlmin

> 174 2\*

194 3\*,\*\*

> 173 3\*

189 2\*,\*\*

> 159 1\*

179 2\*,\*\*

> 223 6

> 215 3\*

> 221 4

> 255 3

> > 3

Ascorbic acid, mkmol/L

> 37.26 1.60\*,\*\*

> 41.79 0.44\*,\*\*

> 40.08 1.72\*,\*\*

> 54.30 0.97\*,\*\*

> > 35.46 1.10\*

48.59 0.96\*,\*\*

> 75.95 0.74\*

> 74.79 0.34\*

67.04 1.40

70.61 2.31

79.49 2.18\*

**Table 25.** The effect of grape polyphenol concentrates and grape wines on the plasma oxidant/antioxidant status in rats with the emotional-painful stress (in each group n=10).

Parameter

Isolated double bonds in ApоВ-LP

> 2.07 0.09\*

> 2.35 0.15\*

> 1.92 0.25\*

2.47 0.03\*,\*\*

> 1.52 0.06\*

> 2.21 0.07\*

> 5.43 0.12\*

> 5.55 0.15\*

4.83 0.11

5.39 0.04\*

5.75 0.17\*

Diene conjugates in ApоВ-LP

> 27.03 1.30\*,#

26.86 0.43\*

18.45 0.18\*,\*\*

14.50 0.35\*,\*\*

25.73 0.68\*,\*\*

27.68 0.60\*,\*\*

> 15.67 0.15\*

> 15.81 0.15\*

16.87 0.21

15.47 0.34\*

17.82 0.24

> 7.77 0.21\*

8.25 0.66\*,\*\*

8.91 0.29\*\*

10.42 0.16\*,\*\*

> 7.04 0.07\*

8.03 0.06\*,\*\*

> 10.03 1.72

11.52 0.24\*

9.49 0.24

9.97 0.16

10.22 0.36

Group

Str.+Con. Podarok Magaracha (AD)

Str.+Con. Krasen (АD)

Srt.+Wine Podarok Magaracha

Srt.+Wine Krasen

Str.+Con. Podarok Magaracha (DW)

Str.+Con. Krasen (DW)

Wine Podarok Magaracha

> Wine Krasen

Con. Podarok Magaracha (AD)

Con. Krasen (AD)

Ethanol 236

The data presented as mean±SD \* - р0.05 versus intact animals \*\* - р0.05 versus stressed animals Based on our findings, it is possible to state that antioxidant complexes, particularly polyphenol extracts and the concentrates obtained from *Vitis Vinifera*, which are safe and reveal the potent antioxidant and stress-protective activity, should be used for reduction of proaterogenic states consequences in the complex prophylactic and treatment of atherosclerosis as effective stress-protective remedies.

Thus*,* administration of *Vitis Vinifera* substances can prevent the increase of the total lipoprotein and ApoB-LP content in the blood, and prevent the free radical process activation in the plasma lipoprotein particles, and, in general, normalize the lipid and lipoprotein metabolism in the liver in metabolic syndrome. These results have proven the ability of the investigated complexes to reduce such negative consequence of metabolic syndrome as development of atherosclerosis.

In addition, according to obtained research data the polyphenolic concentrates possess a potent protective activity both in acute and chronical neurogenic stress.

Our studies suggest that multicomponent active substances with antioxidant properties are more effective in correction of the proatherogenic states caused by stress and metabolic syndrome negative effects in comparison with individual antioxidants (particularly, αtocopherol). The research data suggest that the increased plasma antioxidant activity alone does not result in decreased foam cell formation, at least in the studied animal model. Moreover, *in vitro* studies have shown that α-tocopherol can be pro-oxidative rather than protective for lipids in isolated LDL. Similarly with vitamin E, vitamin C additives do not offer consistent benefit against atherosclerosis in animals.

The occurrence of tocopherol-mediated peroxidation and the mode of its prevention predicts that the balance of α-tocopherol and available coantioxidants, rather than αtocopherol alone, determines whether LDL lipid peroxidation occurs in biological systems. Inhibition of the free radical process with the polyphenolic complexes administration can be associated with their ability to increase the level of antioxidants – α-tocopherol, ascorbic acid and reduced glutathione in the test animal liver tissue compared with the group of the stressed animals. The complexes obtained from *Vitis Vinifera*, in particular, polyphenolic concentrates "Enoant" and "Polyphen", as well as grape wines (particularly "Cabernet") with their moderate use revealed the potent antioxidant activity. The preliminary results also suggest that coantioxidants inhibit lipoprotein lipid peroxidation *in vivo*. Thus, if LDL oxidation causes atherosclerosis, the requirement for coantioxidants may explain why supplementation with individual antioxidants, particular vitamin E alone, overall has yielded inconclusive results in the controlled human and animal intervention studies.

In conclusion, our research results may be used for the development of the atherosclerosis prophylaxis strategy, and treatment of diabetes mellitus and metabolic syndrome because recent studies proved insufficient effectiveness of α-tocopherol and advantages of multicomponent antioxidant complexes administration. The high effectiveness of the polyphenolic complexes obtained from *Vitis Vinifera,* including polyphenolic concentrates

"Enoant" (from grape of "Cabernet" grade) and "Polyphen" (from grape of "Rkatsiteli" grade) produced by the National Institute of Grape and Wine "Magarach" has been proven. Our results also confirmed the high effectiveness of the antioxidant complexes from grapes in the correction of the endothelial dysfunction, thus, including these extracts in the treatment schemes is very reasonable. As it would be expected from our observations, increasing the antioxidant oxidant defense by antioxidant supplementation has the ability to restore the endothelial vasomotor function.

An important question to be asked is whether the polyphenol antioxidants exerted their inhibitory effect on lesion progression only because of their antioxidant properties or, possibly, because of additional biological properties, in particular – the phytoestrogen activity.

However, further studies, especially in humans, are required to validate the role of these antioxidants in inhibiting LDL oxidation.

Nevertheless, there are some limitations in the use of the concentrates produced from red grade grapes because of the uric acid content changes.


*2.3. Mechanism*

Conversion of free FAs into CoA esters constitutes an initial activation step before peroxisomal β-oxidation. This reaction is catalyzed by specific acyl-CoA synthetase connected to the cytosolic side of the peroxisomal membrane [69]. It was proved, using protease protection assays, that acyl-CoAs but not free FAs bind to the TMD of the transporter [70]. It is therefore only after activation that the fatty acyl-CoAs are transported to the peroxisomal matrix through peroxisomal ABC transporters. Fatty acyl-CoA are captured on the cytosolic side by the TMD, enhancing the affinity of NBD for ATP. ATP molecules are then hydrolyzed, thus producing the energy needed to switch the conformation of TMD and eventually allowing the translocation of substrates from the cytosol into the peroxisomal

matrix [22,71]. However, the exact mechanism of transport remains controversial. Two models are commonly considered. The first implies that esterified FAs are delivered directly to the peroxisomal matrix, whereas in the other model, free FAs are transported into the peroxisomal matrix after the hydrolysis of acyl-CoAs, which are re-esterified by acyl-CoA synthetase when in the peroxisomal lumen.

Although the process of cleavage and reactivation of acyl-CoAs seems to be a waste of energy as two ATP molecules are needed for the activation reaction, such a mechanism is crucial for the specific permeabilization of the substrates of β-oxidation [18]. Several studies have been done in an attempt to figure out the correct model for the transportation mechanism. Early studies on yeast models have demonstrated that fatty acyl-CoAs are hydrolyzed before being transported. This hydrolysis occurs when acyl-CoAs interacts with the heterodimer Pxa1p-Pxa2p at the cytosolic side of the peroxisomal membrane [72]. The peroxisomal ABC transporters would release a free fatty acid that should be re-esterified inside the peroxisome before its catabolic processing. In addition, intrinsic acyl-CoA thioesterase activity has been found in COMATOSE (CTS), a homolog of human ABCD1 in *Arabidopsis thaliana*, proving again that VLCFA-CoA is hydrolyzed prior to transport [73]. Very recently, the work of Kawaguchi et al. provided further proof of the transport mechanism [74]. After expressing human His-tagged ABCD1 in methylotrophic yeast, they directly demonstrated that ABCD1 transports the FA moiety after the hydrolysis of VLCFA-CoA and that acyl-CoA synthetase is required before the β-oxidation of VLCFA-CoA within the peroxisomes. When it comes to the fate of the free CoA, they are released in the peroxisomal lumen, as revealed using isolated peroxisomes from *Saccharomyces cerevisiae* [75].

Finally, after their re-esterification, substrates are directly delivered to specific acyl-CoA oxidases to initiate the β-oxidation process. Peroxisomal acyl-coenzyme A oxidase 1 (ACOX1) catalyzes the first and rate-limiting step of the β-oxidation pathway dedicated to straight-chain fatty acids, which includes LCFAs, VLCFAs, PUFAs, and dicarboxylic acids [76]. Other acyl-CoA oxidases also exist, ACOX2 and ACOX3. ACOX2 is specific to bile acid intermediates [76] whereas the oxidation of branched-chain FAs depends on both ACOX2 and ACOX3 enzymes [77,78]. Of note, mitochondria catalyze the β-oxidation of the majority of short, medium, and long chain FAs but not that of VLCFAs [79]. In yeast and plants, this process of FA β-oxidation occurs exclusively in peroxisomes, whereas in higher eukaryotes, the catabolism of VLCFAs is initiated solely in the peroxisomes [7,80].


**1. Introduction**

Disaster management remains a crucial challenge for the international community. National and supranational institutions and non-governmental organizations have to deal with high impact disasters such as the 2010 Haiti earthquake, the 2013 Typhoon Haiyan or the 2018 Sulawesi earthquake and tsunami. This global concern is reflected in the scientific literature and specifically in the Operations Research community, where the increasing interest on humanitarian logistics and disaster management entails the publication of surveys [**? ?** ], books [**? ?** ], special issues [**? ?** ] and journals (such as the Journal of Humanitarian Logistics and Supply Chain Management).

Within the disaster management cycle, the response phase includes the actions to be performed once the disaster has occurred. Research on the response phase may focus on specific types of disasters, such as earthquakes [**? ?** ] or floods [**?** ], or can address particular problems for generic disasters, such as location-allocation [**?** ], evacuation [**?** ], or debris management [**?** ]. Among the operations to be carried out after a disaster strikes, the distribution of aid to the affected people plays an important role. Along with the usual objectives in commercial logistics, mainly the operation cost, in humanitarian logistics other attributes have to be considered [**? ?** ], such as the response time, the equity of the distribution, etc. Therefore, multi-criteria models are necessary to face many realistic problems in this context. **?** ] provide a review on multi-criteria decision making in humanitarian aid classifying the existing models as deterministic or stochastic and pre-disaster or post-disaster.

Many times the distribution operations must be performed under unsecure conditions [**?** ]. Assaults to relief vehicles are sadly frequent and this complicates the work of organizations that execute the aid missions. Avoiding as much as possible the most dangerous roads and forcing the vehicles to travel together forming convoys are two ways of improving the security of the operations.

In this paper, we develop an algorithm based on Ant Colony Optimization (ACO) metaheuristic to solve the multi-criteria last mile distribution problem in an unsecure environment stated in [**?** ]. The complexity of this problem makes it necessary to use a heuristic approach to solve realistic cases. Besides, due to the strong links between the elements of a solution, it is difficult to develop heuristics that perform progressive alterations of a solution while ensuring feasibility, as it is the case of evolutionary algorithms or local search heuristics. Instead, algorithms based on construction of solutions, such as Greedy Randomized Adaptive Search Algorithm (GRASP) or Ant Colony Optimization metaheuristics, fit much better.

The ACO methodology that we propose in this work presents some novelties that make it suitable for routing and distribution problems in which the vehicles must travel in convoys. The use of convoys is appropriate not only in humanitarian operations, but in any logistics operation performed in an environment of insecurity, in which attacks on delivery vehicles to steal the commodities are possible. The algorithm considers two types of ants—vehicles and aid kits—four types of standard pheromones, and especially the use of effective pheromones, updatable in the course of each solution building process. This new concept allows a better balance between the elements of a solution and a higher diversification in the set of solutions that can be built. Effective pheromones may be applied to a wide variety of problems, especially vehicle routing problems involving a large number of vehicles. The results obtained show that the new ACO Algorithm proposed in this paper is a useful decision aid tool. The ACO Algorithm provides good quality solutions in a short time in realistic test cases, frequently improving the solutions obtained with the GRASP Algorithm developed in [**?** ].

The rest of this paper is organized as follows. Section **??** reviews the relevant literature on multi-criteria optimization for the last mile distribution problem, as well as some applications of Ant Colony Optimization in disaster management. Section **??** presents the problem of transportation for last mile distribution, including the data notation and the description of the six attributes considered. Section **??** introduces the heuristic method proposed to solve the problem. An Ant Colony Optimization algorithm is developed together with the subordinate procedure that guides the construction of feasible solutions. Section **??** is devoted to showing and analyzing the results obtained when applying the proposed metaheuristic to solve two test cases based on real disasters. Finally, Section **??** draws some conclusions derived from this work.

#### **2. Literature Review**

Some of the most significant papers that apply multi-criteria techniques in humanitarian last mile distribution are mentioned below and summarized in Table **??**. **?** ] develop a multicommodity mixed integer programming model in order to maximize the total population covered and minimize the total travel time for earthquake response. **?** ] present a multicommodity, multimodal and multiperiod last mile distribution problem and establish a scalar objective function that aggregates minimization of transportation costs and maximization of demand covered. **?** ] consider three attributes (reliability of the transporting tours, coverage and total travel time) in a distribution problem where infrastructure may be damaged in a post-disaster situation. **?** ] use the scalarization approach to deal with a three-objective (total unsatisfied demand, total travel time and equity of the distribution) integer programming model that allows split deliveries. Two heuristic methods are proposed to solve the model. **?** ] propose a robust stochastic compromise programming model for disaster relief logistics, minimizing the sum of the expected value and the variance of the total cost and minimizing the sum of the maximum unsatisfied demands. **?** ] consider four attributes (amount of aid distributed, time of the operation, equity and cost) in a dynamic flow model that is solved via lexicographical goal programming. **?** ] develop a quadratic multi-objective programming model to face an allocation and

distribution problem. Three attributes are considered: lifesaving utility, delay cost and equity of the distribution. This model allows vehicles to visit a node more than once and provides the individual itineraries of the vehicles. **?** ] aggregate the distribution time, the penalty cost of unsatisfied demand and the costs of opening distribution centers into a single objective function for a location-routing problem considering network failure, random travel times and multiple uses of vehicles. **?** ] apply a multi-objective evolutionary algorithm to solve a multiperiod and multicommodity distribution model which aims to minimize the unsatisfied demand and to minimize the risk of choosing damaged roads. **?** ] present a multimodal model considering dynamic and stochastic demands in order to minimize the travel time and the unmet demand. **?** ] develop a biobjective multiperiod stochastic location-transportation model that allows multiple trips to the vehicles not only in different periods, but also in the same period. The authors propose different heuristic techniques to solve the model. Outside the scope of aid distribution, **?** ] address a post-disaster waste management problem through a possibilistic programming model that minimizes both the cost of waste processing and the carbon emissions and maximizes the job opportunities.

Despite the frequency with which aid distribution is to be carried out in unsafe conditions, only a few papers consider the security as an attribute. **?** ] integrate seven attributes (quantity distributed, time of response, cost, equity, priority, reliability and security) into a linear lexicographical goal programming flow model. This model is multimodal, multidepot, allows split deliveries and vehicles travel in convoys so that they can be escorted. An alternative version of this model is established in **?** ], where the authors set a fixed amount to distribute and combine the six remaining attributes into a non-lexicographical model. Nevertheless, these models do not provide the individual route of each vehicle, so it is not possible to easily implement plans for the relief distribution. **?** ] propose a GRASP metaheuristic to solve a nonlinear multi-objective model that allows vehicles to visit a node more than once and provides the individual itineraries of the vehicles. The mathematical formulation of this model is established in **?** ], where an extensive comparison of the characteristics of the models is made with other studies that deal with the last mile distribution problem.

**Table 1.** The literature related to multi-criteria techniques in humanitarian last mile distribution.

In this study, we propose an Ant Colony Optimization methodology to solve the model developed in [**?** ] and [**?** ]. Even though ACO has been applied frequently to different routing and distribution problems [**????** ], it has barely been used in the context of disaster management. The most remarkable work is **?** ], where the authors apply this metaheuristic to solve a distribution and evacuation problem. More recently, **?** ] combine different techniques, including ACO, to face a location-allocation problem

of relief centers after an earthquake. Therefore, this paper presents the first application of Ant Colony Optimization to a multi-criteria last mile distribution problem. Furthermore, the proposed ACO methodology introduces novelties such as specialized ants, representing both vehicles and aid kits, and effective pheromones, that conveniently balance the elements of the built solutions.

#### **3. Problem Description**

The last mile distribution problem addressed in this paper consists in deciding how to deal out a predetermined amount of humanitarian aid from distribution centers to populations in need by moving a set of vehicles through a transportation network. We consider that the distribution is performed in an insecure and uncertain environment and, as a consequence, the convoys transporting the aid may be attacked while traveling and some roads may be unavailable due to damages produced by the disaster.

Table **??** lists the notation used to represent the main elements of the problem. The transportation network is formed by nodes and arcs. Each node is of one of three types: supply, demand or connection. A supply node represents a distribution center where a certain known amount of aid is available, while the demand nodes represent populations in need or help centers. Each demand node is associated with a certain amount of aid that is desirable to receive and a priority level. A connection node represents an intermediate point, where there is neither available nor demanded aid, but may be used by the vehicles for aid trans-shipment.

**Table 2.** Problem notation.

The arcs represent links (roads, pathways, streets, etc.) connecting the nodes. The links that can be traversed in both directions have two associated arcs, one for each direction. Each arc has associated a length, a speed, an assault probability and a probability of being available (the arc can be used).

The vehicles are initially located at some nodes, which are not necessarily supply nodes. We will call depots to the nodes where there are vehicles at the beginning of the operation. There are different types of vehicles to carry out the distribution. Vehicles of the same type have the same capacity, velocity and costs. For technical reasons, we classify the vehicles in classes: vehicles belonging to a same class are of the same type and are located also in the same depot.

The humanitarian aid is packed in kits composed of water, food, medicines, blankets, etc. The operation to be performed consists of delivering *qglobal* of such kits.

In addition, the design of the routes must be carried out under the following assumptions:

− Vehicles can visit a node more than once. Likewise, they can pick up or deliver humanitarian aid more than once, even in the same node.

As usual in humanitarian logistics, this problem is intrinsically multi-criteria. In this case, we consider six different attributes that will be incorporated into the model as objectives to be minimized. In what follows, a description of the attribute measures used in the model is provided; however, since they are not necessary for the understanding of the paper, we decided to omit the mathematical formulation. The interested reader can see [**?** ] for further details.

$$idp = \frac{q \, glabal}{\sum\_{i \in I\_D} dem\_i} \tag{1}$$

The equity measure is defined as the square root of the quadratic deviation of the proportions of satisfied demand from the ideal demand proportion over all demand nodes.

The objectives can be minimized separately or in an aggregated function as in (**??**), where weights *wg*, *g* ∈ *G* represent the preferences of the decision maker and *fg* is the objective function associated with the attribute *g*.

$$f = \sum\_{\mathcal{S} \in G} w\_{\mathcal{S}} f\_{\mathcal{S}} \tag{2}$$

For obtaining solutions that take all attributes into account in a balanced way, we apply compromise programming (CP), introduced in **?** ]. CP is a multi-criteria decision making technique which consists in minimizing the distance to the ideal point. The components of the ideal point are the optimal values of the attributes when optimized independently. To normalize the attributes, the difference between the ideal and the anti-ideal is used, where the components of the anti-ideal point are the worst values for the attributes on the set of non-dominated solutions.

The CP objective function when using the *p*-norm is denoted as ˜ *fp* and defined as follows:

$$\tilde{f}\_p = \left[ \sum\_{\mathcal{S}' \in G} w\_{\mathcal{S}} \left( \frac{f\_{\mathcal{S}} - i\_{\mathcal{S}}^+}{i\_{\mathcal{S}}^- - i\_{\mathcal{S}}^+} \right)^p \right]^{\frac{1}{p}} \tag{3}$$

Again, *wg*, *g* ∈ *G* represent the preferences of the decision maker, while *i*+*g* and *i*−*g* are the ideal and the anti-ideal components associated with attribute *g*, respectively.

#### **4. Heuristic Approach**

In this section we present an algorithm based on Ant Colony Optimization, that was introduced in **?** ]. The ACO metaheuristic emulates the way in which the ants of a colony manage to find the shortest path from the anthill to a source of food through the pheromone trails deposited by the ants when they move. Extensive information about ACO metaheuristic and its multiple variants can be found in [**???** ]. The metaheuristic algorithm that we propose in this paper takes some elements of the Ant Colony System (ACS) developed in **?** ] to establish a multiple-ant colony system variant. Section **??** is devoted to the main program of the metaheuristic (ACO Algorithm), and Section **??** to the subordinate procedure that performs the construction of solutions (Pheromone Trail Constructive Algorithm or PTC Algorithm).

Table **??** shows the main parameters and variables that are used in the metaheuristic. According to the ACS, two types of pheromone updates are made. The global update is performed every time *m* consecutive solutions are obtained and is intended to improve the construction of future solutions. The local update is performed after each single solution is obtained, and its purpose is to diversify the construction. Parameters *ρg* and *ρl* are the corresponding evaporation rates. Our algorithm uses two types of ants, representing vehicles and humanitarian aid kits, and four types of pheromones, that are deposited by the vehicles or by the aid kits in different situations. First type (respectively second type) pheromones are deposited by the vehicles (aid kits) when moving through the arcs, third type pheromones are associated with the amounts of aid leaving supply nodes, and fourth type pheromones correspond to the amounts of aid remaining at the nodes at the end of the operation. *τ*1*c*,*k*, *τ*2*k* , *τ*3*i* and *τ*4*i* are the variables that indicate the quantities of pheromone of the four types that have been deposited at an arc *k* or at a node *i*.

The construction of solutions takes into account both pheromones and visibility. Pheromones indicate how often each element candidate to be added to the current partial solution was selected by previous ants, while visibility, established as greedy functions, indicates how good an element could be according to each of the objectives considered.

The PTC Algorithm builds new solutions by adding elements iteratively. As a result, the values of the pheromones associated with the added elements should be decreased to prevent an excessive concentration of these elements in the solution. This is achieved through what we called *effective* pheromones, a new concept that refers to variables that at the beginning of the execution of the PTC Algorithm take the same values as the corresponding pheromones variables, but which will be updated during the construction of the solution.

#### *4.1. Ant Colony Optimization Algorithm*

In this subsection we introduce the ACO Algorithm, whose pseudocode is given in Algorithm **??**. Its main steps are explained below.

#### **Algorithm 1** ACO algorithm

1: Perform preprocess.

```
2: Create a solution S∗ with the PTC Algorithm.
```

```
3: for it = 2, ··· , n do
```

```
14: Output S∗.
```
**Preprocess.** It comprises several tasks needed before starting the construction of solutions:

**Pheromone update.** The local update increases the pheromones of the elements that appear less frequently in the previous iteration in order to diversify the construction, as stated in Equations (**??**)–(**??**).

$$
\tau\_{c,k}^1 = (1 - \rho\_l)\tau\_{c,k}^1 + \frac{\rho\_l}{1 + vt\_{c,k}} \quad \forall c, \; \forall k \tag{4}
$$

$$
\pi\_k^2 = (1 - \rho\_l)\tau\_k^2 + \frac{\rho\_l}{1 + lt\_k} \quad \forall k \tag{5}
$$

$$
\pi\_i^3 = (1 - \rho\_l)\pi\_i^3 + \frac{\rho\_l}{1 + qav\_i - fs\_i} \quad \forall i \in I\_S \tag{6}
$$

$$
\tau\_i^4 = (1 - \rho\_l)\tau\_i^4 + \frac{\rho\_l}{1 + fs\_i} \quad \forall i \in I\_D \cup I\_S \tag{7}
$$

In the same way, the global update increases the pheromones of the elements that appear more in the current best known solution in order to guide the construction towards potentially better solutions, as stated in Equations (**??**)–(**??**). The auxiliary parameter *clc* in (**??**) denotes the number of vehicles of class *c*, while the asterisk in *vt*∗ *c*,*k*, *lt*∗ *k* and *f s*∗ *i* indicates that these auxiliary variables correspond to the current best known solution.

$$
\tau\_{\mathfrak{c},k}^1 = (1 - \rho\_{\mathfrak{F}}) \tau\_{\mathfrak{c},k}^1 + \frac{v t\_{\mathfrak{c},k}^\*}{c l\_{\mathfrak{c}}} \quad \forall c \; \forall k \tag{8}
$$

$$
\pi\_k^2 = (1 - \rho\_\mathcal{g})\pi\_k^2 + \frac{l t\_k^\*}{q g l
abla l} \quad \forall k \tag{9}
$$

$$
\tau\_i^3 = (1 - \rho\_\mathcal{S})\tau\_i^3 + \frac{qav\_i - fs\_i^\*}{qav\_i} \quad \forall i \in I\_\mathcal{S} \tag{10}
$$

$$\tau\_i^4 = \begin{cases} (1 - \rho\_\mathcal{S}) \tau\_i^4 + \frac{fs\_i^\*}{qav\_i} & \forall i \in I\_\mathcal{S} \\\ (1 - \rho\_\mathcal{S}) \tau\_i^4 + \frac{fs\_i^\*}{dem\_i} & \forall i \in I\_D \end{cases} \tag{11}$$

The elements that are taken into account in both updates are the following: for the first type pheromones, the class-arc pairs (**??**) and (**??**); for the second type pheromones, the amounts of transported aid in the arcs (**??**) and (**??**); for the third type pheromones, the amounts of aid that have come out from the supply nodes (**??**) and (**??**); and for the fourth type pheromones, the amounts of aid that remain in the supply and demand nodes at the end of the distribution (**??**) and (**??**).

**Increase** *λ***.** Variable *λ*, which controls the weight of pheromones and visibility, grows over the course of the ACO Algorithm so that the importance of pheromones in the construction is increasing. The pheromones weight is updated according to (**??**), where *it* is the current iteration. The exponential expression makes the increase of *λ* to be progressively slower.

$$
\lambda = 1 - \sigma\_1 e^{-\sigma\_2 \frac{\mu}{n}} \tag{12}
$$

#### *4.2. Pheromone Trail Constructive Algorithm*

The pseudocode of PTC Algorithm is shown in Algorithm **??**. Lines 1 and 2 are the constructive steps, where pheromone and greedy functions are introduced in order to guide the solution building process. Lines 3 to 12 coordinate the elements of the current solution and determine if it is feasible. It must be taken into account that, as all vehicles travel in convoys through the arcs, the solution must satisfy the precedence relation induced on the set of arcs. Lines 13 to 17 intend to improve the solution with the help of pheromone functions.

The solution is built iteratively by adding or removing elements, being the nature of the elements dependent on the algorithm step. An element can be, for example, a trip of a vehicle or a demand node where to drop an aid kit. If *ζ* is a candidate element to be added, *phζ* and *grζ* are, respectively, the effective pheromone and the greedy function associated to it. As a result, the score of element *ζ* is calculated as stated in (**??**).

$$sc\_{\mathbb{S}} = \lambda ph\_{\mathbb{S}} + (1 - \lambda)gr\_{\mathbb{S}} \tag{13}$$

The reason why this expression is additive rather than multiplicative—as usual in ACO—is because, in this way, the gradual increase of the pheromone weight can be more accurately controlled.

#### **Algorithm 2** PTC Algorithm

The greedy function of a candidate element *ζ* is stated in (**??**), where *gζ* is the individual greedy function of objective *g*, *wg* is the weight of objective *g* in the objective function and *bg* is an upper bound for *gζ* .

$$\log r\_{\mathbb{S}} = \sum\_{\mathbb{S}^{\xi} \in G} w\_{\mathbb{S}} \left( \frac{b^{\mathbb{S}} - \mathfrak{e}\_{\xi}^{\mathbb{S}}}{b^{\mathbb{S}}} \right) \tag{14}$$

The greedy functions are used only where it is possible to define them in a reasonable way. When there are no greedy functions, the score is simply *scζ* = *phζ* if the decision regards which element should be added to the partial solution (it is convenient to select elements with more pheromone) or as stated in (**??**) if the decision regards which element should be removed from the partial solution (it is convenient to select elements with less pheromone).

$$\text{sc}\_{\%} = \frac{1}{p h\_{\%}} \tag{15}$$

In all steps of the PTC Algorithm, the elements are chosen according to the criterion established in ACS: if *q* ≤ *q*0, where *q* is a random variable uniformly distributed over [0, 1], the candidate element with the highest score is chosen; otherwise, the election is made according to the probabilities stated in (**??**), which are proportional to the scores.

$$pr\_{\tilde{\zeta}} = \frac{\text{sc}\_{\tilde{\zeta}}}{\sum\_{\tilde{\zeta}} \text{sc}\_{\tilde{\zeta}}} \tag{16}$$

The main steps of the PTC Algorithm are described in detail in what follows. **Design itineraries.** The routes to be followed by the vehicles are built iteratively, starting from an empty solution where all the vehicles are parked at their corresponding depots. At each iteration, the candidate elements to be added to the partial solution are the feasible trips of every single vehicle *j* through an available arc *k*. Therefore, each element is a pair *ζ* = (*j*, *k*). The effective pheromone associated to each element is defined in (**??**), where *c* is the vehicle class corresponding to vehicle *j*.

$$
\mu h\_{\zeta} = \left(1 - \frac{\upsilon t\_{c,k}}{c a r\_c}\right) \tau\_{c,k}^1 \tag{17}
$$

On the one hand, the effective pheromone is proportional to the first type pheromone, *τ*1 *c*,*k*, which depends on the vehicles of class *c* traversing arc *k* in the previously obtained solutions. In particular, the more vehicles of class *c* transit through this arc in the current best known solution, the more effective pheromone candidate element (*j*, *k*) has. On the other hand, the effective pheromone depends on the number of vehicles of class *c* traversing arc *k* in the current partial solution, so the greater *vtc*,*k* is, the less effective pheromone candidate element (*j*, *k*) has. If this gradual adjustment of the effective pheromone is not performed, it could lead to an excess of vehicles of a class through an arc for those pairs (*c*, *k*) with higher first type pheromone. This eventual excess of elements in a solution justifies the use of the effective pheromone and it also applies to other steps of the PTC Algorithm.

All the greedy functions used in the PTC Algorithm are defined in the same way as in the Elite Set Constructive Algorithm (ESCA) described in [**?** ]. For example, the individual greedy function of the time objective is the time that takes vehicle *j* to traverse arc *k* (**??**), while the individual greedy function of the reliability objective is the probability that arc *k* is not available (**??**).

$$\varrho\_{\zeta}^{t} = \frac{dist\_{k}}{\min\{vela\_{k'}, velv\_{j}\}} \tag{18}$$

$$
\varrho\_{\mathbb{Q}}^f = 1 - r\_k \tag{19}
$$

The score of each element is calculated as stated in (**??**). Once an element (*j*, *k*) is chosen, arc *k* is added to the route of vehicle *j* and the associated precedence relations are updated, in order to ensure that all vehicles travel together in convoys. The construction process is iterated until the routes cannot be continued any longer.

The augmenting paths from the source to the sink are determined iteratively, node by node, as follows:

• If the current node is the source, the candidate elements are the supply nodes with available aid. The effective pheromone of each supply node, *ζ* = *i*, is stated in (**??**), where *f lsc*,*i* is the current flow from the source to the supply node *i*, i.e., the amount of aid leaving the supply node *i* in the current partial flow solution. The effective pheromone depends positively on the third type pheromone (which is related to the amount of aid leaving supply node *i* in the previously obtained solutions) and negatively on the amount of aid leaving *i* in the current solution.

$$
\mu l \eta\_{\tilde{\varsigma}} = \left( 1 - \frac{f l\_{\kappa \jmath j}}{q a v\_i} \right) \tau\_i^3 \tag{20}
$$

No greedy function is defined in this case, so the score of each supply node *i* is equal to the corresponding effective pheromone.

•If the current node is a node *i*0 different from the source, a candidate element is a node *i* linked to *i*0 through an arc *k* = (*i*0, *i*) with positive residual capacity. The effective pheromone of each candidate node, *ζ* = *i*, is calculated as stated in (**??**). It depends positively on the amount of aid transported through arc *k* in the previously obtained solutions and negatively on the amount of aid transported through arc *k* in the current solution.

$$
\mu \text{pl}\_{\mathbb{V}} = \left(1 - \frac{f l\_{i0,i}}{q \text{g} l \text{obal}}\right) \text{r}\_{k}^{2} \tag{21}
$$

The individual greedy function of the time objective is the arriving time of the convoy associated to arc *k*, the individual greedy function of the cost objective is the average cost of transporting an aid kit by arc *k* among the vehicles of the associated convoy, etc.

If the current node *i*0 is a demand node, in addition to the other candidates mentioned earlier, the sink *sn* must be considered as a candidate element as well. The selection of the sink would finish the augmenting path. Its effective pheromone (**??**) depends positively on the amount of aid remaining at demand node *i*0 at the end of the operation in the previously obtained solutions and negatively on the remaining aid in the current solution.

$$
\mu h\_{sn} = \left(1 - \frac{f l\_{i\_0 sn}}{dem\_{i\_0}}\right) \tau\_{i\_0}^4 \tag{22}
$$

The individual greedy function of the equity objective associated to the sink measures how the deviation of the proportion of satisfied demand at node *i*0 from the ideal demand proportion would vary if the sink is selected; the individual greedy function of the priority objective is obtained in a similar way; and the individual greedy functions of the other objectives take the value 0, because finishing the augmenting path does not worsen those objectives.

Finally, the scores are established according to (**??**).

Each time a node is selected and labeled, the current node is updated. Once an augmenting path is determined, the flow is updated according to the Ford-Fulkerson algorithm. The process is iterated until all the planned amount of aid is distributed from the source to the sink or it is not possible to find an augmenting path, in which case the current solution would be discarded and the algorithm would be restarted at Step 1.

The augmenting paths from the source to the sink are determined in a similar way as in the Aid distribution step:

• If the current node is the source, the candidate elements are the improvable nodes. The effective pheromone of each improvable node, *ζ* = *i*, is stated in (**??**). The effective pheromone depends positively on the final stock at node *i* in the previously obtained solutions and negatively on the final stock at node *i* in the current solution.

$$
\mu h\_{\tilde{\varsigma}} = \left(1 - \frac{fs\_i}{dem\_i}\right) \tau\_i^4 \tag{23}
$$

The score of each improvable node *i* is equal to the corresponding effective pheromone.

• If the current node is a node *i*0 different from the source, a candidate element is a node *i* linked to *i*0 through an arc *k* = (*i*0, *i*) with positive residual capacity. The effective pheromone of each candidate node, *ζ* = *i*, is calculated as stated in (**??**). It depends positively on the amount of aid transported through arc *k* in the previously obtained solutions and negatively on the amount of aid transported through arc *k* in the current solution.

$$
\eta h\_{\tilde{\zeta}} = \left(1 - \frac{l t\_k - f l\_{i\_0, i}}{1 + q g
delta l}\right) \tau\_k^2 \tag{24}
$$

Since we are interested in removing needless flow, the score of each improvable node *i* is inversely proportional to the effective pheromone, as established in (**??**). Note that the unit added in the denominator of the fraction in (**??**) prevents the effective pheromone from being 0.

If the current node *i*0 is a demand node linked to the sink, the sink must be considered as a candidate element as well. Its effective pheromone (**??**) depends positively on the final stock at demand node *i*0 at the end of the operation in the previously obtained solutions and negatively on the final stock at *i*0 in the current solution.

$$
\eta h\_{sn} = \left(1 - \frac{sf\_{i\_0}}{1 + dem\_{i\_0}}\right) \pi\_{i\_0}^4 \tag{25}
$$

The score of the sink is established according to (**??**) as well.

Each time an augmenting path is obtained, the flow is updated according to the Ford-Fulkerson algorithm. The final stock at each demand node *i*, *s fi*, and the amount of aid transported through each arc *k*, *ltk*, are updated as well. The algorithm stops when the current flow cannot be increased further.

**Adjust supply nodes.** The purpose of this step is to eliminate as much incoming flow of aid as possible at supply nodes in which there is aid available at the end of the operation. Once again, a variation of the maximum flow problem is solved. The auxiliary transportation network includes all the nodes and the used arcs, which are now reversed because we try to return flow back from the supply nodes. All supply nodes are linked to the source and to the sink by dummy arcs. The capacity of each dummy arc leaving the source is the minimum between the amount of load arriving to the corresponding supply node and the final stock at it. The capacity of each dummy arc reaching the sink is the amount of aid available initially at the corresponding supply node. The capacity of each remaining arc is equal to the amount of aid transported through it (in opposite direction).

The calculation of the effective pheromones is as follows:

• If the current node is the source, the candidate elements are the supply nodes. The effective pheromone of each supply node, *ζ* = *i*, is stated in (**??**). The effective pheromone depends positively on the final stock at node *i* in the previously obtained solutions and negatively on the final stock at node *i* in the current solution.

$$\text{pl}\_{\mathbb{V}} = \left(1 - \frac{fs\_i}{1 + qav\_i}\right) \tau\_i^4 \tag{26}$$

Since we want part of the available (but unused) aid to leave the supply node and, therefore, the final stock decreases, the score of each supply node *i* is inversely proportional to the corresponding effective pheromone.

$$
\mu \text{pl}\_{\zeta} = \left(1 - \frac{v t\_{c,k}}{1 + c a r\_c} \right) \tau\_{c,k}^{1} \tag{27}
$$

Since we want to trim the routes, the scores are inversely proportional to the corresponding effective pheromones. Once an element (*j*, *k*) is chosen, arc *k* is removed from the route of vehicle *j* and the variable *vtc*,*k* is updated.

#### **5. Computational Study**

#### *5.1. Implementation Details and Calibration*

The algorithms introduced in the previous sections have been implemented in Fortran 95. The results presented in this section have been obtained by executing the algorithms on a personal computer with an Intel Core i5-2450M processor with 2.5 Ghz and 8Gb RAM running Windows 7.

Before applying the ACO Algorithm to our test cases, it is necessary to tune several parameters, namely the number of solutions obtained between each couple of global updates ( *m*), the probability of choosing the element with the highest score in any of the decision situations (*q*0), the local and global evaporation rates (*ρl* and *ρg*) and the two parameters which determine how to update the weight of the pheromone and the greedy functions in the constructive process ( *σ*1 and *σ*2). We randomly generated a set of five small test instances, with 10 to 17 nodes, 29 to 41 arcs and 20 to 30 vehicles, which allow us to do a large number of executions in order to obtain a set of parameters that makes the ACO Algorithm work well on realistic cases.

The calibration consists of two stages. In the first stage, the parameters were taken in three pairs, (*m*, *q*0), (*ρl*, *ρg*) and (*σ*1, *σ*2), and for each one, the ACO Algorithm was executed on all test instances, varying the values of the two chosen parameters in a range of between ten and twenty values, while the other four parameters were fixed to arbitrary values. With the results obtained in this first stage we established a new assignment of values.

In the second stage of the calibration, we took the parameters one by one and executed the ACO Algorithm 100 times on each instance, to obtain 2000 solutions in each execution, varying the value of the chosen parameter and fixing the values of the other five parameters according to the assignment established after the first stage. We collected the average of the best objective function value over the one hundred executions and over the five test instances, that is represented in the vertical axis of the figures included in this subsection, while the values tested for the parameter are represented in the horizontal axis. Each execution required an average time of 1.5 s of computation.

To calibrate the number of solutions obtained between each couple of global updates, we varied *m* = 0, 1, ... , 10. In the literature, *m* = 5 or *m* = 10 have been commonly used and our results did not contradict it, so we decided to set *m* = 5, which requires less computation time.

In the same way, the probability of choosing the element with the highest score is usually set to a high value in the literature; however, we found that *q*0 = 0.3 was the best value in our experiments, as shown in Figure **??** (note that the values higher than 0.5 have been rejected in the first stage of the calibration). This means that in 30% of the decision situations the candidate with the highest score will be selected directly, whereas in the remaining 70% the selection will be made at random, where the probability of each candidate is proportional to its score.

**Figure 1.** Calibration of the probability of choosing the element with the highest score.

Figure **??**a,b shows the results obtained for the evaporation rates. According to those, we set *ρl* = 0.02 and *ρg* = 0.04. Since both are small values, the pheromone trails will be very persistent, especially in the case of local updates, which are performed more frequently than the global ones.

Finally, we varied *σ*1 = 0, 0.05, ... , 0.5 and *σ*2 = 0, 1, ... , 10. The results led us to set *σ*1 = 0.2, meaning that at the start of the metaheuristic the weight of the greedy functions is 0.2, and *σ*2 = 5.

**Figure 2.** Calibration of the evaporation rates.

#### *5.2. Description of the Test Cases*

The proposed ACO Algorithm has been tested on two realistic test cases with different characteristics: Haiti earthquake 2010 and Niger famine 2005.

The first test case is based on the earthquake that devastated Port-au-Prince, Haiti's capital, and its metropolitan area in 2010. It was introduced in [**?** ], and it is very appropriate to illustrate the problem, since many of the characteristics that define it are present here. It develops mainly in an urban environment, but the effects of the disaster produced a great uncertainty about the transitability of many streets. Besides, aid operations were carried out in an environment of insecurity, with frequent attacks to convoys.

The logistic map, illustrated in Figure **??**, corresponds to the situation of the Port-au-Prince area fifteen days after the main earthquake. The network consists of 24 nodes and 84 arcs, corresponding to 42 two-way links. Nodes 1, 2 and 3 are supply nodes, with 60, 80 and 140 tons of humanitarian aid available, respectively. Nodes 10, 12, 13, 16, 17, 18, 20, 21 y 22 are demand nodes, with demands of 30, 40, 30, 30, 10, 30, 40, 20 y 20 tons of aid, respectively. Of these demand nodes, node 13 has special priority. The remaining 12 nodes are connection nodes and can be used to transship load.

The links have been classified by their quality (according to the thickness of the line: the thicker the line, the higher the speed) and by their reliability (according to the color scale shown in the figure). The arcs also have different ransack probabilities.

To perform the distribution operations there are 300 vehicles, classified in three different types according to their velocities, costs and capacities. There are 70 large vehicles, 95 medium size vehicles, and 135 small vehicles, which can carry 3, 2 and 1 tons of aid, respectively. The vehicles are located at the supply nodes and at the connection nodes 4, 6 and 18. The purpose of the operation is to deliver 150 tons of aid among the demand nodes.

The data for the Haiti test case come from various documents available in January 2010: a logistic map provided by the United Nations Office for the Coordination of Humanitarian Affairs (OCHA) [**?** ], a map reporting the satellite-identified IDP concentrations, road and bridge obstacles in central Port-au-Prince, Haiti, created by the United Nations Institute for Training and Research (UNITAR), and a map focused only on the Port-au-Prince road conditions created by the World Food Programme (WFP) Emergency Preparedness and Response [**?** ].

**Figure 3.** Transport network at Port-au-Prince.

The second test case, introduced in **?** ], is based on the 2005 Niger food crisis. The logistic network (see Figure **??**) consists of two supply nodes (Niamey and Gaya cities), two demand nodes (Agadez and Zinder), three connection nodes (Dosso, Tahoua and Maradi) and 11 two-way roads. There are a total of 376 vehicles available of three different types located at the supply and connection nodes to distribute 500 tons of humanitarian aid. The data for this test case come from the report of the real operation developed by the Spanish organization Fondo de Ayuda Humanitaria y de Emergencias de Farmamundi (FAHE) [**?** ] and public websites such as Google Maps.

**Figure 4.** Transport network at Niger.

Both test cases are explained in detail in [**?** ] and their data are available on the website: www.mat.ucm.es/humlog.

## *5.3. Results*

In this section we will compare and analyze the results obtained when applying the ACO Algorithm, proposed in this paper, and the GRASP Algorithm, introduced in [**?** ], to the Haiti and Niger test cases. Furthermore, we compare as well the best solutions obtained by both algorithms with the optimal values for these test cases when they are available. These optimal values were obtained by solving to optimality relaxed flow models [**? ?** ].

Firstly, we present the results obtained in the Haiti test case. Both algorithms—ACO and GRASP—have been executed ten times minimizing each objective independently and ten more times minimizing the aggregated function with equal weights for all attributes. The parameters of the algorithm were set according to the results of Section **??**. A fixed computation time of 1500 s has been set for each of the 140 executions—two algorithms, seven objectives per algorithm and ten runs per objective.

Table **??** shows the computational results obtained with the ACO Algorithm and with the GRASP Algorithm in the Haiti test case. For each algorithm, the objective values of the best and the worst solutions among the ten executions are presented, as well as the average and the standard deviation of the objective values of the ten solutions. In this test case the operation time is given in minutes, the cost in dollars, and the security and the reliability in expected lost tons. We can see that the ACO Algorithm provides a better best solution when minimizing time, cost and security attributes, although for the cost the average is worse than with GRASP. For the reliability attribute, both algorithms provide the same best result, but this best is obtained in all GRASP executions—standard deviation is equal to zero—and only in two ACO executions. Both algorithms provide optimal solutions for the priority attribute in all executions and ACO performs worse than GRASP for the equity and for the aggregated function. The ACO algorithm provides around 10,000 feasible solutions on average while the GRASP algorithm provides 25,000. Nevertheless, this substancial difference does not translate into a proportional improvement in performance.

**Table 4.** Results obtained with ten executions per objective of GRASP and ACO algorithms in Haiti earthquake 2010.

Table **??** shows the best values obtained among all tests performed during the experimentation—not only the ten executions of the previous table—in the Haiti test case, together with the optimal values for the attributes in which they are known. The optimal values for time and cost were obtained by using the static flow model described in [**?** ], which although has some different assumptions, can serve as a reference for these two attributes that are measured in a very similar way. The ACO Algorithm provides the assumed optimal value for the time attribute and stays very close to the optimum for the cost, substantially improving the results obtained with the GRASP Algorithm. For the priority and the reliability attributes, both algorithms provide the same values, and almost the same for the security, whereas for the equity attribute, the ACO Algorithm performs worse. For the aggregated function, the GRASP Algorithm also works better, probably because of the higher influence

of the equity attribute in the aggregated function, due to its tighter bound in comparison to the other attributes.

Figure **??** depicts the itineraries of the best solutions found for the time and the cost attributes. The numbers on the arcs represent the amounts of aid transported by the vehicles traversing them (a "0" on an arc indicates that an empty convoy is traversing the arc). In the best solution for the time (Figure **??**a), the three demand nodes farther from the supply nodes (18, 21 and 22) do not receive any aid and the distribution is mainly carried out by fast arcs, regardless of their dangerousness or reliability. The transportation is done by small and medium-sized vehicles, because they are faster than the large ones. When the cost attribute is optimized (Figure **??**b), only a few arcs are used and the convoys consist mainly of large vehicles, giving a better cost-capacity ratio. In this solution three demand nodes, including the priority node (13), are not visited.

**Table 5.** Best values obtained with several executions of GRASP and ACO algorithms in Haiti earthquake 2010.

(**a**) Optimize time

(**b**) Optimize cost

**Figure 5.** Best itineraries optimizing time and cost in Haiti earthquake 2010.

Figure **??** shows the solution obtained when minimizing ˜ *f*2 with equal weights. In this case, the ACO algorithm has been executed until 25,000 feasible solutions have been generated. This compromise solution provides reasonably good values for all attributes, as it can be observed in Table **??**, all demand nodes receive some aid and, in particular, the priority node receives 26 out of 30 tons demanded.

**Figure 6.** Best itineraries optimizing ˜ *f*2 in Haiti earthquake 2010.

**Table 6.** Best objective values optimizing *f*2 in Haiti earthquake 2010.

˜

Finally, we summarize the results obtained in the Niger test case. As in the Haiti test case, GRASP and ACO have been executed ten times minimizing each objective independently and ten times minimizing the aggregated function, running 1500 s per execution. Again, Table **??** shows the objective values of the best and the worst solutions among the ten executions for each attribute, together with the average and the standard deviation of the objective values of the ten solutions. The operation time is given in hours, the cost in euros, and the security and the reliability in expected lost tons. Both algorithms provide optimal solutions for the equity and priority objectives in the ten executions. The ACO Algorithm performs slightly better than GRASP for the security and significantly better for the remaining attributes—time, cost and reliability—and also for the aggregated function. The best value obtained with ACO for the cost is only 1.6% higher than the assumed optimal value of the cost attribute, which is EUR 63,240.0, obtained using the static flow model [**?** ]. ACO is significantly better for the cost, the time, the reliability and the aggregated function. Despite the higher speed of GRASP generating feasible solutions in this test case—9000 on average per execution with GRASP versus 6400 with ACO—ACO performance is much better than that of GRASP.

**Table 7.** Results obtained with ten executions per objective of GRASP and ACO algorithms in Niger Famine 2005.

It is remarkable that the number of feasible solutions generated using the same computation time is lower in the Niger test case than in the Haiti test case, even though the logistic network is much smaller in the former. This is because there are other factors that influence the computation time as, for example, the structure of the network, the number and the location of the vehicles, the relation between the amount to be distributed and the total capacity of the vehicles, etc.

The proposed compromise solution for the Niger test case is shown in Table **??**. All attribute values remain reasonably close to the best values obtained by optimizing each objective individually, except the value for the priority objective. This is partially due to the high confrontation between the equity and priority attributes. For optimizing the priority attribute, the node with the highest priority level (Zinder) should receive all the demanded aid, which would lead to an inequitable distribution.

**Table 8.** Best objective values optimizing ˜ *f*2 in Niger famine 2005.

In summary, the results obtained in both test cases show that the ACO Algorithm is able to obtain good solutions for all attributes when optimizing each of them individually, together with good compromise solutions as well. In particular, it provides the best known solutions for several attributes, such as the operation time in the Haiti test case or the reliability of the itineraries in the Niger test case. In addition, the computational effort required by the ACO algorithm is quite reasonable, making it suitable to solve realistic cases.

#### *5.4. Managerial Insights*

#### **6. Concluding Remarks**

This paper addresses a humanitarian last mile distribution problem in response to a disaster. The distribution must be done under unsecure conditions, so the vehicles must travel in convoys to try to prevent assaults. Furthermore, there is uncertainty about the state of the roads, that may have been damaged by the disaster. The proposed model is multimodal, multidepot, allows split deliveries and multiple trips and provides the individual itineraries of the vehicles. Besides, the model considers six performance criteria—time, cost, equity, priority, security and reliability—and is solved by a compromise programming approach.

The main contribution of this paper is the proposed methodology to solve the model. We develop an algorithm based on Ant Colony Optimization that considers several types of ants and pheromones. The methodology introduces the concept of effective pheromones to balance the elements of the provided solutions and to diversify the set of obtainable solutions. This methodology is especially appropriate to approach routing problems in which the vehicles are required to travel in convoys, but some of its most novel elements, as it is the case of the effective pheromones, could also be successfully applied to many other complex distribution and routing problems.

The algorithm is tested on two realistic case studies of public domain, and its performance is compared to that of the GRASP Algorithm, the only existing solution method in the literature. The ACO Algorithm performs well for all the objective functions tested in our experiments and provides, for both test cases considered, the best solutions obtained to date for several attributes. In comparison with the available solution method in the literature—the GRASP Algorithm—in the Haiti test case, both algorithms have their own strengths, but in the Niger test case, ACO clearly outperforms GRASP. The computation time required to obtain good solutions is reasonable in all cases, proving that the ACO Algorithm can be perfectly applied in an emergency context. Furthermore, GRASP and ACO algorithms actually complement each other and could even be used in combination in order to improve the quality of the final solutions provided to the decision makers.

One interesting future research line to continue this work comprises the development of a stochastic model in which the possibility of suffering assaults or finding blocked roads could be represented by random variables. Additionally, the design of simulation models to replicate real situations and test different solution strategies could also be worth considering.

**Author Contributions:** Conceptualization, J.M.F., M.T.O. and G.T.; Methodology, J.M.F., M.T.O. and G.T.; Software, J.M.F.; Validation, J.M.F., M.T.O. and G.T.; Formal analysis, J.M.F., M.T.O. and G.T.; Investigation, J.M.F., M.T.O. and G.T.; Data curation, J.M.F., M.T.O. and G.T.; Writing—original draft preparation, J.M.F.; Writing—review and editing, M.T.O. and G.T.; Visualization, J.M.F., M.T.O. and G.T.; Project administration, M.T.O. and G.T. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the European Commission grant H2020 MSCA-RISE 691161 (GEO-SAFE), the Government of Spain grant MTM2015-65803-R and the Government of Madrid, grant S2013/ICE-2845. The APC was funded by the European Commission grant H2020 MSCA-RISE 691161 (GEO-SAFE).

**Conflicts of Interest:** The authors declare no conflict of interest.


**6. Conclusions**

In this paper, we propose a novel electricity theft detection scheme based on TextCNN. We innovatively formulated the electricity data into two-dimensional time-series in order to capture the intraday and daily correlations of electricity consumption data. Then, we discussed the relationship between DNN, CNN and TextCNN, and explained why TextCNN is the most suitable classifier for our purposes, considering both the e fficiency and e ffectiveness. Additionally, in order to balance the electricity consumption dataset, we proposed a data augmentation method. We conducted extensive experiments on di fferent realistic datasets to prove the e ffectiveness of the proposed scheme, including the residential and industrial datasets from a province in China and the public Irish residential dataset. The experimental results show that the proposed method outperforms other methods, such as LR, SVM, DNN and 1D CNN. At the same time, we analyzed the importance and e ffectiveness of data augmentation.

**Author Contributions:** Conceptualization and methodology, X.F. and H.H.; software, Z.L. and Y.Y.; validation, W.G., H.Q. and H.F.; writing, H.H., Z.L., Y.Y. and C.Y.; supervision and project administration, X.F., H.Q., H.F. and Y.D. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by Project Supported by the China Southern Power Grid Corporation, grant number GDKJXM20185800.

**Conflicts of Interest:** The authors declare no conflict of interest.


**7. Materials and methods**

452 patients who underwent aortic valve surgery were retrospectively analyzed by being divided into the following groups: single and multiple aortic valve replacements, aortic

Neurological Complications in Aortic Valve Surgery and Rehabilitation Treatment Used 197

diagnosed or treated with medication, a diet or exercise by a doctor

been diagnosed or treated with medication or chest physiotherapy.

no response to verbal or visual stimuli.

motor functions, seizures and psychiatric disorders

by a doctor with medication, a diet or exercise.

consumers of cigarettes (>2500 cigarettes/person/ year).

doctor with medication, a diet or exercise.

ray, dyspnea, or ventricular gallop.

valves that were operated on.

with a cardiovascular emergency).

immobilization.

reversible.

alcohol.

amputation.

severe <30%.

kinesitherapy.

valve replacement plus coronary artery bypass grafting and aortic valve replacement plus tube graft due to root aneurysm between January 2008 and November 2010.

Fig. 7. Reimplantation using the inclusion (David) technique

An NC was defined as the occurrence of a cerebrovascular accident (STROKE) (ischemic or hemorrhagic), transient ischemic attack (TIA), spinal cord injury, peripheral neuropathy, seizure, stupor, coma, polyneuropathy of critically ill patient, dementia, acute delirium or encephalopathy.

A comparative analysis was carried out on the incidence of NCs according to a series of preoperative and postoperative clinical variables: arterial hypertension, diabetes mellitus, dyslipidemia, chronic obstructive pulmonary disease, heart failure, renal failure, prior STROKE, smoker, drinker, calcified valves, endocarditis, peripheral arterial disease, previous revascularization, aortic atherosclerosis, acute myocardial infarction within 3 months prior to surgery, left ventricle ejection fraction, aortic clamping time and CPB time, postoperative arrhythmia, number of transfusions required, whether resuscitation techniques were required, drugs the patients were taking, surgical priority (scheduled, priority and urgent) and type of surgery performed.

An attempt was made to determine which of these variables was more important statistically in relation to the others.

### **7.1 Definitions used**

valve replacement plus coronary artery bypass grafting and aortic valve replacement plus

An NC was defined as the occurrence of a cerebrovascular accident (STROKE) (ischemic or hemorrhagic), transient ischemic attack (TIA), spinal cord injury, peripheral neuropathy, seizure, stupor, coma, polyneuropathy of critically ill patient, dementia, acute delirium or

A comparative analysis was carried out on the incidence of NCs according to a series of preoperative and postoperative clinical variables: arterial hypertension, diabetes mellitus, dyslipidemia, chronic obstructive pulmonary disease, heart failure, renal failure, prior STROKE, smoker, drinker, calcified valves, endocarditis, peripheral arterial disease, previous revascularization, aortic atherosclerosis, acute myocardial infarction within 3 months prior to surgery, left ventricle ejection fraction, aortic clamping time and CPB time, postoperative arrhythmia, number of transfusions required, whether resuscitation techniques were required, drugs the patients were taking, surgical priority (scheduled,

An attempt was made to determine which of these variables was more important

haemorrhage persisting longer than 24 hours or which leaves side-effects.

tube graft due to root aneurysm between January 2008 and November 2010.

Fig. 7. Reimplantation using the inclusion (David) technique

priority and urgent) and type of surgery performed.

statistically in relation to the others.

motor, speech or sight functions.

hours and with no side-effects.

without loss of consciousness.

**7.1 Definitions used**

mobility

mutism.

encephalopathy.

Neurological Complications in Aortic Valve Surgery and Rehabilitation Treatment Used 199

actinomyces odontolyticus, Granulicatella, pseudomonas aeruginosa, staphylococcus epidermidis, S. pneumoniae, 2 cases reported as Gram strains and two other cases as

**Drugs Anti AHT Sintrom Diuretics Hypolipemiant Antiplatelet OAD/Insulin**  % 67 18.1 34.3 42.5 29.9 17.5

NC STROKE TIA SCI PN Seizure Stupor Coma CIP ACS Encephalopathy % 1.5 0.67 0.4 0.2 1.3 1.3 0.4 1.1 6.6 5 Table 3. Percentages of NCs (neurological complications): TIA (transient ischemic attack), SCI Spinal Cord Injury, PN (peripheral neuropathy), CIP (critical illness polyneuropathy),

Table 3 shows the percentages of different NCs in Ao valve surgery. The rate of strokes was 1.5%, corresponding to 7 cases. Of these 7 cases, 2 came from the 71 patients who underwent valve replacement + CABG AO (2.8%), 3 were from the 227 patients that underwent single Ao valve replacement surgery (1.3%), one case was from the 72 patients who underwent multiple valve replacement (1.4%) and finally there was one case from the 12 patients who underwent a 2nd Ao valve replacement (8.3%). Table 4 shows the three most common NCs that have been seen in this study according to the different surgical techniques and the

Of the 23 cases of encephalopathy, 17.4% had a metabolic cause, in 56.5% the cause was hypoxic-ischemic and 26% had various causes. Of the strokes, 57.1% were right hemisphere and 42.9% were on the left. There was only one case in which the stroke was haemorrhagic. Of the ischemia (6 cases - 85.7%), 4 were of cardioembolic origin (66.6%), 2 border territory

Replacement + CBPG 71 2(2.8) 8(11.3) 5(7) 128.5 105.7 1 Ao valve replacement 227 3(1.3) 12(5.3) 6(2.6) 96.5 75.4 Multiple vv replacement 72 1(1,4) 4(5.5) 5(6,9) 151.2 126.8 Ao root repl. (Bentall) 58 0 5(8.6) 5(8.6) 148.2 118.2 Ao root repl. (David) 10 0 0 0 133.6 129.2 2nd Ao valve repl. 12 1(8,3) 3(25) 2(16.6) 166.2 127.3 Endoscopic valve repl. 2 0 0 0 - - Table 4. Most common NC percentages according to the different techniques used. Ao (aorta), vv (valve), repl (replacement), n: total number of patients, ACS (acute confusional syndrome). T. CPB (cardiopulmonary bypass time measured in minutes), T. Clamp (aortic clamping

There are only 2 cases described of spinal cord injury. One was a side effect of Ao valve replacement surgery and the other came after surgery for aortic arch replacement using the Bentall technique. The two cases were side effects to a spinal cord ischemia. One was a ASIA

**CPB**

**T. Clamp**

Table 2. Drugs used by patients with operation on the Ao valve

average times for CPB and Ao clamping measured in minutes.

(33.3%) (Man in the Barrel Syndrome) and one was lacunar (16.6%).

time measured in minutes). CBPG (coronary artery bypass grafting)

**Type of surgery n Stroke ACS Encephalopathy T.**

undetermined.

ACS (confusional syndrome)


**5. Conclusion**

*Animal Models in Medicine and Biology*

consecutive weekly measurements.

challenging to perform.

vessel lumen (stenosis) and reduction of the blood flow.

Mechanical endothelium injury is a standard strategy to induce hyperplasia in various animal models such as mice, rabbits or minipigs [9, 10, 17]. In the present study, the endothelial injury was induced by two different approaches in rabbits: first by overexpansion of a BMS in the right iliac artery and second with an inflated balloon in the contralateral left iliac artery. Stenting and balloon injury resulted in a respective permanent and transient overexpansion of the vessel inducing a well-recognised vascular response [11]. An acute inflammation rapidly follows the induced endothelium injury or denudation and peaks after a few days; then, the inflammation temporally declines during the resolution phase. Acute inflammation is accompanied by the proliferation of smooth muscle cells that results in hyperplasia. Consequent thickening of the wall (hyperplasia) may induce a narrowing of the

In the present study, we showed that the rabbit iliac arteries and the distal abdominal aorta could be successfully monitored using high-resolution ultrasound for longitudinal and non-invasive investigation. The quality of the images acquired allowed rigorous measurements of the wall and vessel sizes. It is important to note that for some measurements, the possible presence of oedema could impair the quality of the image acquisition and may explain variations of the parameters in

Besides, the metallic structures of the BMS were visualised as shadows. Although the presence of the metallic stent did not impair the wall thickness and vessel diameters' quantification, the blood flow velocity measurements were often

lumen size followed by a reduction suggesting a weakening of the artery. Notably, the overexpansion of the stent affected the wall thickness that increased rapidly and remained elevated. Histology analysis revealed the formation of neointimal hyperplasia. In contrast, balloon injury results in a transient wall thickening recorded 2 weeks post-intervention with a successive return to the initial dimensions. Transient inflammatory response to the balloon injury may explain the wall. Accordingly, Welt et al. [12] reported a transient inflammation following a balloon injury model as compared with a sustained accumulation of inflammatory

cells such as monocytes in stented iliac arteries of rabbits up to 14 days.

Furthermore, Virmani et al. [13] reviewed the temporal vascular response to BMS implanted in the rabbit iliac artery. Between 7 and 14 days, the intimal thickness increased due to inflammation and cell proliferation and then peaked at 1 month with a return of the cell proliferation to the basal level. Further shrinkage has been reported from 3 to 6 months due to the extracellular matrix remodelling. In agreement, we also report the shrinkage of the wall thickness observed at day

As far as the evaluation of the blood flow is concerned, the VTI and the calculated changes in flow revealed notable individual variations. Nevertheless, the longitudinal evaluation showed that in all animals, the VTI decreased overtime in the left and right segments but remained stable in the proximal segments. Meanwhile, the calculated blood flow was significantly reduced at the site of the balloon injury

and the distal segment but maintained unchanged in the stented segment.

Comparing the acute balloon injury with the chronic injury associated with stent implantation, we reported that the response of the vessel differed with the type of intervention. The primary function of the BMS stent is to provide mechanical support. As expected, the stent allowed the maintenance of increased vessel diameter in all animals. In contrast, balloon inflation induced a short-term increase in the

**4. Discussion**

**190**

28 in the left artery.

The non-invasive, real-time imaging of the rabbit iliac arteries and the distal abdominal aorta for the quantification of lumen diameter and wall thickness using high-resolution ultrasound permit the monitoring of the progression of the wall and vessel following balloon angioplasty and endovascular stent implantation. Conveniently combined with the blood flow analysis, this methodological approach would be essential to evaluate novel therapeutic approaches to prevent hyperplasia.


*3.1. Principle of Electricity Measurement*

Most industrial customers equipped with three-phase smart meters. According to the different wiring modes, it can be divided into two types: three-phase four-wire and three-phase three-wire. In this paper, we primary introduce our approach with three-phase four-wire as an example. For single phase scenario, the electricity energy is calculated by following equation:

$$E = \sum\_{n=0}^{N-1} P\_n \cdot \Delta t \tag{1}$$

where, Δ*t* is the cycle of computation, *E* is the active electricity energy in a certain time period *N* · Δ*t*, *Pn* is the average active power at time n. It can be further calculated by voltage and current:

$$P\_n = \mathcal{U}\_n \cdot I\_n \cdot \cos \phi\_n \tag{2}$$

where, *Un* and *In* are average voltage and current at time *n*, *φn* is the phase difference between voltage and current at the same time stamp and *cosφn* is called power factor. Further, the equation of active energy are rewritten by:

$$E = \sum\_{n=0}^{N-1} \left[ \mathcal{U}\_n \cdot I\_n \cdot \cos \phi\_n \right] \cdot \Delta t \tag{3}$$

For three-phase four-wire situation, the total active electricity energy is calculated by:

$$E\_{total} = \sum\_{n=0}^{N-1} P\_{total} \cdot \Delta t = \sum\_{n=0}^{N-1} \left[ P\_A^\text{\textquotedbl} + P\_B^\text{\textquotedbl} + P\_C^\text{\textquotedbl} \right] \cdot \Delta t \tag{4}$$

Generally, electricity energy measured by voltage and current, and the relationship between them is very important to NTL detection.


The effect of scattering feedback of individual scattering part on random lasing is demonstrated by measuring the emission spectra of the random lasers under different pump power densities in Figure 2. A blue random laser was obtained by coupling the C1 solution and PDLC film, and named as RLB. The emission spectra of RLB at different pump power densities when the external electric field is removed from the PDLC film in OFF state are shown in Figure 2a. Only broadband fluorescence spectra can be observed as the pump power density is less than 0.30 MW cm−2. However, when the pump power density is greater than 0.30 MW cm−2, many narrow linewidth peaks can be observed around 457.5 nm. The results show that the emission light from the C1 solution is randomly scattered by the PDLC film and the coherent resonance is established in the coupling system. Meanwhile, the variations of integral intensity of the emission spectra with pump power density are shown in Figure 2c. When the pump power density is greater than 0.30 MW cm−2, the integral intensity increases sharply with increasing pump power density, which indicates that the threshold of the random lasing is about 0.30 MW cm−2 [27,28]. The emission spectra of RLB at different pump power densities when the electric field is applied to the PDLC film in ON state are shown in Figure 2b. Coherent feedback random lasing can also be obtained when the pump power density is greater than 0.20 MW cm−2, and the intensity is much higher than that of PDLC film in OFF state. The corresponding threshold of RLB is about 0.20 MW cm−2 (Figure 2c), which is smaller than that of the PDLC film in OFF state. The increase in intensity and decrease in threshold are due to the weak scattering feedback of the PDLC film in ON state.

**Figure 2.** Performance characterization of the random lasers. (**a**,**b**) The emission spectra of RLB at different pump power densities when the electric field is removed (**a**) and applied (**b**). (**c**) The integral intensity varies with pump power density for RLB. (**d**,**e**) The emission spectra of RLG at different pump power densities when the electric field is removed (**d**) and applied (**e**). (**f**) The integral intensity varies with pump power density for RLG. (**g**,**h**) The emission spectra of RLR at different pump power densities when the electric field is removed (**g**) and applied (**h**). (**i**) The integral intensity varies with pump power density for RLR. The corresponding spectral integral range is from 440 nm to 470 nm for RLB, from 510 nm to 570 nm for RLG and from 625 nm to 700 nm for RLR.

Similarly, the effect of scattering feedback of individual scattering part on the emission spectra of green random laser was also studied. A green random laser named as RLG was obtained by coupling the C153 solution and PDLC film. When the electric field of PDLC film is removed (OFF), coherent feedback random lasing with a central wavelength of about 542.7 nm and a threshold of about 1.15 MW cm−2 can be observed (Figure 2d,f). Meanwhile, when the electric field is applied to the PDLC film (ON), coherent feedback random lasing with a threshold of about 0.87 MW cm−2 can be obtained (Figure 2e–f). The regulation effect of scattering feedback on emission spectra of RLG is similar to that of RLB. In addition, the effect of scattering feedback of individual scattering part on the emission spectra of red random laser is demonstrated. A red random laser named as RLR was achieved by coupling the DCJTB solution and PDLC film. When the electric field of

PDLC film is removed (OFF), coherent feedback random lasing with a threshold of about 1.07 MW cm−2 can be observed around 662.9 nm in Figure 2g,i. When the PDLC film is applied electric field and is in ON state, coherent feedback random lasing with a threshold of about 0.87 MW cm−2 is obtained (Figure 2h,i). The experimental results show that the random lasing can always be achieved based on the separate coupling configuration with pure gain region and PDLC scattering part. Moreover, the threshold of the blue, green and red random lasing can be effectively modified through regulating the scattering feedback of individual scattering part.

The effect of scattering feedback of individual scattering part on the intensity of random lasing is demonstrated by continuously pumping 100 times at a fixed pump power density in Figure 3. The integral intensities of the emission spectra are calculated when the PDLC film is with external electric field (ON, black sphere) and without electric field (OFF, blue sphere) at a pump power density of 0.62 MW cm−2 in Figure 3a. The emission intensity of RLB in ON state is much higher than that in OFF state, and the average integral intensity of 100 spectra in ON state is about 109 times of that in OFF state. Similarly, RLG and RLR are continuously pumped for 100 times under the pump power densities of 1.39 MW cm−2 and 1.46 MW cm−2, respectively. Their integral intensities of the emission spectra are demonstrated through when the electric field is applied (ON, black sphere) and removed (OFF, blue sphere) in Figure 3b,c. The scattering feedback of individual scattering part has a similar regulation effect on the emission intensity of RLG and RLR as RLB. The average integral intensities of RLG and RLR in ON state are about 79 and 48 times higher than those in OFF state, respectively. Therefore, the intensity of the random lasing can be dynamically manipulated by applying or removing an external electric field to the PDLC film.

**Figure 3.** Regulation effect of scattering feedback on the intensity of random lasing. (**a**) The integral intensity evolution of the emission spectra when the RLB is in ON state (black) or OFF state (blue) under a pump power density of 0.62 MW cm−2. (**b**) The integral intensity evolution of the emission spectra when the RLG is in ON state (black) or OFF state (blue) at 1.39 MW cm−2. (**c**) The integral intensity evolution of the emission spectra when the RLR is in ON state (black) or OFF state (blue) at 1.46 MW cm−2. The corresponding integral range is (440, 470) nm for RLB, (510, 570) nm for RLG and (625, 700) nm for RLR.

The effect of scattering feedback of individual scattering part on the spatial coherence of random lasing is also demonstrated by calculating the visibility of interference fringe. The Young's double-slit interference optical path is used to study the spatial coherence of the obtained random lasing. Two slits with a width of 50 μm are separated by 160 μm, and a CCD is positioned behind the slit plane to measure the far-field interference pattern. The visibility (*V*) of the interference fringe in the central region is calculated by *V* = (*I*max − *I*min)/(*I*max + *I*min), where *I*max and *I*min are the local maximum and minimum values of the interference intensity near the observation point, respectively [21]. For example, the optical photo and the corresponding intensity distribution of the interference fringe generated by the blue random lasing at the pump power density of 0.35 MW cm−2 are shown in Figure 4a,b. The *V* is about 0.35 for the random lasing related to ON state, and about 0.16 for the random lasing under OFF state. The results show that the spatial coherence of blue random lasing in OFF state is much lower than that in ON state. Similarly, the optical photos and intensity distribution of the interference fringes obtained by the

green random lasing passing through the double-slit interference optical path at the pump power density of 1.32 MW cm−2 are shown in Figure 4c,d. The *V* are about 0.34 (ON) and 0.16 (OFF), respectively. In addition, the spatial coherence of the red random lasing is also characterized at 1.32 MW cm−2. The *V* of the interference fringes of red random lasing are 0.31 (ON, Figure 4e) and 0.14 (OFF, Figure 4f), respectively. The experimental results show that the visibility of interference fringes of the obtained random lasing in ON state is greater than that in OFF state. Moreover, the spatial coherence of random lasing emitting from another side are also measured, and the results prove that there exists the weak coupling between the gain region and the scattering region of the PDLC film supplying tunable backscattering. Therefore, the spatial coherence of the obtained random lasing can be dynamically manipulated by electrically regulating the scattering coefficient of the PDLC scattering part.

**Figure 4.** Regulation effect of scattering feedback on the spatial coherence of random lasing. (**a**,**b**) The intensity distribution and visibility of the interference fringes obtained by the double slit interference of RLB under ON state (**a**) or OFF state (**b**). (**c**,**d**) The intensity distribution and visibility of the interference fringes for the RLG under ON state (**c**) or OFF state (**d**). (**e**,**f**) The intensity distribution and visibility of the interference fringes for the RLR under ON state (**e**) or OFF state (**f**). Illustrations are the optical photos of the interference fringes, respectively.

#### **4. Conclusions**

In summary, a random lasing with tunable spatial coherence is designed based on the separate coupling between gain medium and scattering part. The proposed random laser has the advantages of flexible controllability, easy fabrication and recyclability. First, when the PDLC film is separately coupled with different color gain dyes, random lasing with different wavelengths is realized under the excitation of pulsed laser. Then, the scattering intensity of PDLC film is manipulated by applying or removing electric field, and the threshold and intensity of random lasing are regulated. Finally, the spatial coherence of the random lasing can also be dynamically regulated. In short, the intensity and spatial coherence of random lasing can be manipulated by controlling the electric field. The proposed random lasing is expected to be used in interference measurement, security, optical imaging and display.

**Author Contributions:** Conceptualization, Z.W. and W.G.; Investigation, Y.B. and J.Z.; Methodology, Y.B. and D.L.; Project administration, Z.W.; Supervision, Z.W.; Validation, Z.W.; Visualization, Y.B. and Z.W.; Writing, Y.B., H.Y. and Z.W. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by National Natural Science Foundation of China (grant Nos. 92150109, 11574033 and 61975018).

**Acknowledgments:** The authors acknowledge support from the National Natural Science Foundation of China and Beijing Higher Education Young Elite Teacher Project.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

MDPI St. Alban-Anlage 66 4052 Basel Switzerland Tel. +41 61 683 77 34 Fax +41 61 302 89 18 www.mdpi.com

*Crystals* Editorial Office E-mail<EMAIL_ADDRESS>www.mdpi.com/journal/crystals

Academic Open Access Publishing

www.mdpi.com ISBN 978-3-0365-7673-2