MATLAB的遗传算法解决车辆路径规划问题

%% 1. 数据准备（以20客户单仓库为例）
clear; clc;
num_customers = 20; % 客户数量
num_vehicles = 4;   % 车辆数量
vehicle_capacity = 100; % 车辆容量
[coordinates, demands] = generate_vrp_data(num_customers); % 生成测试数据
​
%% 2. 遗传算法参数设置
options = optimoptions('ga',...
    'PopulationSize', 100,...
    'MaxGenerations', 500,...
    'CrossoverFcn', @cxOrdered,...
    'MutationFcn', @mutSwap,...
    'SelectionFcn', @selectiontournament,...
    'PlotFcn', {@gaplotbestf,@gaplotdistance});
​
%% 3. 染色体编码（基于仓库分隔符）
chromosome = ; % 示例染色体
​
%% 4. 适应度函数（含容量约束）
function fitness = calculateFitness(chromosome, coordinates, demands, vehicle_capacity)
    num_vehicles = sum(chromosome == 0) - 1; % 车辆数=分隔符数-1
    total_distance = 0;
    valid = true;
    
    % 路径分割
    routes = strsplit(num2str(chromosome), '0');
    for i = 1:numel(routes)
        route = routes{i};
        if isempty(route)
            continue;
        end
        route = route(route ~= ' '); % 去除空格
        route = str2num(route); 
        
        % 容量检查
        load_sum = sum(demands(route));
        if load_sum > vehicle_capacity
            valid = false;
            break;
        end
        
        % 距离计算
        current_pos = ; % 仓库坐标
        for j = 1:numel(route)
            next_pos = coordinates(route(j),:);
            total_distance = total_distance + norm(current_pos - next_pos);
            current_pos = next_pos;
        end
        total_distance = total_distance + norm(current_pos, ); % 返回仓库
    end
    
    if valid
        fitness = 1 / total_distance; % 适应度与距离成反比
    else
        fitness = 0; % 非法解适应度为0
    end
end
​
%% 5. 遗传算法主循环
nvars = num_customers + 1; % 染色体长度（客户数+仓库分隔符）
lb = ; % 下限
ub = [num_customers,num_customers,num_customers](@ref); % 上限
​
% 运行遗传算法
[best_solution, best_fitness] = ga(@(x) calculateFitness(x, coordinates, demands, vehicle_capacity),...
    nvars, [], [], [], [], lb, ub, [], options);
​
%% 6. 结果解析与可视化
best_chromosome = best_solution;
disp('最优路径:');
disp(best_chromosome);
​
% 绘制路径
figure;
hold on;
plot(coordinates(:,2), coordinates(:,3),'ko'); % 客户点
plot(coordinates(1,2), coordinates(1,3),'ro'); % 仓库
routes = strsplit(num2str(best_chromosome), '0');
for i = 1:numel(routes)
    route = routes{i};
    route = route(route ~= ' ');
    route = str2num(route);
    if isempty(route)
        continue;
    end
    route_coords = coordinates(route,:);
    plot(route_coords(:,2), route_coords(:,3),'b-');
end
hold off;
title(sprintf('最优路径（总距离=%.2f）', 1/best_fitness));

% 容量约束检查函数
function valid = check_capacity(chromosome, demands, capacity)
    index = find(chromosome == 0);
    index = [0,index,length(chromosome)+1](@ref);
    for i = 2:numel(index)
        route = chromosome(index(i-1)+1:index(i)-1);
        if sum(demands(route)) > capacity
            valid = false;
            return;
        end
    end
    valid = true;
end

% 启用并行计算池
parpool('local');
​
% 并行适应度评估
fitness = gcp().JobManager.evaluate(@(x) calculateFitness(x, coordinates, demands, vehicle_capacity), population);

% 2-opt局部优化
function improved = two_opt(route)
    improved = route;
    n = length(route);
    for i = 1:n-2
        for j = i+2:n
            new_route = route;
            new_route(i:j) = route(j:-1:i);
            if calculateRouteDistance(new_route) < calculateRouteDistance(improved)
                improved = new_route;
            end
        end
    end
end

% 保留前5%最优个体
elite_ratio = 0.05;
num_elites = round(elite_ratio * pop_size);
[~, idx] = mink(fitness, num_elites);
elite_population = population(idx,:);

% 多仓库路径规划流程
function multi_depot_vrp()
    % 仓库分配阶段
    depot_assign = kmeans(customer_coords, num_depots);
    
    % 车辆路径规划
    for i = 1:num_depots
        sub_customers = customers(depot_assign == i);
        sub_route = ga_solve(sub_customers);
        global_route = [global_route, sub_route, '0'];
    end
end

% 时间窗约束处理
function feasible = check_time_window(route, time_windows)
    current_time = 0;
    for i = 1:length(route)
        customer = route(i);
        ready_time = time_windows(customer,1);
        due_time = time_windows(customer,2);
        
        if current_time < ready_time
            current_time = ready_time;
        end
        service_time = 10; % 假设服务时间10分钟
        current_time = current_time + service_time;
        
        if current_time > due_time
            feasible = false;
            return;
        end
    end
    feasible = true;
end