Files
2026_DesignAI/dofile/examples/02-spatial-intelligence/spatial_optimization.py
T
pengxiao 219232de74 refactor: 重组项目目录结构
以讲义内容为骨架迁移到标准目录格式:
- officefile/ 主内容(12章 + 附录 + CC4SI补充)
- dofile/ 代码示例(11个Python脚本)
- data/ 图片资源
- output/ 生成输出(忽略)
- Archive/ 归档旧目录(忽略)
- .claude/skills/ 保留markdown-to-docx工具链
- .pandoc/ 保留CSL和本地化配置

Co-Authored-By: Claude Opus 4.7 <noreply@anthropic.com>
2026-05-25 14:00:56 +08:00

962 lines
29 KiB
Python

"""
空间优化示例 (Spatial Optimization Example)
==========================================
本示例展示空间智能系统中的空间优化方法。
空间优化是在空间约束下寻找最优解的过程。
核心概念:
1. 目标函数 - 需要最大化或最小化的目标
2. 约束条件 - 空间和非空间限制
3. 决策变量 - 可控制的空间变量
4. 优化算法 - 求解最优解的方法
5. 帕累托前沿 - 多目标优化中的解集
应用场景:
- 选址优化
- 路径优化
- 空间配置优化
- 资源分配优化
作者: CC4SI 项目组
"""
import math
import random
from typing import List, Dict, Tuple, Optional, Any, Callable
from dataclasses import dataclass, field
from enum import Enum
from abc import ABC, abstractmethod
import heapq
# ============================================================================
# 基础数据结构
# ============================================================================
@dataclass
class Point:
"""二维点"""
x: float
y: float
def distance_to(self, other: 'Point') -> float:
return math.sqrt((self.x - other.x)**2 + (self.y - other.y)**2)
def __repr__(self) -> str:
return f"({self.x:.2f}, {self.y:.2f})"
@dataclass
class DemandPoint:
"""需求点"""
id: str
location: Point
demand: float # 需求量
population: int = 0
def __repr__(self) -> str:
return f"Demand({self.id}, demand={self.demand})"
@dataclass
class Facility:
"""设施"""
id: str
location: Point
capacity: float # 服务能力
fixed_cost: float = 0 # 固定成本
variable_cost: float = 0 # 单位可变成本
def __repr__(self) -> str:
return f"Facility({self.id}, loc={self.location})"
@dataclass
class OptimizationResult:
"""优化结果"""
success: bool
objective_value: float
solution: Any
iterations: int = 0
convergence_history: List[float] = field(default_factory=list)
metadata: Dict[str, Any] = field(default_factory=dict)
# ============================================================================
# 约束条件
# ============================================================================
class ConstraintType(Enum):
"""约束类型"""
EQUALITY = "equality" # 等式约束
INEQUALITY = "inequality" # 不等式约束
BOUNDS = "bounds" # 边界约束
@dataclass
class Constraint:
"""约束条件"""
name: str
constraint_type: ConstraintType
rhs: float # 右端值
lhs_function: Optional[Callable[[Any], float]] = None # 左端函数
def is_satisfied(self, variables: Any, tolerance: float = 1e-6) -> bool:
"""检查约束是否满足"""
if self.lhs_function is None:
return True
lhs = self.lhs_function(variables)
if self.constraint_type == ConstraintType.EQUALITY:
return abs(lhs - self.rhs) < tolerance
elif self.constraint_type == ConstraintType.INEQUALITY:
return lhs <= self.rhs + tolerance
return True
# ============================================================================
# 优化问题定义
# ============================================================================
class OptimizationProblem(ABC):
"""优化问题抽象基类"""
def __init__(self, name: str = ""):
self.name = name
self.constraints: List[Constraint] = []
self.objective_calls = 0
@abstractmethod
def objective(self, variables: Any) -> float:
"""目标函数"""
pass
@abstractmethod
def get_initial_solution(self) -> Any:
"""获取初始解"""
pass
def add_constraint(self, constraint: Constraint) -> None:
"""添加约束"""
self.constraints.append(constraint)
def is_feasible(self, variables: Any) -> bool:
"""检查解是否可行"""
return all(c.is_satisfied(variables) for c in self.constraints)
# ============================================================================
# 贪心算法
# ============================================================================
class GreedyOptimizer:
"""
贪心优化器
每一步选择当前最优的选项。
"""
def __init__(self, problem: OptimizationProblem):
self.problem = problem
def optimize(self, max_iterations: int = 1000) -> OptimizationResult:
"""
执行贪心优化
Args:
max_iterations: 最大迭代次数
Returns:
优化结果
"""
current_solution = self.problem.get_initial_solution()
current_value = self.problem.objective(current_solution)
history = [current_value]
for iteration in range(max_iterations):
# 生成邻居解
neighbors = self._generate_neighbors(current_solution)
# 找最优邻居
best_neighbor = None
best_neighbor_value = float('inf')
for neighbor in neighbors:
if self.problem.is_feasible(neighbor):
value = self.problem.objective(neighbor)
if value < best_neighbor_value:
best_neighbor = neighbor
best_neighbor_value = value
# 如果没有改进,停止
if best_neighbor is None or best_neighbor_value >= current_value:
break
current_solution = best_neighbor
current_value = best_neighbor_value
history.append(current_value)
return OptimizationResult(
success=True,
objective_value=current_value,
solution=current_solution,
iterations=len(history),
convergence_history=history
)
def _generate_neighbors(self, solution: Any) -> List[Any]:
"""生成邻居解 (需要根据具体问题实现)"""
return []
# ============================================================================
# 模拟退火算法
# ============================================================================
class SimulatedAnnealing:
"""
模拟退火算法
一种概率性全局优化算法,能够跳出局部最优。
"""
def __init__(self, problem: OptimizationProblem):
self.problem = problem
def optimize(self,
initial_temp: float = 1000.0,
cooling_rate: float = 0.95,
min_temp: float = 0.01,
max_iterations: int = 10000) -> OptimizationResult:
"""
执行模拟退火优化
Args:
initial_temp: 初始温度
cooling_rate: 冷却速率
min_temp: 最小温度
max_iterations: 最大迭代次数
Returns:
优化结果
"""
current_solution = self.problem.get_initial_solution()
current_value = self.problem.objective(current_solution)
best_solution = current_solution
best_value = current_value
temperature = initial_temp
history = [current_value]
iteration = 0
while temperature > min_temp and iteration < max_iterations:
# 生成邻居解
neighbor = self._generate_neighbor(current_solution)
if self.problem.is_feasible(neighbor):
neighbor_value = self.problem.objective(neighbor)
# 决定是否接受新解
delta = neighbor_value - current_value
if delta < 0 or random.random() < math.exp(-delta / temperature):
current_solution = neighbor
current_value = neighbor_value
# 更新最优解
if current_value < best_value:
best_solution = current_solution
best_value = current_value
history.append(best_value)
temperature *= cooling_rate
iteration += 1
return OptimizationResult(
success=True,
objective_value=best_value,
solution=best_solution,
iterations=iteration,
convergence_history=history
)
def _generate_neighbor(self, solution: Any) -> Any:
"""生成邻居解 (需要根据具体问题实现)"""
return solution
# ============================================================================
# 遗传算法
# =============================================================================
class GeneticAlgorithm:
"""
遗传算法
模拟自然进化的全局优化算法。
"""
def __init__(self, problem: OptimizationProblem):
self.problem = problem
def optimize(self,
population_size: int = 50,
generations: int = 100,
mutation_rate: float = 0.1,
crossover_rate: float = 0.8,
elitism_count: int = 2) -> OptimizationResult:
"""
执行遗传算法优化
Args:
population_size: 种群大小
generations: 迭代代数
mutation_rate: 变异率
crossover_rate: 交叉率
elitism_count: 精英保留数量
Returns:
优化结果
"""
# 初始化种群
population = self._initialize_population(population_size)
history = []
for generation in range(generations):
# 评估适应度
fitness = []
for individual in population:
value = self.problem.objective(individual)
fitness.append(1.0 / (1.0 + value)) # 转换为适应度 (越小越好 -> 越大越好)
# 记录最优值
best_idx = max(range(len(fitness)), key=lambda i: fitness[i])
best_value = self.problem.objective(population[best_idx])
history.append(best_value)
# 选择
selected = self._selection(population, fitness)
# 交叉
offspring = self._crossover(selected, crossover_rate)
# 变异
offspring = self._mutation(offspring, mutation_rate)
# 精英保留
if elitism_count > 0:
elite_indices = sorted(range(len(fitness)),
key=lambda i: fitness[i], reverse=True)[:elitism_count]
for i, idx in enumerate(elite_indices):
offspring[i] = population[idx]
population = offspring
# 返回最优解
final_fitness = [self.problem.objective(ind) for ind in population]
best_idx = min(range(len(final_fitness)), key=lambda i: final_fitness[i])
return OptimizationResult(
success=True,
objective_value=final_fitness[best_idx],
solution=population[best_idx],
iterations=generations,
convergence_history=history
)
def _initialize_population(self, size: int) -> List[Any]:
"""初始化种群"""
return [self.problem.get_initial_solution() for _ in range(size)]
def _selection(self, population: List[Any], fitness: List[float]) -> List[Any]:
"""锦标赛选择"""
selected = []
tournament_size = max(3, len(population) // 10)
for _ in range(len(population)):
# 随机选择tournament_size个个体
contestants = random.sample(list(zip(population, fitness)), tournament_size)
# 选择适应度最高的
winner = max(contestants, key=lambda x: x[1])[0]
selected.append(winner)
return selected
def _crossover(self, population: List[Any], rate: float) -> List[Any]:
"""交叉操作"""
offspring = []
for i in range(0, len(population), 2):
parent1 = population[i]
parent2 = population[i + 1] if i + 1 < len(population) else population[0]
if random.random() < rate:
child1, child2 = self._crossover_operators(parent1, parent2)
else:
child1, child2 = parent1, parent2
offspring.extend([child1, child2])
return offspring[:len(population)]
def _crossover_operators(self, parent1: Any, parent2: Any) -> Tuple[Any, Any]:
"""交叉算子 (需要根据具体问题实现)"""
return parent1, parent2
def _mutation(self, population: List[Any], rate: float) -> List[Any]:
"""变异操作"""
mutated = []
for individual in population:
if random.random() < rate:
mutated.append(self._mutate(individual))
else:
mutated.append(individual)
return mutated
def _mutate(self, individual: Any) -> Any:
"""变异算子 (需要根据具体问题实现)"""
return individual
# ============================================================================
# 选址优化问题
# ============================================================================
class LocationProblem(OptimizationProblem):
"""
设施选址问题
在给定候选位置中选择最优的设施位置组合。
"""
def __init__(self,
demand_points: List[DemandPoint],
candidate_locations: List[Point],
num_facilities: int,
fixed_costs: List[float] = None,
transportation_cost: float = 1.0):
"""
初始化选址问题
Args:
demand_points: 需求点列表
candidate_locations: 候选位置列表
num_facilities: 设施数量
fixed_costs: 各候选位置的固定成本
transportation_cost: 单位运输成本
"""
super().__init__("设施选址问题")
self.demand_points = demand_points
self.candidate_locations = candidate_locations
self.num_facilities = num_facilities
self.fixed_costs = fixed_costs or [0] * len(candidate_locations)
self.transportation_cost = transportation_cost
def objective(self, solution: List[int]) -> float:
"""
计算目标函数值
Args:
solution: 选中的候选位置索引列表
Returns:
总成本
"""
self.objective_calls += 1
if not solution or len(solution) != self.num_facilities:
return float('inf')
total_cost = 0.0
# 固定成本
for idx in solution:
if 0 <= idx < len(self.fixed_costs):
total_cost += self.fixed_costs[idx]
# 运输成本 (每个需求点分配到最近的设施)
for demand in self.demand_points:
min_dist = float('inf')
for facility_idx in solution:
if 0 <= facility_idx < len(self.candidate_locations):
facility_loc = self.candidate_locations[facility_idx]
dist = demand.location.distance_to(facility_loc)
min_dist = min(min_dist, dist)
total_cost += min_dist * demand.demand * self.transportation_cost
return total_cost
def get_initial_solution(self) -> List[int]:
"""获取初始解 (随机选择)"""
n = len(self.candidate_locations)
if n <= self.num_facilities:
return list(range(n))
return random.sample(range(n), self.num_facilities)
class LocationSAOptimizer(SimulatedAnnealing):
"""针对选址问题的模拟退火优化器"""
def _generate_neighbor(self, solution: List[int]) -> List[int]:
"""生成邻居解"""
if not solution:
return solution
neighbor = solution.copy()
n = len(self.problem.candidate_locations)
# 随机选择一个操作
operation = random.choice(['replace', 'swap'])
if operation == 'replace':
# 替换一个设施
idx = random.randint(0, len(neighbor) - 1)
available = [i for i in range(n) if i not in neighbor]
if available:
neighbor[idx] = random.choice(available)
elif operation == 'swap' and len(neighbor) >= 1:
# 交换一个设施
idx = random.randint(0, len(neighbor) - 1)
available = [i for i in range(n) if i not in neighbor]
if available:
neighbor[idx] = random.choice(available)
return neighbor
# ============================================================================
# P-中值问题
# ============================================================================
class PMedianProblem(LocationProblem):
"""
P-中值问题
选择P个设施位置,使需求点到最近设施的总加权距离最小。
"""
def __init__(self,
demand_points: List[DemandPoint],
candidate_locations: List[Point],
p: int):
super().__init__(demand_points, candidate_locations, p, [0] * len(candidate_locations))
# ============================================================================
# 覆盖问题
# ============================================================================
class MaxCoverageProblem(OptimizationProblem):
"""
最大覆盖问题
在给定设施数量限制下,最大化覆盖的需求量。
"""
def __init__(self,
demand_points: List[DemandPoint],
candidate_locations: List[Point],
num_facilities: int,
coverage_radius: float):
"""
初始化最大覆盖问题
Args:
demand_points: 需求点列表
candidate_locations: 候选位置列表
num_facilities: 设施数量
coverage_radius: 覆盖半径
"""
super().__init__("最大覆盖问题")
self.demand_points = demand_points
self.candidate_locations = candidate_locations
self.num_facilities = num_facilities
self.coverage_radius = coverage_radius
# 预计算覆盖关系
self.coverage_matrix = self._compute_coverage()
def _compute_coverage(self) -> List[List[bool]]:
"""计算覆盖矩阵"""
matrix = []
for loc in self.candidate_locations:
coverage = []
for demand in self.demand_points:
covered = loc.distance_to(demand.location) <= self.coverage_radius
coverage.append(covered)
matrix.append(coverage)
return matrix
def objective(self, solution: List[int]) -> float:
"""
计算覆盖的需求量 (负值,因为算法最小化)
Args:
solution: 选中的候选位置索引列表
Returns:
负的覆盖需求量
"""
covered = [False] * len(self.demand_points)
for facility_idx in solution:
if 0 <= facility_idx < len(self.coverage_matrix):
for j, is_covered in enumerate(self.coverage_matrix[facility_idx]):
if is_covered:
covered[j] = True
total_demand = sum(
self.demand_points[j].demand
for j, c in enumerate(covered) if c
)
return -total_demand # 负值用于最小化
def get_initial_solution(self) -> List[int]:
"""获取初始解"""
n = len(self.candidate_locations)
if n <= self.num_facilities:
return list(range(n))
return random.sample(range(n), self.num_facilities)
class MaxCoverageGAOptimizer(GeneticAlgorithm):
"""针对最大覆盖问题的遗传算法优化器"""
def _crossover_operators(self, parent1: List[int], parent2: List[int]) -> Tuple[List[int], List[int]]:
"""单点交叉"""
if not parent1 or not parent2:
return parent1, parent2
size = min(len(parent1), len(parent2))
if size < 2:
return parent1, parent2
point = random.randint(1, size - 1)
child1 = parent1[:point] + [x for x in parent2[point:] if x not in parent1[:point]]
child2 = parent2[:point] + [x for x in parent1[point:] if x not in parent2[:point]]
# 补足长度
all_indices = set(range(len(self.problem.candidate_locations)))
while len(child1) < len(parent1):
available = list(all_indices - set(child1))
if available:
child1.append(random.choice(available))
else:
break
while len(child2) < len(parent2):
available = list(all_indices - set(child2))
if available:
child2.append(random.choice(available))
else:
break
return child1[:len(parent1)], child2[:len(parent2)]
def _mutate(self, individual: List[int]) -> List[int]:
"""变异操作"""
if not individual:
return individual
mutated = individual.copy()
n = len(self.problem.candidate_locations)
# 随机替换一个基因
idx = random.randint(0, len(mutated) - 1)
available = [i for i in range(n) if i not in mutated]
if available:
mutated[idx] = random.choice(available)
return mutated
# ============================================================================
# 路径优化 (TSP)
# ============================================================================
class TSPProblem(OptimizationProblem):
"""
旅行商问题 (TSP)
寻找访问所有城市的最短路径。
"""
def __init__(self, cities: List[Point]):
"""
初始化TSP问题
Args:
cities: 城市位置列表
"""
super().__init__("旅行商问题")
self.cities = cities
self.n = len(cities)
def objective(self, solution: List[int]) -> float:
"""
计算路径总长度
Args:
solution: 城市访问顺序列表
Returns:
路径总长度
"""
if not solution or len(solution) != self.n:
return float('inf')
total = 0.0
for i in range(len(solution)):
from_idx = solution[i]
to_idx = solution[(i + 1) % len(solution)]
total += self.cities[from_idx].distance_to(self.cities[to_idx])
return total
def get_initial_solution(self) -> List[int]:
"""获取初始解 (随机顺序)"""
solution = list(range(self.n))
random.shuffle(solution)
return solution
class TSPSAOptimizer(SimulatedAnnealing):
"""针对TSP的模拟退火优化器"""
def _generate_neighbor(self, solution: List[int]) -> List[int]:
"""生成邻居解 (2-opt交换)"""
if len(solution) < 2:
return solution
neighbor = solution.copy()
# 随机选择两个位置并交换
i, j = random.sample(range(len(neighbor)), 2)
neighbor[i], neighbor[j] = neighbor[j], neighbor[i]
return neighbor
# ============================================================================
# 主程序
# ========================================================================
def main():
"""主程序 - 演示空间优化的使用"""
print("="*70)
print("空间优化示例演示")
print("="*70)
random.seed(42)
# ========================================================================
# 1. P-中值问题 (选址优化)
# ========================================================================
print("\n[部分 1] P-中值问题 - 商场选址优化")
print("-" * 50)
# 生成需求点
demand_points = []
for i in range(20):
demand_points.append(DemandPoint(
id=f"D{i}",
location=Point(random.uniform(0, 100), random.uniform(0, 100)),
demand=random.uniform(50, 200)
))
print(f"\n需求点数量: {len(demand_points)}")
print(f"总需求量: {sum(d.demand for d in demand_points):.1f}")
# 候选位置
candidate_locations = [
Point(20, 20), Point(50, 20), Point(80, 20),
Point(20, 50), Point(50, 50), Point(80, 50),
Point(20, 80), Point(50, 80), Point(80, 80)
]
print(f"\n候选位置数量: {len(candidate_locations)}")
print("候选位置:", [str(loc) for loc in candidate_locations])
# 创建P-中值问题 (选择3个位置)
p_median = PMedianProblem(demand_points, candidate_locations, p=3)
# 使用模拟退火求解
sa_optimizer = LocationSAOptimizer(p_median)
sa_result = sa_optimizer.optimize(
initial_temp=100,
cooling_rate=0.95,
min_temp=0.1,
max_iterations=1000
)
print(f"\n模拟退火结果:")
print(f" 目标函数值: {sa_result.objective_value:.2f}")
print(f" 迭代次数: {sa_result.iterations}")
print(f" 选中的位置: {[str(candidate_locations[i]) for i in sa_result.solution]}")
# 贪心算法对比
greedy_optimizer = GreedyOptimizer(p_median)
# 简单的贪心: 逐步添加最优位置
best_solution = None
best_value = float('inf')
for _ in range(100):
solution = p_median.get_initial_solution()
value = p_median.objective(solution)
if value < best_value:
best_value = value
best_solution = solution
print(f"\n随机搜索对比:")
print(f" 目标函数值: {best_value:.2f}")
print(f" 选中的位置: {[str(candidate_locations[i]) for i in best_solution]}")
# ========================================================================
# 2. 最大覆盖问题
# ========================================================================
print("\n\n[部分 2] 最大覆盖问题 - 5G基站选址")
print("-" * 50)
# 创建更大规模的需求点
coverage_demands = []
for i in range(50):
coverage_demands.append(DemandPoint(
id=f"C{i}",
location=Point(random.uniform(0, 100), random.uniform(0, 100)),
demand=random.uniform(10, 100)
))
coverage_radius = 25
num_bases = 5
max_coverage = MaxCoverageProblem(
coverage_demands,
candidate_locations,
num_bases,
coverage_radius
)
# 使用遗传算法求解
ga_optimizer = MaxCoverageGAOptimizer(max_coverage)
ga_result = ga_optimizer.optimize(
population_size=50,
generations=100,
mutation_rate=0.1
)
covered_demand = -ga_result.objective_value
total_demand = sum(d.demand for d in coverage_demands)
coverage_ratio = covered_demand / total_demand * 100
print(f"\n遗传算法结果:")
print(f" 覆盖需求量: {covered_demand:.1f} / {total_demand:.1f}")
print(f" 覆盖率: {coverage_ratio:.1f}%")
print(f" 迭代次数: {ga_result.iterations}")
print(f" 选中的位置: {[str(candidate_locations[i]) for i in ga_result.solution]}")
# ========================================================================
# 3. 旅行商问题 (TSP)
# ========================================================================
print("\n\n[部分 3] 旅行商问题 - 配送路线优化")
print("-" * 50)
# 生成城市
cities = []
for i in range(15):
cities.append(Point(random.uniform(0, 100), random.uniform(0, 100)))
print(f"\n城市数量: {len(cities)}")
print(f"城市位置: {[str(city) for city in cities[:5]]}...")
tsp = TSPProblem(cities)
tsp_optimizer = TSPSAOptimizer(tsp)
tsp_result = tsp_optimizer.optimize(
initial_temp=1000,
cooling_rate=0.99,
min_temp=0.01,
max_iterations=5000
)
print(f"\n模拟退火结果:")
print(f" 最短路径长度: {tsp_result.objective_value:.2f}")
print(f" 访问顺序: {[cities[i].__repr__() for i in tsp_result.solution[:5]]}...")
# 对比随机解
random_solution = list(range(len(cities)))
random.shuffle(random_solution)
random_length = tsp.objective(random_solution)
print(f"\n随机解对比:")
print(f" 路径长度: {random_length:.2f}")
print(f" 改进: {(1 - tsp_result.objective_value / random_length) * 100:.1f}%")
# ========================================================================
# 4. 多目标优化讨论
# ========================================================================
print("\n\n[部分 4] 多目标优化说明")
print("-" * 50)
print("""
在实际应用中,空间优化往往涉及多个目标:
1. 成本最小化
- 设施建设成本
- 运营成本
- 运输成本
2. 服务最大化
- 覆盖范围
- 服务质量
- 响应时间
3. 公平性
- 服务均等化
- 负载均衡
4. 环境影响
- 最小化污染
- 保护生态
处理方法:
- 加权求和法 (将多目标转为单目标)
- 帕累托优化 (寻找非劣解集)
- 约束法 (将部分目标转为约束)
- 目标规划 (设定目标满意水平)
""")
# ========================================================================
# 5. 收敛过程可视化 (文本)
# ========================================================================
print("\n[部分 5] 优化过程")
print("-" * 50)
if len(ga_result.convergence_history) > 0:
print("\n遗传算法收敛过程:")
steps = min(10, len(ga_result.convergence_history))
step_size = len(ga_result.convergence_history) // steps
for i in range(0, len(ga_result.convergence_history), step_size):
iteration = i
value = -ga_result.convergence_history[i] # 转回正值
print(f" 迭代 {iteration:4d}: 覆盖需求量 = {value:.1f}")
print("\n" + "="*70)
print("演示完成!")
print("="*70)
if __name__ == "__main__":
main()