clear all;
clc;
format short g;
global D;
global q;
global q1;
global ss;
global E;
global L;
global ELL;
% test=xlsread('test.xlsx',1);
%
% position=test(:,2:3);
% data=[
% 100 100 0 0 0 40
% 60 100 10 0.2 0 3
% 20 110 15 0.3 0 4
% 160 150 21 0.3 1 6
% 160 110 16 0.4 0 4
% 30 20 25 0.3 1 4
% 100 60 22 0.2 0 6
% 100 170 15 0.1 0 5
% 30 10 12 0.4 0 6
% 60 50 15 0.3 1 7
% 160 160 20 0.2 0 6
% 50 140 23 0.3 0 6];
position=[18.70 15.29
16.47 8.45
20.07 10.14
19.39 13.37
25.27 14.24
22.00 10.04
25.47 17.02
15.79 15.10
16.60 12.38
14.05 18.12
17.53 17.38
23.52 13.45
19.41 18.13
22.11 12.51
11.25 11.04
14.17 9.76
24.00 19.89
12.21 14.50];
% Matrice de distance entre le Centre de distribution et chaque point de demande D
D=squareform(pdist(position(:,1:2),'euclidean'));
% Coût unitaire de déplacement du véhicule : h
% Utilisation des coûts fixes par véhicule ：R.
% Vitesse de déplacement du véhicule v;
% Charge maximale du véhicule Qmax.
h=1;R=10;v=1;Qmax=200;epsilon=0.001;
R0=2;R1=0.8;
% Matrice de temps de déplacement du véhicule entre le Centre de distribution et chaque point de demande T
T=D/v;
% Itérations maximales
Iter=30;
% Nombre de gouttelettes N
N=size(position,1);
% Paramètres de mise à jour de la vitesse ：
av=1;bv=0.01;cv=1;
% Paramètres de mise à jour de la quantité de sol ：
as=1;bs=0.01;cs=1;
% Coefficient de renouvellement local des sols
alpha=0.9;
% Coefficient global de renouvellement des sols
beta=0.9;
% Quantité initiale de sol entre deux points Initsoil
Initsoil=1000;
% Matrice initiale des quantités de sol W ???
W=ones(N)*1000;
% Vitesse initiale par goutte InitVel;
% InitVel=100;
InitVel=randperm(100,1);
% Valeur de la fonction objective optimale globale
TotalZ=1000000;
% Chemin Global optimal
TotalRoute=[];
% Procédure principale
t=0;
WaterDrop(1,N)=struct('SW',[],... % Matrice initiale de quantité de sol pour les gouttelettes d'eau
'Source',[],... % Le point de départ des gouttelettes d'eau est 1; Centre de distribution
'Target',[],... %
'VisitNode',[],... % Séquence de points visités (Chemin d'accès)
'UnvisitNode',[],... % Collection de points non visités
'Vel',[],... % Vitesse initiale des gouttelettes
'Soil',[],... % Quantité initiale de sol transportée par les gouttelettes
'Q',[],... % Charge au départ des gouttelettes
'S',[],... % Les gouttelettes arrivent. source(k) Heure du point
'ZZ',[],... % Valeur initiale de la fonction objective pour les gouttelettes
'FK',[]); % Desource(k) L'ensemble des points de demande suivants où vous pouvez commencer FV
% q=test(:,4); % Demande par point de demande
% q1=TestData(:,5);
q=[0 6 5 11 6 3 8 5 6 4 5 7 6 10 9 4 7 8 ]';
% q1=[0 3.6 2 4.6 3.6 2.4 4.8 3 3.6 2.4 3 4.2 3.6 4 5.4 2.6 2.4 3]';
% ss=test(:,7); % Temps de service par point de demande
ss=[0 1.8 1.0 2.3 1.8 1.2 2.4 1.5 1.8 1.2 1.5 2.1 1.8 2.0 2.7 1.3 1.2 1.5]';
E=[0 5.0 4 1 2.0 5 2 1 3 1 2 2.0 2 3 2 3 2 1]';
L=[40 20 15 20 20 15 18 24 27 20 16 20 10 25 28 24.0 20 23]';
%
% E=test(:,5); % Fenêtre de temps inférieure par point de demande
% L=test(:,6); % Limite supérieure de la fenêtre de temps par point de demande
ELL=L-E;
% EP=test(:,8);
% LP=test(:,9);
I=1:N;
while t < Iter
clc
fprintf('No%d Sous - évolution \n',t+1) ;
% Définir la variable dynamique initiale ：
% Valeur de la fonction objective optimale pour cette itération
Z=1000000;
% Chemin optimal pour cette itération
Route=[];
for k=1:N
WaterDrop(k).SW=W;
WaterDrop(k).Source=1;
WaterDrop(k).VisitNode=[1];
WaterDrop(k).UnvisitNode=I;
WaterDrop(k).Vel= InitVel;
% Quantité initiale de sol transportée par les gouttelettes
WaterDrop(k).Soil=0;
% Charge au départ des gouttelettes
WaterDrop(k).Q(WaterDrop(k).Source)=0;
% Les gouttelettes arrivent. source(k) Heure du point ;
WaterDrop(k).S(WaterDrop(k).Source)=0;
% Nok Valeur initiale de la fonction objective pour les gouttelettes
WaterDrop(k).ZZ=0;
end
for k=1:N
while WaterDrop(k).Source~=1 || ~isequal(WaterDrop(k).UnvisitNode,[1])
WaterDrop(k).FK=[];
if WaterDrop(k).Source==1
alt=2;
else
alt=1;
end
for i=alt:length(WaterDrop(k).UnvisitNode)
WaterDrop(k).Q(WaterDrop(k).UnvisitNode(i))=WaterDrop(k).Q(WaterDrop(k).Source)+q(WaterDrop(k).UnvisitNode(i));
% Jugementi La charge au point est - elle inférieure à la charge maximale du véhicule et est - elle atteinte? i Si l'heure du point est i Dans la fenêtre de temps requise par le point
if WaterDrop(k).Q(WaterDrop(k).UnvisitNode(i))<Qmax;
WaterDrop(k).FK=[WaterDrop(k).FK,WaterDrop(k).UnvisitNode(i)];
end
end
% for i=alt:N
% if ismember(i,WaterDrop(k).UnvisitNode)
% % if sum(ismember(WaterDrop(k).UnvisitNode,i))
% WaterDrop(k).Q(i)=WaterDrop(k).Q(WaterDrop(k).Source)+q(i);
% WaterDrop(k).S(i)=WaterDrop(k).S(WaterDrop(k).Source)+ss(WaterDrop(k).Source)+T(WaterDrop(k).Source,i);
% % Jugementi La charge au point est - elle inférieure à la charge maximale du véhicule et est - elle atteinte? i Si l'heure du point est i Dans la fenêtre de temps requise par le point
% if WaterDrop(k).Q(i)<Qmax && WaterDrop(k).S(i)>=E(i) && WaterDrop(k).S(i)<=L(i)
% WaterDrop(k).FK=[WaterDrop(k).FK,i];
% end
% end
% end
if isempty(WaterDrop(k).FK)
WaterDrop(k).Target=1;
else
% Calculer la probabilité d'atteindre le prochain point de demande disponible
% De source(k)Est - ce que la quantité minimale de terre sur le chemin vers le point de demande suivant qui peut être desservi est inférieure à0
Minsoil=0;
for u=1:length(WaterDrop(k).FK)
if WaterDrop(k).SW(WaterDrop(k).Source,WaterDrop(k).FK(u))<Minsoil
Minsoil=WaterDrop(k).SW(WaterDrop(k).Source,WaterDrop(k).FK(u));
end
end
% Trouver la fonction suivante correspondant au point de demande de service fLa somme de（ Peut ajouter des améliorations , Ajuster la probabilité de sélection des noeuds ）
SumF=0;
for u=1:length(WaterDrop(k).FK)
g(WaterDrop(k).Source,WaterDrop(k).FK(u))=WaterDrop(k).SW(WaterDrop(k).Source,WaterDrop(k).FK(u))-Minsoil;
f(WaterDrop(k).Source,WaterDrop(k).FK(u))=1/(0.01+g(WaterDrop(k).Source,WaterDrop(k).FK(u)));
SumF=SumF+f(WaterDrop(k).Source,WaterDrop(k).FK(u));
end
% Trouver la probabilité correspondant au point de demande de service
P=[];
U=[];
for u=1:length(WaterDrop(k).FK)
P=[P,f(WaterDrop(k).Source,WaterDrop(k).FK(u))/SumF];
U=[U,WaterDrop(k).FK(u)];
end
% Selon la probabilité P(source(k),u), De la roulette FV(k) Sélectionner un point au hasard dans l'ensemble
WaterDrop(k).Target=U(RoutedGame(P));
end
% Mettre à jour la vitesse des gouttelettes
WaterDrop(k).Vel=WaterDrop(k).Vel+av/(bv+cv*WaterDrop(k).SW(WaterDrop(k).Source,WaterDrop(k).Target)^2);
Maxtime=max(epsilon,WaterDrop(k).Vel);
% Calculer le temps nécessaire pour que les gouttelettes atteignent le point cible
TT(WaterDrop(k).Source,WaterDrop(k).Target)=D(WaterDrop(k).Source,WaterDrop(k).Target)/Maxtime;
% Calculer l'accroissement du sol transporté par les gouttelettes d'eau
DeltaSW(WaterDrop(k).Source,WaterDrop(k).Target)= as/(bs+cs*TT(WaterDrop(k).Source,WaterDrop(k).Target)^2);
% Mettre à jour la quantité de sol sur le chemin que les gouttelettes traversent
WaterDrop(k).SW(WaterDrop(k).Source,WaterDrop(k).Target)=...
(1-alpha)*WaterDrop(k).SW(WaterDrop(k).Source,WaterDrop(k).Target)-alpha*DeltaSW(WaterDrop(k).Source,WaterDrop(k).Target);
% Renouveler la quantité de terre transportée par les gouttelettes
WaterDrop(k).Soil=WaterDrop(k).Soil+DeltaSW(WaterDrop(k).Source,WaterDrop(k).Target);
% % Facteur de pénalité
% cfe1=0;
% cfe2=0;
% cfe0=0;
WaterDrop(k).S(WaterDrop(k).Target)=WaterDrop(k).S(WaterDrop(k).Source)+ss(WaterDrop(k).Source)+T(WaterDrop(k).Source,WaterDrop(k).Target);
% if WaterDrop(k).S(WaterDrop(k).Target)<E
% cfe1=EP(WaterDrop(k).Target)*abs(E-WaterDrop(k).Target);
% elseif WaterDrop(k).S(WaterDrop(k).Target)>L
% cfe1=LP(WaterDrop(k).Target)*abs(WaterDrop(k).Target-L);
% end
% cfe0=cfe1+cfe2;
if WaterDrop(k).Target==1
% （Rassemblement ordonné） Retour du véhicule de distribution au centre de distribution
WaterDrop(k).VisitNode=[WaterDrop(k).VisitNode,WaterDrop(k).Target];
% Augmentation du coût fixe de l'utilisation d'un véhicule dans la fonction objective （ Indique une boucle qui forme un véhicule ）
WaterDrop(k).ZZ=WaterDrop(k).ZZ+h*D(WaterDrop(k).Source,WaterDrop(k).Target)+D(WaterDrop(k).Source,WaterDrop(k).Target)*(R0+R1*(q(WaterDrop(k).Target)))+R;
WaterDrop(k).Source=WaterDrop(k).Target;
WaterDrop(k).Q(WaterDrop(k).Source)=0;
WaterDrop(k).S(WaterDrop(k).Source)=0;
else
% Accès àTarget(k)Point, Augmentation du nombre de points visités 1- Oui.
WaterDrop(k).VisitNode=[WaterDrop(k).VisitNode,WaterDrop(k).Target];
% Réduction d'un point de non - accès
WaterDrop(k).UnvisitNode(WaterDrop(k).UnvisitNode==WaterDrop(k).Target)=[]; %%%%%？？？？？
% La fonction objective n'augmente que le coût de fonctionnement du véhicule
WaterDrop(k).ZZ=WaterDrop(k).ZZ+h*D(WaterDrop(k).Source,WaterDrop(k).Target)+D(WaterDrop(k).Source,WaterDrop(k).Target)*(R0+R1*(q(WaterDrop(k).Target)));
% Calcul par k Charge du véhicule au point suivant
WaterDrop(k).Q(WaterDrop(k).Target)=WaterDrop(k).Q(WaterDrop(k).Source)+q(WaterDrop(k).Target);
% Calcul par k Quand le véhicule arrive au point suivant
WaterDrop(k).S(WaterDrop(k).Target)=WaterDrop(k).S(WaterDrop(k).Source)+...
ss(WaterDrop(k).Source)+T(WaterDrop(k).Source,WaterDrop(k).Target);
% Prenez le point suivant comme point de départ
WaterDrop(k).Source=WaterDrop(k).Target;
end
%Sik Le chemin des gouttelettes est terminé （ C'est - à - dire que tous les points de demande sont servis ）
end
if WaterDrop(k).ZZ<Z
Z=WaterDrop(k).ZZ;
Route=WaterDrop(k).VisitNode;
Soil=WaterDrop(k).Soil;
end
end
% Chemin correspondant à la solution optimale de ce cycle Route Mise à jour de la quantité de sol sur
M=size(Route,2);
% La quantité de sol sur les bords au - delà du chemin est constante .
for i=1:M-1
W(Route(i),Route(i+1))=(1+beta)* W(Route(i),Route(i+1))-beta*Soil/(M-1);
end
% Mettre à jour les Itérations
t=t+1;
end
TotalZ
TotalRoute
LL
DrawPath(TotalRoute,position);
% Diagramme itératif
figure
x=1:1:Iter;
y=bestZ;
plot(x,y,'r--')
%plot(x,y);
hold on
%plot(x,y1,'r--');% Adaptation moyenne
title(' Processus d'optimisation ')
xlabel('Nombre d'Itérations')
ylabel(' Adaptation optimale ')
axis([0,Iter,500,1200])


1.
2.
3.
4.
5.
6.
7.
8.
9.
10.
11.
12.
13.
14.
15.
16.
17.
18.
19.
20.
21.
22.
23.
24.
25.
26.
27.
28.
29.
30.
31.
32.
33.
34.
35.
36.
37.
38.
39.
40.
41.
42.
43.
44.
45.
46.
47.
48.
49.
50.
51.
52.
53.
54.
55.
56.
57.
58.
59.
60.
61.
62.
63.
64.
65.
66.
67.
68.
69.
70.
71.
72.
73.
74.
75.
76.
77.
78.
79.
80.
81.
82.
83.
84.
85.
86.
87.
88.
89.
90.
91.
92.
93.
94.
95.
96.
97.
98.
99.
100.
101.
102.
103.
104.
105.
106.
107.
108.
109.
110.
111.
112.
113.
114.
115.
116.
117.
118.
119.
120.
121.
122.
123.
124.
125.
126.
127.
128.
129.
130.
131.
132.
133.
134.
135.
136.
137.
138.
139.
140.
141.
142.
143.
144.
145.
146.
147.
148.
149.
150.
151.
152.
153.
154.
155.
156.
157.
158.
159.
160.
161.
162.
163.
164.
165.
166.
167.
168.
169.
170.
171.
172.
173.
174.
175.
176.
177.
178.
179.
180.
181.
182.
183.
184.
185.
186.
187.
188.
189.
190.
191.
192.
193.
194.
195.
196.
197.
198.
199.
200.
201.
202.
203.
204.
205.
206.
207.
208.
209.
210.
211.
212.
213.
214.
215.
216.
217.
218.
219.
220.
221.
222.
223.
224.
225.
226.
227.
228.
229.
230.
231.
232.
233.
234.
235.
236.
237.
238.
239.
240.
241.
242.
243.
244.
245.
246.
247.
248.
249.
250.
251.
252.
253.
254.
255.
256.
257.
258.
259.
260.
261.
262.
263.
264.
265.
266.
267.
268.