Triển khai thuật toán tối ưu bầy đàn (PSO) trong MATLAB để giải quyết bài toán định tuyến phương tiện (VRP). Phương pháp này được sử dụng để giải bài toán định tuyến phương tiện có cửa sổ thời gian (VRPTW).

Triển khai mã nguồn MATLAB

1. Hàm thuật toán PSO

function [optimalPath, minimalCost] = psAlgorithm_vrp(customerCount, vehicleCount, demandList, vehicleCapacity, distanceMatrix, timeConstraints, iterationLimit, swarmSize)
    % Tham số thuật toán
    cognitiveFactor = 1.5; % Hệ số học tập cá nhân
    socialFactor = 1.5; % Hệ số học tập xã hội
    inertiaWeight = 0.8; % Trọng số quán tính
    totalNodes = customerCount + 1; % Bao gồm kho trung tâm
    dimension = totalNodes * 2; % Kích thước không gian tìm kiếm

    % Khởi tạo bầy đàn
    swarmPositions = rand(swarmSize, dimension); % Vị trí ban đầu ngẫu nhiên
    velocityMatrix = zeros(swarmSize, dimension); % Ma trận vận tốc
    personalBest = swarmPositions; % Vị trí tốt nhất cá thể
    personalBestValue = inf(swarmSize, 1); % Giá trị tốt nhất cá thể
    globalBest = swarmPositions(1, :); % Vị trí tốt nhất toàn cục
    globalBestValue = inf; % Giá trị tốt nhất toàn cục

    % Vòng lặp chính
    for epoch = 1:iterationLimit
        for idx = 1:swarmSize
            % Giải mã vị trí hạt
            solution = decodeSwarmPosition(swarmPositions(idx, :), totalNodes, vehicleCount, demandList, vehicleCapacity, timeConstraints);
            fitness = calculateSolutionFitness(solution, distanceMatrix, demandList, vehicleCapacity);

            % Cập nhật cá nhân tốt nhất
            if fitness < personalBestValue(idx)
                personalBest(idx, :) = swarmPositions(idx, :);
                personalBestValue(idx) = fitness;
            end

            % Cập nhật toàn cục tốt nhất
            if fitness < globalBestValue
                globalBest = swarmPositions(idx, :);
                globalBestValue = fitness;
            end
        end

        % Cập nhật vận tốc và vị trí
        for idx = 1:swarmSize
            r1 = rand();
            r2 = rand();
            velocityMatrix(idx, :) = inertiaWeight * velocityMatrix(idx, :) + ...
                cognitiveFactor * r1 * (personalBest(idx, :) - swarmPositions(idx, :)) + ...
                socialFactor * r2 * (globalBest - swarmPositions(idx, :));
            swarmPositions(idx, :) = swarmPositions(idx, :) + velocityMatrix(idx, :);
        end

        fprintf('Epoch %d: Global Best Value = %.2f\n', epoch, globalBestValue);
    end

    % Kết quả cuối cùng
    optimalPath = decodeSwarmPosition(globalBest, totalNodes, vehicleCount, demandList, vehicleCapacity, timeConstraints);
    minimalCost = globalBestValue;
end

2. Hàm giải mã hạt

function solutionPath = decodeSwarmPosition(positionVector, nodeCount, fleetSize, requirementList, maxCapacity, timeLimits)
    % Chuyển đổi vị trí hạt thành lộ trình phương tiện
    vehicleSchedules = cell(1, fleetSize);
    for vehicle = 1:fleetSize
        vehicleSchedules{vehicle} = [1]; % Bắt đầu từ kho
    end

    % Phân bổ nhiệm vụ cho phương tiện
    for index = 1:nodeCount-1
        vehicleAssignment = round(positionVector(index));
        taskAssignment = round(positionVector(index + nodeCount - 1));
        if taskAssignment > 1 && requirementList(taskAssignment) <= maxCapacity
            vehicleSchedules{vehicleAssignment} = [vehicleSchedules{vehicleAssignment}, taskAssignment];
        end
    end

    % Kết thúc quay về kho
    for vehicle = 1:fleetSize
        vehicleSchedules{vehicle} = [vehicleSchedules{vehicle}, 1];
    end

    solutionPath = vehicleSchedules;
end

3. Hàm đánh giá lộ trình

function fitnessScore = calculateSolutionFitness(pathSolution, distanceData, demandInfo, capacityLimit)
    % Tính toán chi phí của lộ trình
    fitnessScore = 0;
    for vehicleIndex = 1:length(pathSolution)
        currentPath = pathSolution{vehicleIndex};
        currentLoad = 0;
        for step = 1:length(currentPath)-1
            fitnessScore = fitnessScore + distanceData(currentPath(step), currentPath(step+1));
            currentLoad = currentLoad + demandInfo(currentPath(step));
            if currentLoad > capacityLimit
                fitnessScore = inf; % Phạt do quá tải
                break;
            end
        end
    end
end

4. Chương trình chính

% Chương trình chính
clc;
clear;

% Thiết lập tham số
customerNumber = 10; % Số lượng khách hàng
vehicleNumber = 3; % Số lượng xe
demandArray = [0, 1, 1, 2, 2, 1, 1, 2, 1, 1, 1]; % Yêu cầu của từng điểm, kho là 0
truckCapacity = 3; % Sức chứa mỗi xe
distanceTable = randi(100, customerNumber+1, customerNumber+1); % Ma trận khoảng cách ngẫu nhiên
timeFrame = [0, 100; 10, 20; 15, 30; 20, 40; 25, 50; 30, 60; 35, 70; 40, 80; 45, 90; 50, 100; 55, 110]; % Cửa sổ thời gian
maxIterations = 100; % Số lần lặp tối đa
particleCount = 50; % Số lượng hạt

% Giải bài toán VRP
[optimalPath, minimalCost] = psAlgorithm_vrp(customerNumber, vehicleNumber, demandArray, truckCapacity, distanceTable, timeFrame, maxIterations, particleCount);

% Hiển thị kết quả
disp('Lộ trình tối ưu:');
disp(optimalPath);
fprintf('Chi phí thấp nhất: %.2f\n', minimalCost);

Mô tả chi tiết

1. Thuật toán PSO:

• Mỗi hạt đại diện cho một giải pháp tiềm năng, vị trí và vận tốc được cập nhật dựa trên vị trí tốt nhất cá nhân và toàn cục.
• Sử dụng hàm giải mã để chuyển đổi vị trí hạt thành lộ trình phương tiện cụ thể.

2. Đánh giá lộ trình:

• Dựa vào ma trận khoảng cách và sức chứa xe để đánh giá chi phí lộ trình.
• Vi phạm tải trọng sẽ bị phạt để đảm bảo tính khả thi của nghiệm.

3. Chương trình chính:

• Thiết lập các tham số bài toán, khởi tạo bầy đàn và thực hiện quá trình tối ưu hóa lặp lại.

Lưu ý quan trọng

• Hiệu chỉnh tham số: Điều chỉnh các thông số PSO như trọng số quán tính, hệ số học tập để tối ưu hiệu suất.
• Xử lý ràng buộc: Đối với các ràng buộc như cửa sổ thời gian, có thể sử dụng hàm phạt để cải thiện đánh giá độ phù hợp.
• Nâng cấp thuật toán: Kết hợp với các kỹ thuật khác như luyện kim mô phỏng có thể cải thiện thêm hiệu suất thuật toán.