% 初始化参数
Np = 100; % 种群大小
Nc = 5; % 染色体数量
Nv = 2; % 变量数量
MaxIt = 300; % 最大迭代次数
INF = 1e6; % 适应度初始值
% 初始化变量范围
x_min = [0 0 0 0 0];
x_max = [6 6 6 6 6];
% 初始化需求点和仓库信息
demandF = 4;
demandG = 4;
warehouses = ["A","B","C","D","E"]; % 仓库名称
stock = [6, 6, 0, 4, 6]; % 仓库库存
distanceF = [1360, 1181, 2388, 1368, 1840]; % 需求点F到仓库距离
speedF = [98*0.8, 101, 87*0.8, 110*0.6, 93*0.8]; % 需求点F到仓库运输速度
costF = [15, 9, 19, 11, 10]; % 需求点F到仓库单位运输成本
distanceG = [2708, 2526, 3509, 548, 3001]; % 需求点G到仓库距离
speedG = [105, 80, 85*0.6, 112, 54]; % 需求点G到仓库运输速度
costG = [17, 9, 16, 6, 12]; % 需求点G到仓库单位运输成本
% 初始化种群
P = zeros(Np, Nc, Nv);
for i = 1:Np
for j = 1:Nc
P(i, j, :) = unifrnd(x_min(j), x_max(j), [1,Nv]);
end
end
% 差分进化算法主程序
F = zeros(1, Np); % 初始化适应度
bestF = NaN(MaxIt, 1); % 记录每次迭代的最优适应度
bestX = NaN(MaxIt, Nv); % 记录每次迭代的最优解
for t = 1:MaxIt
for i = 1:Np
% 随机选择3个个体
idx = randsample(1:Np, 3);
a = P(idx(1), :, :);
b = P(idx(2), :, :);
c = P(idx(3), :, :);
% 随机选择一个变量维度
j = randi(Nc);
% 差分变异操作
v = a(j,:) + rand*(b(j,:) - c(j,:));
v = min(v, x_max);
v = max(v, x_min);
% 交叉操作
if rand < 0.9
u = zeros(1, Nc, Nv);
k = randi(Nc, [1,Nc]);
for jj = 1:Nc
if jj == j || rand < 0.1
u(1, k(jj), :) = v;
else
u(1, k(jj), :) = P(i, k(jj), :);
end
end
else
u = P(i,:,:);
end
% 适应度评价
v_idx = randi(Nv); % 随机选择一个变量维度
cost = 0; % 总成本
timeF = 0; % 运输需求点F货物总时间
timeG = 0; % 运输需求点G货物总时间
for k = 1:Nc
warehouse = warehouses(k);
demandF_k = u(1, k, 1);
demandG_k = u(1, k, 2);
% 运输需求点F货物
if demandF_k > 0
time = distanceF(k)/speedF(k);
cost = cost + demandF_k*costF(k);
timeF = timeF + time;
% 库存不足，惩罚适应度
if demandF_k > stock(k)
cost = cost + INF;
end
end
% 运输需求点G货物
if demandG_k > 0
time = distanceG(k)/speedG(k);
cost = cost + demandG_k*costG(k);
timeG = timeG + time;
% 库存不足，惩罚适应度
if demandG_k > stock(k)
cost = cost + INF;
end
end
% 变异
if v_idx == 1
u(1, k, 1) = v; % 修改第一个变量
else
u(1, k, 2) = v; % 修改第二个变量
end
end
F(i) = cost;
P(i, :, :) = u;
end
% 记录每次迭代的最优解和适应度
[bestF(t), idx] = min(F);
bestX(t, :) = P(idx, :, :);
end
% 绘制适应度变化曲线
figure;
plot(bestF);
title('Fitness Curve');
xlabel('Generation');
ylabel('Fitness');
% 输出最优解
disp("===============");
disp("Optimization Results:");
disp("===============");
disp("Minimum Cost:");
disp(bestF(end));
disp("Optimal Solution:");
disp(bestX(end, :));
% 绘制桶图
figure;
clf;
hold on;
% 需求点F
bar(1:5, squeeze(bestX(end, :, 1)), 'r');
% 需求点G
bar(1:5, squeeze(bestX(end, :, 2)), 'b');
% 仓库库存
plot(1:5, stock, 'k--');
hold off;
title('Optimal Solution');
xlabel('Warehouse');
ylabel('Amount');
legend('Demand Point F', 'Demand Point G', 'Stock');
麻烦帮我看看 v = a(j,:) + rand*(b(j,:) - c(j,:));说这句索引超出数组边界
就是数组不正确呀呀呀