import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
import random
import math
from collections import defaultdict
# 三种地方的位置
file_path = "C:\\Users\\Dd\\Desktop\\工作簿1.xlsx"
sheet_name = 'Sheet3'
df = pd.read_excel(file_path, sheet_name=sheet_name)
# 提取坐标数据并分组
a = [tuple(row) for row in df[["栅格x坐标", "栅格y坐标", "高程(米)"]].values]
location_z = a[0:9]
name_z = ['Z1', 'Z2', 'Z3', 'Z4', 'Z5', 'Z6', 'Z7', 'Z8', 'Z9']
location_c = a[9:18]
name_c = ['C1', 'C2', 'C3', 'C4', 'C5', 'C6', 'C7', 'C8', 'C9']
location_h = [a[18]]
name_h = ['H']
# 参数设置
POPULATION_SIZE = 50 # 种族数
GENERATION = 100 # 遗传代数
MUTATION_RATE = 0.1 # 变异概率
ELITE_RATE = 0.2 # 精英比例
CROSSOVER_RATE = 0.8 # 交叉概率
one_day = 7.5 # 人均一天7.5kg
# 阵地驻扎人数
z_people = [68, 60, 38, 60, 90, 60, 38, 90, 60]
# 阵地初始物资
z_start = [1020.0, 900.0, 570.0, 900.0, 1350.0, 900.0, 570.0, 1350.0, 900.0]
# a的装卸时间
a_time = 0.5 # 装卸时间(小时)
# b的装卸时间
b_time = 1 # 装卸时间(小时)
# 中转仓库物资容量上限
c_max = [600, 800, 800, 800, 800, 800, 600, 800, 500]
# 飞行速度
a_v = 10 # km/h
b_v = 15 # km/h
# 飞机最大载重
a_m = 300 # kg
b_m = 500 # kg
# 飞机每公里消耗电量百分比
a_consumption = 0.01 # 1%/km
b_consumption = 0.01 # 1%/km
# 模拟天数
DAYS = 15
# 每天工作时间(小时)
WORK_HOURS = 24
# 计算距离函数
def calculate_distance(loc1, loc2):
return math.sqrt((loc1[0] - loc2[0]) ** 2 + (loc1[1] - loc2[1]) ** 2 + ((loc1[2] - loc2[2]) / 5) ** 2)
# 计算飞行时间
def flight_time(distance, speed):
return distance / speed
# 计算电量消耗
def battery_consumption(distance, consumption_rate):
return distance * consumption_rate
# 充电函数
def charge_battery(current_battery, hours):
if hours >= 2:
if current_battery < 0.2:
return 0.2
elif current_battery < 0.4:
return 0.4
elif current_battery < 0.6:
return 0.6
elif current_battery < 0.8:
return 0.8
else:
return 1.0
return current_battery
# 初始化距离矩阵
distances_c_to_h = [calculate_distance(loc, location_h[0]) for loc in location_c]
distances_c_to_z = [[calculate_distance(c_loc, z_loc) for z_loc in location_z] for c_loc in location_c]
distances_z_to_z = [[calculate_distance(z_loc1, z_loc2) for z_loc2 in location_z] for z_loc1 in location_z]
# 个体表示:仓库启用状态 + B飞机位置 + A飞机位置
def create_individual():
# 仓库启用状态 (9个仓库)
ck_active = [random.choice([0, 1]) for _ in range(9)]
# 确保至少有一个仓库启用
if sum(ck_active) == 0:
ck_active[random.randint(0, 8)] = 1
# B飞机位置 (10-50架)
b_count = random.randint(1, 50)
b_locations = [random.choice(name_c) for _ in range(b_count)]
# A飞机位置 (10-50架)
a_count = random.randint(1, 50)
a_locations = [random.choice(name_c) for _ in range(a_count)]
return [ck_active, b_locations, a_locations]
# 创建初始种群
def initialize_population():
return [create_individual() for _ in range(POPULATION_SIZE)]
# 适应度函数
def fitness(individual):
ck_active, b_locations, a_locations = individual
penalty = 0
# 1. 仓库启用惩罚
active_warehouses = sum(ck_active)
penalty += active_warehouses * 1000 # 启用仓库惩罚
# 2. 飞机数量惩罚
total_planes = len(b_locations) + len(a_locations)
penalty += total_planes * 100 # 飞机数量惩罚
# 3. 模拟15天物流
# 初始化库存
z_inventory = z_start.copy()
c_inventory = [0] * 9
# 飞机状态,位置,电量,结束时间
a_planes = {i: {'location': loc, 'battery': 1.0, 'busy_until': 0}
for i, loc in enumerate(a_locations)}
b_planes = {i: {'location': loc, 'battery': 1.0, 'busy_until': 0}
for i, loc in enumerate(b_locations)}
# 每天模拟
for day in range(DAYS):
# 消耗阵地物资
for i in range(9):
z_inventory[i] -= z_people[i] * one_day
if z_inventory[i] < 0:
penalty += abs(z_inventory[i]) * 10000 # 物资短缺惩罚
# 模拟24小时
for hour in range(WORK_HOURS):
current_time = day * WORK_HOURS + hour
# B飞机任务调度 (从H到C)
for plane_id, plane in b_planes.items():
if current_time >= plane['busy_until']:
c_index = name_c.index(plane['location'])
# 检查仓库是否启用
if ck_active[c_index]:
# 计算飞行时间
distance = distances_c_to_h[c_index]
flight_t = flight_time(distance, b_v)
total_t = flight_t * 2 + b_time # 往返时间+装卸时间
# 检查电量
consumption = battery_consumption(distance * 2, b_consumption)
if plane['battery'] >= consumption:
# 执行飞行
plane['busy_until'] = current_time + total_t
plane['battery'] -= consumption
# 运输物资
capacity = min(b_m, c_max[c_index] - c_inventory[c_index])
if capacity > 0:
c_inventory[c_index] += capacity
else:
# 充电
plane['busy_until'] = current_time + 1
plane['battery'] = charge_battery(plane['battery'], 1)
# A飞机任务调度 (从C到Z)
for plane_id, plane in a_planes.items():
if current_time >= plane['busy_until']:
c_index = name_c.index(plane['location'])
# 检查仓库是否启用
if ck_active[c_index] and c_inventory[c_index] > 0:
# 找到最近的Z点
z_index = np.argmin(distances_c_to_z[c_index])
distance = distances_c_to_z[c_index][z_index]
flight_t = flight_time(distance, a_v)
total_t = flight_t * 2 + a_time * 2 # 往返时间+装卸时间
# 检查电量
consumption = battery_consumption(distance * 2, a_consumption)
if plane['battery'] >= consumption:
# 执行飞行
plane['busy_until'] = current_time + total_t
plane['battery'] -= consumption
# 运输物资
capacity = min(a_m, c_inventory[c_index],
z_inventory[z_index] + z_people[z_index] * one_day)
if capacity > 0:
c_inventory[c_index] -= capacity
z_inventory[z_index] += capacity
else:
# 充电
plane['busy_until'] = current_time + 1
plane['battery'] = charge_battery(plane['battery'], 1)
# 检查仓库容量
for i in range(9):
if c_inventory[i] > c_max[i]:
penalty += (c_inventory[i] - c_max[i]) * 100 # 仓库超容惩罚
# 最终检查阵地物资
for i in range(9):
if z_inventory[i] < 0:
penalty += abs(z_inventory[i]) * 10000 # 物资短缺惩罚
# 适应度 = 惩罚的倒数 (最小化惩罚)
return 1 / (1 + penalty)
# 选择操作 (锦标赛选择)
def tournament_selection(population, fitnesses, tournament_size=3):
selected = []
for _ in range(len(population)):
tournament = random.sample(list(zip(population, fitnesses)), tournament_size)
winner = max(tournament, key=lambda x: x[1])[0]
selected.append(winner)
return selected
# 交叉操作
def crossover(parent1, parent2):
if random.random() < CROSSOVER_RATE:
# 仓库启用状态交叉
ck1, b1, a1 = parent1
ck2, b2, a2 = parent2
# 单点交叉
crossover_point = random.randint(1, 8)
new_ck = ck1[:crossover_point] + ck2[crossover_point:]
# 飞机位置交叉
min_b = min(len(b1), len(b2))
if min_b > 0:
b_crossover = random.randint(1, min_b - 1)
new_b = b1[:b_crossover] + b2[b_crossover:]
else:
new_b = b1 + b2
min_a = min(len(a1), len(a2))
if min_a > 0:
a_crossover = random.randint(1, min_a - 1)
new_a = a1[:a_crossover] + a2[a_crossover:]
else:
new_a = a1 + a2
return [new_ck, new_b, new_a]
return parent1
# 变异操作
def mutate(individual):
ck_active, b_locations, a_locations = individual
# 仓库启用状态变异
for i in range(len(ck_active)):
if random.random() < MUTATION_RATE:
ck_active[i] = 1 - ck_active[i]
# 确保至少有一个仓库启用
if sum(ck_active) == 0:
ck_active[random.randint(0, 8)] = 1
# B飞机位置变异
for i in range(len(b_locations)):
if random.random() < MUTATION_RATE:
b_locations[i] = random.choice(name_c)
# A飞机位置变异
for i in range(len(a_locations)):
if random.random() < MUTATION_RATE:
a_locations[i] = random.choice(name_c)
# 随机增加/减少飞机
if random.random() < MUTATION_RATE:
if random.choice([True, False]) and len(b_locations) < 50:
b_locations.append(random.choice(name_c))
elif len(b_locations) > 10:
b_locations.pop(random.randint(0, len(b_locations) - 1))
if random.random() < MUTATION_RATE:
if random.choice([True, False]) and len(a_locations) < 50:
a_locations.append(random.choice(name_c))
elif len(a_locations) > 10:
a_locations.pop(random.randint(0, len(a_locations) - 1))
return [ck_active, b_locations, a_locations]
# 遗传算法主函数
def genetic_algorithm():
population = initialize_population()
best_fitness = -float('inf')
best_individual = None
fitness_history = []
for generation in range(GENERATION):
# 计算适应度
fitnesses = [fitness(ind) for ind in population]
# 更新最佳个体
current_best_fitness = max(fitnesses)
current_best_index = fitnesses.index(current_best_fitness)
if current_best_fitness > best_fitness:
best_fitness = current_best_fitness
best_individual = population[current_best_index]
fitness_history.append(best_fitness)
# 选择
selected = tournament_selection(population, fitnesses)
# 交叉
new_population = []
for i in range(0, len(selected), 2):
if i + 1 < len(selected):
child1 = crossover(selected[i], selected[i + 1])
child2 = crossover(selected[i + 1], selected[i])
new_population.extend([child1, child2])
else:
new_population.append(selected[i])
# 变异
population = [mutate(ind) for ind in new_population]
# 精英保留
elite_size = int(POPULATION_SIZE * ELITE_RATE)
elite_indices = np.argsort(fitnesses)[-elite_size:]
for i in elite_indices:
population.append(population[i])
# 保持种群大小
population = population[:POPULATION_SIZE]
print(f"Generation {generation + 1}, Best Fitness: {best_fitness:.6f}")
return best_individual, best_fitness, fitness_history
# 运行遗传算法
best_solution, best_fitness, fitness_history = genetic_algorithm()
# 打印最佳解决方案
print("\nBest Solution Found:")
print(f"Fitness: {best_fitness}")
print(f"Active Warehouses: {[name_c[i] for i, active in enumerate(best_solution[0]) if active]}")
print(f"Number of B Planes: {len(best_solution[1])}")
print(f"Number of A Planes: {len(best_solution[2])}")
# 可视化适应度变化
plt.figure(figsize=(10, 6))
plt.plot(fitness_history)
plt.title('Genetic Algorithm Fitness Progression')
plt.xlabel('Generation')
plt.ylabel('Fitness')
plt.grid(True)
plt.show()
# 模拟最佳解决方案的详细物流
def simulate_solution(individual):
ck_active, b_locations, a_locations = individual
# 初始化库存
z_inventory = z_start.copy()
c_inventory = [0] * 9
# 飞机状态
a_planes = {i: {'location': loc, 'battery': 1.0, 'busy_until': 0, 'log': []}
for i, loc in enumerate(a_locations)}
b_planes = {i: {'location': loc, 'battery': 1.0, 'busy_until': 0, 'log': []}
for i, loc in enumerate(b_locations)}
# 日志记录
daily_log = []
# 每天模拟
for day in range(DAYS):
day_log = {
'day': day + 1,
'z_inventory': z_inventory.copy(),
'c_inventory': c_inventory.copy(),
'b_flights': 0,
'a_flights': 0,
'b_cargo': 0,
'a_cargo': 0
}
# 消耗阵地物资
for i in range(9):
z_inventory[i] -= z_people[i] * one_day
# 模拟24小时
for hour in range(WORK_HOURS):
current_time = day * WORK_HOURS + hour
# B飞机任务调度 (从H到C)
for plane_id, plane in b_planes.items():
if current_time >= plane['busy_until']:
c_index = name_c.index(plane['location'])
# 检查仓库是否启用
if ck_active[c_index]:
# 计算飞行时间
distance = distances_c_to_h[c_index]
flight_t = flight_time(distance, b_v)
total_t = flight_t * 2 + b_time # 往返时间+装卸时间
# 检查电量
consumption = battery_consumption(distance * 2, b_consumption)
if plane['battery'] >= consumption:
# 执行飞行
plane['busy_until'] = current_time + total_t
plane['battery'] -= consumption
# 运输物资
capacity = min(b_m, c_max[c_index] - c_inventory[c_index])
if capacity > 0:
c_inventory[c_index] += capacity
day_log['b_flights'] += 1
day_log['b_cargo'] += capacity
# 记录日志
plane['log'].append({
'day': day + 1,
'time': hour,
'from': 'H',
'to': name_c[c_index],
'cargo': capacity,
'type': 'B'
})
else:
# 充电
plane['busy_until'] = current_time + 1
plane['battery'] = charge_battery(plane['battery'], 1)
# A飞机任务调度 (从C到Z)
for plane_id, plane in a_planes.items():
if current_time >= plane['busy_until']:
c_index = name_c.index(plane['location'])
# 检查仓库是否启用
if ck_active[c_index] and c_inventory[c_index] > 0:
# 找到最近的Z点
z_index = np.argmin(distances_c_to_z[c_index])
distance = distances_c_to_z[c_index][z_index]
flight_t = flight_time(distance, a_v)
total_t = flight_t * 2 + a_time * 2 # 往返时间+装卸时间
# 检查电量
consumption = battery_consumption(distance * 2, a_consumption)
if plane['battery'] >= consumption:
# 执行飞行
plane['busy_until'] = current_time + total_t
plane['battery'] -= consumption
# 运输物资
capacity = min(a_m, c_inventory[c_index],
z_inventory[z_index] + z_people[z_index] * one_day)
if capacity > 0:
c_inventory[c_index] -= capacity
z_inventory[z_index] += capacity
day_log['a_flights'] += 1
day_log['a_cargo'] += capacity
# 记录日志
plane['log'].append({
'day': day + 1,
'time': hour,
'from': name_c[c_index],
'to': name_z[z_index],
'cargo': capacity,
'type': 'A'
})
else:
# 充电
plane['busy_until'] = current_time + 1
plane['battery'] = charge_battery(plane['battery'], 1)
daily_log.append(day_log)
return daily_log, a_planes, b_planes
# 模拟最佳解决方案
daily_logs, a_planes, b_planes = simulate_solution(best_solution)
# 打印每日日志
print("\nDaily Logistics Summary:")
for log in daily_logs:
print(f"Day {log['day']}:")
print(f" Z Inventory: {[f'{inv:.1f}' for inv in log['z_inventory']]}")
print(f" C Inventory: {[f'{inv:.1f}' for inv in log['c_inventory']]}")
print(f" B Flights: {log['b_flights']}, Cargo: {log['b_cargo']}kg")
print(f" A Flights: {log['a_flights']}, Cargo: {log['a_cargo']}kg")
# 打印飞机飞行日志
print("\nPlane Flight Logs (Sample):")
for plane_id, plane in list(a_planes.items())[:2]: # 只打印前两架A飞机
print(f"\nA Plane {plane_id} Log:")
for entry in plane['log'][:5]: # 只打印前5条记录
print(
f" Day {entry['day']}, Time {entry['time']}: {entry['from']} -> {entry['to']}, Cargo: {entry['cargo']}kg")
for plane_id, plane in list(b_planes.items())[:2]: # 只打印前两架B飞机
print(f"\nB Plane {plane_id} Log:")
for entry in plane['log'][:5]: # 只打印前5条记录
print(
f" Day {entry['day']}, Time {entry['time']}: {entry['from']} -> {entry['to']}, Cargo: {entry['cargo']}kg"))进行修改不是每天计算而是每小时计算是每天的1/24,物资消耗速度是每日按小时等比例进行,前沿征地最多可以储存每日需求量的3.5倍,各前沿征地相邻两次得到补给时间不超过3天,且每次得到补给时,日常物资持有量不低于每日需求量的1/4,且可以从c点到了z点之后如果还有剩余可以给另外的z点进行补充,飞机不用返回原来的起始C点,可以在任意c点进行装和充电