Solving Vehicle Routing Problem Using Artificial Fish Swarm Algorithm in MATLAB

  1. Algorithm Fundamentals and VRP Problem Modeling

1.1 VRP Problem Definition

The Vehicle Routing Problem (VRP) can be described as follows: m vehicles depart from a distribution center to serve n customer locations, where each customer has a specific demand quantity and each vehicle has a capacity constraint. The objective is to find the route assignment with minimum total distance.

1.2 Artificial Fish Swarm Algorithm (AFSA) Enhancement

In this implementation, the position of each artificial fish represents a VRP solution. The algorithm optimizes solution quality through four fish behaviors: prey, swarm, follow, and random movement.

  1. Complete MATLAB Implementation

2.1 Main Program Structure

%% AFSA-based VRP Solver
clear; clc; close all;

%% 1. Problem Configuration
config.custCount = 31;
config.vehCapacity = 100;
config.maxVehicles = 5;
config.demandBounds = [5, 25];

%% 2. Algorithm Parameters
algo.fishCount = 50;
algo.maxGenerations = 200;
algo perceptionRange = 0.3;
algo.moveStep = 0.2;
algo.congestionFactor = 0.6;
algo.maxAttempts = 5;
algo.randProb = 0.2;

%% 3. Initialize Problem Data
rng(2024);
[locationData, distMatrix] = generateVRPData(config);

%% 4. Initialize Fish Population
school = initFishSchool(algo, config, locationData);

%% 5. Main Loop - AFSA Optimization
bestFish = school(1);
convergenceCurve = zeros(algo.maxGenerations, 1);

for gen = 1:algo.maxGenerations
    fprintf('Generation %d/%d: ', gen, algo.maxGenerations);
    
    for i = 1:algo.fishCount
        randVal = rand();
        
        if randVal < 0.35
            school(i) = executeFollow(school, i, algo, distMatrix, locationData, config);
        elseif randVal < 0.7
            school(i) = executeSwarm(school, i, algo, distMatrix, locationData, config);
        else
            school(i) = executePrey(school(i), algo, distMatrix, locationData, config);
        end
        
        if school(i).travelDist < bestFish.travelDist
            bestFish = school(i);
        end
    end
    
    convergenceCurve(gen) = bestFish.travelDist;
    
    % Adaptive parameter adjustment
    algo.moveStep = algo.moveStep * 0.995;
    algo.perceptionRange = algo.perceptionRange * 0.997;
    
    fprintf('Best Distance: %.2f, Vehicles Used: %d\n', ...
        bestFish.travelDist, numel(bestFish.routes));
    
    if mod(gen, 20) == 0
        visualizeRoutes(bestFish, locationData, sprintf('Generation %d Best', gen));
    end
end

%% 6. Results Display
displayResults(bestFish, convergenceCurve, locationData);
showSolutionInfo(bestFish, locationData);

2.2 Data Initialization Function

function [locationData, distMatrix] = generateVRPData(cfg)
    n = cfg.custCount;
    
    % Sample city coordinates (31 locations)
    cityCoords = [
        1304, 2312; 3639, 1315; 4177, 2244; 3712, 1399; 3488, 1535;
        3326, 1556; 3238, 1229; 4196, 1004; 4312, 790; 4386, 570;
        3007, 1970; 2562, 1756; 2788, 1491; 2381, 1676; 1332, 695;
        3715, 1678; 3918, 2179; 4061, 2370; 3780, 2212; 3676, 2578;
        4029, 2838; 4263, 2931; 3429, 1908; 3507, 2367; 3394, 2643;
        3439, 3201; 2935, 3240; 3140, 3550; 2545, 2357; 2778, 2826;
        2370, 2975
    ];
    
    if n > size(cityCoords, 1)
        cityCoords = repmat(cityCoords, ceil(n/size(cityCoords,1)), 1);
    end
    locationData.coords = cityCoords(1:n, :);
    
    % Generate random demands
    locationData.demands = zeros(n, 1);
    locationData.demands(2:end) = randi(cfg.demandBounds, n-1, 1);
    
    % Compute distance matrix (Euclidean)
    distMatrix = zeros(n, n);
    for i = 1:n
        for j = 1:n
            if i ~= j
                dx = locationData.coords(i, 1) - locationData.coords(j, 1);
                dy = locationData.coords(i, 2) - locationData.coords(j, 2);
                distMatrix(i, j) = sqrt(dx^2 + dy^2);
            end
        end
    end
    
    distMatrix = distMatrix / 100;
end

