需求
假设你是小王,你的毕业论文是仿照这个写一个“面相交通走廊的氢能站设计”(选址优化),想的是选几个城市,然后在这几个城市之间选址。就是加氢站选址优化,在几个城市之间选址。最后得出来一个选址的地点就行。现在你的任务就是写matlab代码,生成matlab代码之后,将各个函数的结构用中括号表明层次关系。需要的数据另外作为一个.inp文件放在主函数同一层
各个函数结构
mainHydrogenStationLocation
├── loadInputData
├── initializePopulation
├── calculateFitness
│ ├── calculateConstructionCost
│ ├── calculateCoverageRate
│ └── calculateAvgDistance
├── updatePopulationGA
│ ├── selection
│ ├── crossover
│ └── mutate
└── visualizeResults
输入数据说明:
hydrogen_station_input.inp
文件包含了模型所需的所有输入数据:- 城市信息(名称、坐标)
- 需求点信息(坐标、需求量)
- 模型参数(加氢站数量、种群大小、迭代次数等)
模型特点:
- 专注于城市间交通走廊的氢能站选址
- 考虑建设成本、需求覆盖率和服务距离三个主要因素
- 使用遗传算法进行优化求解
- 提供可视化结果展示
- 所有数据都从外部输入文件读取,便于修改和扩展
代码
main_hydrogen_station_location.m
% 主函数:面相交通走廊的氢能站选址优化模型
function [bestSolution, bestFitness, convergenceCurve] = mainHydrogenStationLocation()
% 验证文件路径
inputFile = fullfile(pwd, 'hydrogen_station_input.inp');
if ~exist(inputFile, 'file')
error('输入文件必须命名为hydrogen_station_input.inp并放在当前目录');
end
% 读取数据
[cities, demandPoints, params] = loadInputData(inputFile);
% 参数验证
if params.nStations <= 0 || isempty(cities)
error('参数读取错误,请检查输入文件格式');
end
dim = params.nStations * 2;
% 初始化种群
population = initializePopulation(params.popSize, dim, cities);
% 优化循环
convergenceCurve = zeros(params.maxGen, 1);
gbest = population(1,:);
gbestFitness = inf;
for gen = 1:params.maxGen
fitness = calculateFitness(population, cities, demandPoints, params);
% 更新最优解
[minFitness, idx] = min(fitness);
if minFitness < gbestFitness
gbest = population(idx,:);
gbestFitness = minFitness;
end
convergenceCurve(gen) = gbestFitness;
% 种群更新
population = updatePopulationGA(population, fitness, cities, params);
% 进度显示
if mod(gen,50) == 0
fprintf('迭代 %d: 最优适应度=%.2f\n', gen, gbestFitness);
end
end
% 结果可视化
visualizeResults(gbest, cities, demandPoints, params, convergenceCurve);
% 输出结果
bestSolution = reshape(gbest, 2, params.nStations)';
bestFitness = gbestFitness;
disp('最优加氢站坐标:');
disp(bestSolution);
end
load_input_data.m
function [cities, demandPoints, params] = loadInputData(inputFile)
% 增加文件存在性检查
if ~exist(inputFile, 'file')
error('输入文件不存在: %s', inputFile);
end
fid = fopen(inputFile, 'r');
if fid == -1
error('无法打开文件: %s', inputFile);
end
% 读取城市数据
nCities = fscanf(fid, '%d', 1);
cities = struct('name', cell(1,nCities), 'x', [], 'y', []);
for i = 1:nCities
cities(i).name = fscanf(fid, '%s', 1);
cities(i).x = fscanf(fid, '%f', 1);
cities(i).y = fscanf(fid, '%f', 1);
fscanf(fid, '\n'); % 跳过换行符
end
% 读取需求点数据
nDemandPoints = fscanf(fid, '%d', 1);
demandPoints = struct('x', [], 'y', [], 'demand', []);
for i = 1:nDemandPoints
demandPoints(i).x = fscanf(fid, '%f', 1);
demandPoints(i).y = fscanf(fid, '%f', 1);
demandPoints(i).demand = fscanf(fid, '%f', 1);
fscanf(fid, '\n');
end
% 读取算法参数
params.nStations = fscanf(fid, '%d', 1);
params.popSize = fscanf(fid, '%d', 1);
params.maxGen = fscanf(fid, '%d', 1);
params.costWeight = fscanf(fid, '%f', 1);
params.coverageWeight = fscanf(fid, '%f', 1);
params.distanceWeight = fscanf(fid, '%f', 1);
params.serviceRadius = fscanf(fid, '%f', 1);
fclose(fid);
% 打印验证参数
fprintf('成功读取参数:\n');
fprintf('加氢站数量: %d\n', params.nStations);
fprintf('城市数量: %d\n', length(cities));
end
initialize_population.m
% 初始化种群
function population = initializePopulation(popSize, dim, cities)
% 在城市之间的交通走廊上生成初始解
nCities = length(cities); % 修改此处,直接获取结构体数组长度
population = zeros(popSize, dim);
for i = 1:popSize
for j = 1:(dim/2)
% 随机选择两个城市
city1 = randi(nCities);
city2 = randi(nCities);
% 在这两个城市之间的连线上随机生成一个点
r = rand;
population(i, 2*j-1) = cities(city1).x + r * (cities(city2).x - cities(city1).x); % 修改此处,通过结构体数组索引访问
population(i, 2*j) = cities(city1).y + r * (cities(city2).y - cities(city1).y); % 修改此处,通过结构体数组索引访问
end
end
end
calculate_fitness.m
% 计算适应度
function fitness = calculateFitness(population, cities, demandPoints, params)
popSize = size(population, 1);
fitness = zeros(popSize, 1);
for i = 1:popSize
locations = reshape(population(i,:), 2, params.nStations)';
constructionCost = calculateConstructionCost(locations, cities, params);
coverageRate = calculateCoverageRate(locations, demandPoints, params);
avgDistance = calculateAvgDistance(locations, demandPoints);
% 修正后的适应度公式(全部正向加权)
fitness(i) = params.costWeight * constructionCost ...
+ (1 - coverageRate) * params.coverageWeight ...
+ avgDistance * params.distanceWeight;
end
end
calculate_construction_cost.m
% 计算建设成本
function cost = calculateConstructionCost(locations, cities, params)
nStations = size(locations, 2);
baseCost = 1000;
cost = 0;
for i = 1:nStations
% 计算到所有城市的距离
distToCities = arrayfun(@(c) norm([c.x, c.y] - locations(:,i)'), cities);
minDist = min(distToCities);
cost = cost + baseCost + minDist * 10;
end
end
calculate_coverage_rate.m
% 计算需求覆盖率
function coverage = calculateCoverageRate(locations, demandPoints, params)
% 修正结构体数组访问方式
nPoints = length(demandPoints); % 正确获取需求点数量
totalDemand = sum([demandPoints.demand]); % 结构体数组字段需要转换为数值数组
coveredDemand = 0;
for i = 1:nPoints
% 获取当前需求点坐标和需求
currentX = demandPoints(i).x;
currentY = demandPoints(i).y;
currentDemand = demandPoints(i).demand;
% 检查所有加氢站
for j = 1:size(locations, 2)
stationX = locations(1,j);
stationY = locations(2,j);
dist = norm([stationX-currentX, stationY-currentY]);
if dist <= params.serviceRadius
coveredDemand = coveredDemand + currentDemand;
break; % 只要有一个覆盖就跳出循环
end
end
end
coverage = coveredDemand / totalDemand;
end
calculate_avg_distance.m
% 计算需求点到最近加氢站的平均距离
function avgDist = calculateAvgDistance(locations, demandPoints)
% 修正结构体数组访问方式
nPoints = length(demandPoints);
totalDist = 0;
for i = 1:nPoints
demandX = demandPoints(i).x;
demandY = demandPoints(i).y;
minDist = inf;
for j = 1:size(locations, 2)
stationX = locations(1,j);
stationY = locations(2,j);
dist = norm([stationX-demandX, stationY-demandY]);
if dist < minDist
minDist = dist;
end
end
totalDist = totalDist + minDist;
end
avgDist = totalDist / nPoints;
end
update_population_ga.m
% 使用遗传算法更新种群
function newPopulation = updatePopulationGA(population, fitness, cities, params)
popSize = size(population, 1);
newPopulation = zeros(popSize, size(population,2));
% 锦标赛选择
tournamentSize = 3;
for i = 1:popSize
% 选择
parent1 = selection(population, fitness, tournamentSize);
parent2 = selection(population, fitness, tournamentSize);
% 交叉 (80%概率)
if rand < 0.8
[child1, child2] = crossover(parent1, parent2);
newPopulation(i,:) = child1;
else
newPopulation(i,:) = parent1;
end
% 变异 (10%概率)
if rand < 0.1
newPopulation(i,:) = mutate(newPopulation(i,:), cities, params);
end
end
end
selection.m
% 选择操作(锦标赛选择)
function selected = selection(population, fitness, tournamentSize)
popSize = size(population, 1);
% 随机选择几个个体进行锦标赛
tournamentIdx = randperm(popSize, tournamentSize);
tournamentFitness = fitness(tournamentIdx);
% 找到最优个体
[~, winnerIdx] = min(tournamentFitness);
selected = population(tournamentIdx(winnerIdx), :);
end
crossover.m
% 交叉操作
function [child1, child2] = crossover(parent1, parent2)
dim = length(parent1);
% 单点交叉
crossoverPoint = randi([1, dim-1]);
child1 = [parent1(1:crossoverPoint), parent2(crossoverPoint+1:end)];
child2 = [parent2(1:crossoverPoint), parent1(crossoverPoint+1:end)];
end
mutate.m
% 变异操作
function mutant = mutate(individual, cities, params)
dim = length(individual);
mutant = individual;
% 随机选择变异点
mutationPoint = randi(dim);
% 确保访问城市坐标的正确方式
cityCoords = [[cities.x]; [cities.y]]';
if mod(mutationPoint,2) == 1 % X坐标
c1 = randi(length(cities));
c2 = randi(length(cities));
mutant(mutationPoint) = cityCoords(c1,1) + rand*(cityCoords(c2,1)-cityCoords(c1,1));
else % Y坐标
c1 = randi(length(cities));
c2 = randi(length(cities));
mutant(mutationPoint) = cityCoords(c1,2) + rand*(cityCoords(c2,2)-cityCoords(c1,2));
end
end
visualize_results.m
% 可视化结果
function visualizeResults(bestSolution, cities, demandPoints, params, convergenceCurve)
% 数据完整性检查
if length(cities) < 5
warning('城市数据不完整,只找到%d个城市', length(cities));
end
figure;
hold on;
% 绘制城市
cityCoords = [[cities.x]; [cities.y]]';
scatter(cityCoords(:,1), cityCoords(:,2), 100, 'r', 'filled', 'DisplayName','城市');
for i = 1:length(cities)
text(cities(i).x+2, cities(i).y, cities(i).name, 'FontSize',10);
end
% 绘制需求点(大小与需求量成正比)
demandSizes = [demandPoints.demand] * 20;
scatter([demandPoints.x], [demandPoints.y], demandSizes, 'b', 'filled', 'DisplayName','需求点');
% 绘制加氢站
stationCoords = reshape(bestSolution, 2, [])';
scatter(stationCoords(:,1), stationCoords(:,2), 150, 'g', 'p', 'filled', 'DisplayName','加氢站');
% 绘制服务范围
theta = linspace(0,2*pi,100);
for i = 1:size(stationCoords,1)
x = stationCoords(i,1) + params.serviceRadius*cos(theta);
y = stationCoords(i,2) + params.serviceRadius*sin(theta);
plot(x, y, 'g--', 'LineWidth',0.5);
end
% 收敛曲线
figure;
plot(convergenceCurve, 'LineWidth',2);
xlabel('迭代次数');
ylabel('适应度值');
title('收敛曲线');
grid on;
legend('show');
axis equal;
hold off;
end
数据
1. 城市信息表(5个)
城市名称 | X坐标 (km) | Y坐标 (km) |
---|---|---|
北京 | 10 | 50 |
天津 | 30 | 60 |
济南 | 60 | 40 |
青岛 | 80 | 70 |
石家庄 | 40 | 20 |
2. 需求点信息表(20个)
X坐标 (km) | Y坐标 (km) | 需求量 (吨/天) |
---|---|---|
15 | 45 | 3.2 |
25 | 55 | 2.8 |
40 | 50 | 4.1 |
50 | 60 | 3.5 |
65 | 65 | 2.9 |
75 | 55 | 3.8 |
20 | 30 | 2.5 |
30 | 35 | 3.0 |
45 | 30 | 3.7 |
55 | 25 | 4.2 |
60 | 35 | 3.3 |
70 | 25 | 2.7 |
10 | 20 | 1.8 |
25 | 15 | 2.3 |
35 | 10 | 2.0 |
45 | 15 | 2.6 |
55 | 10 | 3.1 |
65 | 15 | 2.4 |
75 | 20 | 3.9 |
80 | 30 | 2.2 |
3. 算法参数表
参数名称 | 参数值 |
---|---|
要建设的加氢站数量 | 2 |
种群大小 | 100 |
最大迭代次数 | 500 |
成本权重 | 1.0 |
覆盖率权重 | 100.0 |
距离权重 | 1.0 |
服务半径 (km) | 20.0 |
表格说明:
- 城市坐标代表候选站点的潜在建设位置
- 需求点坐标和需求量表示需要服务的交通走廊节点
- 权重参数用于多目标优化(建设成本、覆盖率、服务距离)的平衡
4. 文件
hydrogen_station_input.inp(和main_hydrogen_station_location.m放在同一文件夹下)
5
北京 10 50
天津 30 60
济南 60 40
青岛 80 70
石家庄 40 20
20
15 45 3.2
25 55 2.8
40 50 4.1
50 60 3.5
65 65 2.9
75 55 3.8
20 30 2.5
30 35 3.0
45 30 3.7
55 25 4.2
60 35 3.3
70 25 2.7
10 20 1.8
25 15 2.3
35 10 2.0
45 15 2.6
55 10 3.1
65 15 2.4
75 20 3.9
80 30 2.2
2
100
500
1.0
100.0
1.0
20.0