2.3 Fish Population Initialization

function school = initFishSchool(algo, cfg, locData)
    n = cfg.custCount;
    
    school(algo.fishCount) = struct();
    
    for f = 1:algo.fishCount
        perm = randperm(n-1) + 1;
        
        fish = interpretChromosome(perm, locData, cfg, algo);
        fish.travelDist = computeRouteLength(fish.routes, locData.coords);
        fish.chromosome = perm;
        fish.routes = fish.routes;
        fish.vehCount = numel(fish.routes);
        
        school(f) = fish;
    end
end

function fish = interpretChromosome(perm, locData, cfg, algo)
    routes = {};
    currentPath = [1];
    currentLoad = 0;
    
    for i = 1:numel(perm)
        custIdx = perm(i);
        demand = locData.demands(custIdx);
        
        if currentLoad + demand <= cfg.vehCapacity
            currentPath = [currentPath, custIdx];
            currentLoad = currentLoad + demand;
        else
            currentPath = [currentPath, 1];
            routes{end+1} = currentPath;
            
            currentPath = [1, custIdx];
            currentLoad = demand;
        end
    end
    
    if numel(currentPath) > 1
        currentPath = [currentPath, 1];
        routes{end+1} = currentPath;
    end
    
    fish.routes = routes;
end

2.4 Fish Behavior Implementation

function newFish = executePrey(fish, algo, distMatrix, locData, cfg)
    bestCandidate = fish;
    foundBetter = false;
    
    for attempt = 1:algo.maxAttempts
        neighbor = createNeighbor(fish.chromosome, algo, distMatrix);
        
        newFish = interpretChromosome(neighbor, locData, cfg, algo);
        newFish.travelDist = computeRouteLength(newFish.routes, locData.coords);
        newFish.chromosome = neighbor;
        newFish.vehCount = numel(newFish.routes);
        
        if newFish.travelDist < bestCandidate.travelDist
            bestCandidate = newFish;
            foundBetter = true;
        end
    end
    
    if ~foundBetter
        newFish = executeRandom(fish, algo, distMatrix, locData, cfg);
    else
        newFish = bestCandidate;
    end
end

function newFish = executeSwarm(school, idx, algo, distMatrix, locData, cfg)
    n = algo.fishCount;
    currentFish = school(idx);
    
    neighborIndices = [];
    for i = 1:n
        if i ~= idx
            diff = solutionDifference(currentFish.chromosome, school(i).chromosome);
            if diff < algo.perceptionRange * numel(currentFish.chromosome)
                neighborIndices = [neighborIndices, i];
            end
        end
    end
    
    if ~isempty(neighborIndices)
        centerChrom = computeCentroid(school(neighborIndices));
        
        neighborCount = numel(neighborIndices);
        if (neighborCount / n) < algo.congestionFactor
            newChrom = interpolatePosition(currentFish.chromosome, centerChrom, algo.moveStep);
            newFish = interpretChromosome(newChrom, locData, cfg, algo);
            newFish.travelDist = computeRouteLength(newFish.routes, locData.coords);
            newFish.chromosome = newChrom;
            newFish.vehCount = numel(newFish.routes);
        else
            newFish = executePrey(currentFish, algo, distMatrix, locData, cfg);
        end
    else
        newFish = executePrey(currentFish, algo, distMatrix, locData, cfg);
    end
end

function newFish = executeFollow(school, idx, algo, distMatrix, locData, cfg)
    n = algo.fishCount;
    currentFish = school(idx);
    
    bestNeighbor = currentFish;
    bestDist = currentFish.travelDist;
    
    for i = 1:n
        if i ~= idx
            diff = solutionDifference(currentFish.chromosome, school(i).chromosome);
            if diff < algo.perceptionRange * numel(currentFish.chromosome)
                if school(i).travelDist < bestDist
                    bestNeighbor = school(i);
                    bestDist = school(i).travelDist;
                end
            end
        end
    end
    
    if bestNeighbor.travelDist < currentFish.travelDist
        nearbyCount = 0;
        for i = 1:n
            if i ~= idx
                diff = solutionDifference(bestNeighbor.chromosome, school(i).chromosome);
                if diff < algo.perceptionRange * numel(bestNeighbor.chromosome)
                    nearbyCount = nearbyCount + 1;
                end
            end
        end
        
        if (nearbyCount / n) < algo.congestionFactor
            newChrom = interpolatePosition(currentFish.chromosome, bestNeighbor.chromosome, algo.moveStep);
            newFish = interpretChromosome(newChrom, locData, cfg, algo);
            newFish.travelDist = computeRouteLength(newFish.routes, locData.coords);
            newFish.chromosome = newChrom;
            newFish.vehCount = numel(newFish.routes);
        else
            newFish = executePrey(currentFish, algo, distMatrix, locData, cfg);
        end
    else
        newFish = executePrey(currentFish, algo, distMatrix, locData, cfg);
    end
end

2.5 Helper Functions

function dist = computeRouteLength(routes, coords)
    dist = 0;
    for r = 1:numel(routes)
        path = routes{r};
        for i = 1:numel(path)-1
            fromNode = path(i);
            toNode = path(i+1);
            dx = coords(fromNode, 1) - coords(toNode, 1);
            dy = coords(fromNode, 2) - coords(toNode, 2);
            dist = dist + sqrt(dx^2 + dy^2);
        end
    end
end

function newChrom = createNeighbor(chrom, algo, distMatrix)
    n = numel(chrom);
    newChrom = chrom;
    
    opType = randi(4);
    
    switch opType
        case 1
            i = randi(n);
            j = randi(n);
            while i == j
                j = randi(n);
            end
            temp = newChrom(i);
            newChrom(i) = newChrom(j);
            newChrom(j) = temp;
            
        case 2
            i = randi(n-1);
            j = randi(n);
            while i >= j
                i = randi(n-1);
                j = randi(n);
            end
            newChrom(i:j) = fliplr(newChrom(i:j));
            
        case 3
            i = randi(n);
            j = randi(n);
            while i == j
                j = randi(n);
            end
            if i < j
                newChrom = [newChrom(1:i-1), newChrom(i+1:j), newChrom(i), newChrom(j+1:end)];
            else
                newChrom = [newChrom(1:j), newChrom(i), newChrom(j+1:i-1), newChrom(i+1:end)];
            end
            
        case 4
            if n >= 4
                idx1 = randi(n);
                candidates = 1:n;
                candidates(idx1) = [];
                probs = ones(1, n-1) / (n-1);
                idx2 = randsample(candidates, 1, true, probs);
                
                temp = newChrom(idx1);
                newChrom(idx1) = newChrom(idx2);
                newChrom(idx2) = temp;
            end
    end
end

function diff = solutionDifference(chrom1, chrom2)
    n = numel(chrom1);
    diff = 0;
    for i = 1:n
        if chrom1(i) ~= chrom2(i)
            diff = diff + 1;
        end
    end
    diff = diff / n;
end

function centerChrom = computeCentroid(popMembers)
    nFishes = numel(popMembers);
    n = numel(popMembers(1).chromosome);
    
    centerChrom = zeros(1, n);
    positions = zeros(nFishes, n);
    
    for f = 1:nFishes
        positions(f, :) = popMembers(f).chromosome;
    end
    
    for pos = 1:n
        clients = positions(:, pos);
        uniqueClients = unique(clients);
        counts = histcounts(clients, [uniqueClients; max(uniqueClients)+1]);
        
        [~, maxIdx] = max(counts);
        centerChrom(pos) = uniqueClients(maxIdx);
    end
    
    missing = setdiff(2:max(centerChrom), centerChrom);
    if ~isempty(missing)
        duplicates = find(diff(sort(centerChrom)) == 0);
        if ~isempty(duplicates)
            centerChrom(duplicates(1)) = missing(1);
        end
    end
end

2.6 Visualization Functions

function visualizeRoutes(fish, locData, titleStr)
    figure(1);
    clf;
    hold on;
    
    routes = fish.routes;
    colorMap = lines(numel(routes));
    
    scatter(locData.coords(1,1), locData.coords(1,2), 200, 'r', 'p', 'LineWidth', 2);
    scatter(locData.coords(2:end,1), locData.coords(2:end,2), 100, 'b', 'o', 'filled');
    
    for r = 1:numel(routes)
        path = routes{r};
        x = locData.coords(path, 1);
        y = locData.coords(path, 2);
        
        plot(x, y, '-', 'Color', colorMap(r,:), 'LineWidth', 2);
        plot(x, y, 'o', 'Color', colorMap(r,:), 'MarkerFaceColor', colorMap(r,:));
        
        for i = 1:numel(path)-1
            midX = (x(i) + x(i+1)) / 2;
            midY = (y(i) + y(i+1)) / 2;
            text(midX, midY, sprintf('%d', r), 'FontSize', 8);
        end
    end
    
    text(locData.coords(1,1), locData.coords(1,2), ' Depot', ...
        'FontSize', 10, 'FontWeight', 'bold');
    
    xlabel('X Coordinate');
    ylabel('Y Coordinate');
    title(sprintf('%s\nTotal Distance: %.2f, Vehicles: %d', titleStr, fish.travelDist, numel(routes)));
    grid on;
    axis equal;
    hold off;
    
    drawnow;
end

function displayResults(bestFish, convergenceHistory, locData)
    figure('Position', [100, 100, 1200, 400]);
    
    subplot(1, 3, 1);
    plot(convergenceHistory, 'b-', 'LineWidth', 2);
    xlabel('Generation');
    ylabel('Total Distance');
    title('AFSA Convergence');
    grid on;
    
    avgCost = mean(convergenceHistory);
    stdCost = std(convergenceHistory);
    text(0.05, 0.95, sprintf('Best: %.2f\nMean: %.2f\nStd: %.2f', ...
        convergenceHistory(end), avgCost, stdCost), ...
        'Units', 'normalized', 'FontSize', 9, ...
        'VerticalAlignment', 'top', 'BackgroundColor', 'w');
    
    subplot(1, 3, [2, 3]);
    visualizeRoutes(bestFish, locData, 'Final Solution');
end

function showSolutionInfo(fish, locData)
    fprintf('\n========== Solution Details ==========\n');
    fprintf('Total Distance: %.2f\n', fish.travelDist);
    fprintf('Vehicles Used: %d\n', numel(fish.routes));
    fprintf('\nRoute Breakdown:\n');
    
    for v = 1:numel(fish.routes)
        path = fish.routes{v};
        pathStr = sprintf('%d', path(1));
        totalDemand = 0;
        
        for i = 2:numel(path)
            pathStr = sprintf('%s -> %d', pathStr, path(i));
            if path(i) ~= 1
                totalDemand = totalDemand + locData.demands(path(i));
            end
        end
        
        fprintf('Vehicle %d: %s\n', v, pathStr);
        fprintf('  Customers: %d, Demand: %d\n', ...
            numel(path)-2, totalDemand);
        
        pathDist = 0;
        for i = 1:numel(path)-1
            fromNode = path(i);
            toNode = path(i+1);
            dx = locData.coords(fromNode, 1) - locData.coords(toNode, 1);
            dy = locData.coords(fromNode, 2) - locData.coords(toNode, 2);
            pathDist = pathDist + sqrt(dx^2 + dy^2);
        end
        fprintf('  Route Distance: %.2f\n\n', pathDist);
    end
end
  1. Parameter Tuning Guidelines

3.1 Key Parameter Configuration

%经验参数配置（针对不同规模问题）
if cfg.custCount <= 20
    algo.fishCount = 30;
    algo.maxGenerations = 100;
    algo.perceptionRange = 0.4;
    algo.moveStep = 0.3;
elseif cfg.custCount <= 50
    algo.fishCount = 50;
    algo.maxGenerations = 200;
    algo.perceptionRange = 0.3;
    algo.moveStep = 0.2;
else
    algo.fishCount = 80;
    algo.maxGenerations = 300;
    algo.perceptionRange = 0.2;
    algo.moveStep = 0.15;
end

3.2 Adaptive Parameter Adjustment

if mod(gen, 20) == 0
    improvementRatio = (convergenceHistory(max(1, gen-10)) - convergenceHistory(gen)) / convergenceHistory(max(1, gen-10));
    
    if improvementRatio < 0.01
        algo.moveStep = algo.moveStep * 1.2;
        algo.perceptionRange = min(0.5, algo.perceptionRange * 1.1);
    else
        algo.moveStep = max(0.05, algo.moveStep * 0.95);
        algo.perceptionRange = max(0.1, algo.perceptionRange * 0.98);
    end
end
  1. Algorithm Extensions and Variants

4.1 VRP with Time Windows (VRPTW)

function fish = decodeWithTimeWindows(perm, locData, cfg, algo)
    routes = {};
    currentPath = [1];
    currentLoad = 0;
    currentTime = 0;
    
    for i = 1:numel(perm)
        custIdx = perm(i);
        demand = locData.demands(custIdx);
        svcTime = locData.serviceTime(custIdx);
        timeWin = locData.timeWindows(custIdx, :);
        
        travelTime = computeTravelTime(currentPath(end), custIdx, locData);
        arrivalTime = currentTime + travelTime;
        
        if currentLoad + demand <= cfg.vehCapacity && ...
           arrivalTime <= timeWin(2)
            
            if arrivalTime < timeWin(1)
                waitTime = timeWin(1) - arrivalTime;
                arrivalTime = timeWin(1);
            else
                waitTime = 0;
            end
            
            currentPath = [currentPath, custIdx];
            currentLoad = currentLoad + demand;
            currentTime = arrivalTime + svcTime + waitTime;
        else
            if numel(currentPath) > 1
                returnTime = computeTravelTime(currentPath(end), 1, locData);
                currentPath = [currentPath, 1];
                routes{end+1} = currentPath;
            end
            
            currentPath = [1, custIdx];
            currentLoad = demand;
            travelTime = computeTravelTime(1, custIdx, locData);
            arrivalTime = travelTime;
            
            if arrivalTime < timeWin(1)
                waitTime = timeWin(1) - arrivalTime;
                currentTime = timeWin(1) + svcTime;
            else
                currentTime = arrivalTime + svcTime;
            end
        end
    end
    
    if numel(currentPath) > 1
        currentPath = [currentPath, 1];
        routes{end+1} = currentPath;
    end
    
    fish.routes = routes;
end

4.2 Hybrid AFSA Approach

function newFish = hybridBehavior(fish, algo, distMatrix, locData, cfg)
    if rand() < 0.5
        newFish = swarmWrapper(fish, algo, distMatrix, locData, cfg);
    else
        newFish = followWrapper(fish, algo, distMatrix, locData, cfg);
    end
    
    if rand() < 0.3
        newFish = localRefinement(newFish, algo, locData, cfg);
    end
end

function improved = localRefinement(fish, algo, locData, cfg)
    improved = fish;
    
    for iter = 1:10
        opType = randi(3);
        newChrom = improved.chromosome;
        
        switch opType
            case 1
                i = randi(numel(newChrom)-1);
                j = randi(numel(newChrom));
                while i >= j
                    i = randi(numel(newChrom)-1);
                    j = randi(numel(newChrom));
                end
                newChrom(i:j) = fliplr(newChrom(i:j));
                
            case 2
            case 3
        end
        
        tempFish = interpretChromosome(newChrom, locData, cfg, algo);
        tempFish.travelDist = computeRouteLength(tempFish.routes, locData.coords);
        
        if tempFish.travelDist < improved.travelDist
            improved = tempFish;
            improved.chromosome = newChrom;
        end
    end
end
  1. Usage Instructions and Recommendations

5.1 Running the Program

  1. Save all functions as separate .m files or combine in one main file
  2. Adjust problem parameetrs as needed (customer count, vehicle capacity, etc.)
  3. Execute the main script
  4. Review visualization outputs and detailed results

5.2 Performance Optimization Tips

  1. Parallel Computing: Use parfor for parallel fish behavior evaluation
  2. Memory Archive: Maintain an external archive for historical best solutions
  3. Multi-Population Strategy: Employ multiple subpopulations for parallel search with preiodic information exchange
  4. Heuristic Initialization: Use nearest neighbor or savings algorithm for initial solution generation

5.3 Real-World Application

function locData = importRealData(filename)
    data = readtable(filename);
    locData.coords = [data.X, data.Y];
    locData.demands = data.Demand;
    locData.serviceTime = data.ServiceTime;
    
    if ismember('TimeWindow', data.Properties.VariableNames)
        locData.timeWindows = [data.EarliestTime, data.LatestTime];
    end
end

function results = benchmarkAFSA(problemSizes, numRuns)
    results = struct();
    
    for p = 1:numel(problemSizes)
        cfg.custCount = problemSizes(p);
        cfg.vehCapacity = 100;
        
        execTimes = zeros(numRuns, 1);
        bestCosts = zeros(numRuns, 1);
        
        for r = 1:numRuns
            tic;
            [bestFish, ~] = runAFS(cfg);
            execTimes(r) = toc;
            bestCosts(r) = bestFish.travelDist;
        end
        
        results(p).size = problemSizes(p);
        results(p).avgTime = mean(execTimes);
        results(p).stdTime = std(execTimes);
        results(p).avgCost = mean(bestCosts);
        results(p).bestCost = min(bestCosts);
        results(p).worstCost = max(bestCosts);
    end
end

This implementation provides a complete framework for solving VRP using Artificial Fish Swarm Algorithm, including core algorithm logic, visualization tools, parameter tuning, and extension capabilities.

Tags: MATLAB vehicle-routing-problem artificial-fish-swarm-algorithm Metaheuristics Optimization

Posted on Sun, 25 Oct 2026 16:29:38 +0000 by FSGr33n

Hot Tags

©2026 Freaks City