编辑代码

import numpy as np
import matplotlib.pyplot as plt
from SALib.sample import saltelli
from SALib.analyze import sobol

# ================== 1. 树形网络节点定义 ==================
class ProductionNode:
    def __init__(self, name, node_type, defect_rate, cost, inspect_cost=None, disassemble_cost=None):
        self.name = name
        self.type = node_type  # 'part', 'semi', 'final'
        self.defect_rate = defect_rate  # 次品率
        self.cost = cost  # 购买单价或装配成本
        self.inspect_cost = inspect_cost  # 检测成本
        self.disassemble_cost = disassemble_cost  # 拆解费用
        self.children = []  # 子节点
        self.parent = None  # 父节点（用于拆解溯源）

    def add_child(self, child):
        self.children.append(child)
        child.parent = self

# ================== 2. 树形网络构建函数 ==================
def build_complex_network():
    """根据表2构建两道工序、8零配件的树形网络"""
    # 零配件节点（8个）
    parts = [
        ProductionNode("Part1", "part", 0.10, 2, 1),
        ProductionNode("Part2", "part", 0.10, 8, 1),
        ProductionNode("Part3", "part", 0.10, 12, 2),
        ProductionNode("Part4", "part", 0.10, 2, 1),
        ProductionNode("Part5", "part", 0.10, 8, 1),
        ProductionNode("Part6", "part", 0.10, 12, 2),
        ProductionNode("Part7", "part", 0.10, 8, 1),
        ProductionNode("Part8", "part", 0.10, 12, 2)
    ]
    
    # 半成品节点（3个）
    semi1 = ProductionNode("Semi1", "semi", 0.10, 8, 4, 6)
    semi1.add_child(parts[0])
    semi1.add_child(parts[1])
    
    semi2 = ProductionNode("Semi2", "semi", 0.10, 8, 4, 6)
    semi2.add_child(parts[2])
    semi2.add_child(parts[3])
    
    semi3 = ProductionNode("Semi3", "semi", 0.10, 6, 10, 6)
    semi3.add_child(parts[4])
    semi3.add_child(parts[5])
    semi3.add_child(parts[6])
    semi3.add_child(parts[7])
    
    # 成品节点
    final = ProductionNode("Final", "final", 0.0, 0, None, None)
    final.add_child(semi1)
    final.add_child(semi2)
    final.add_child(semi3)
    
    return final

# ================== 3. 剩余价值计算 ==================
def calculate_residual_value(node, inspect_decisions, disassemble_decisions, depth=0):
    """递归计算节点及其子树的剩余价值（带缩进打印）"""
    indent = "  " * depth
    print(f"{indent}计算节点 [{node.name}] 类型={node.type}")
    
    if node.type == 'part':
        # 零配件剩余价值
        if inspect_decisions.get(node.name, False):
            # 检测后剩余正品价值
            residual = (1 - node.defect_rate) * node.cost - node.inspect_cost
            print(f"{indent}  - 检测后剩余价值: {residual:.2f}")
        else:
            # 不检测则保留全部价值（含次品风险）
            residual = node.cost
            print(f"{indent}  - 不检测剩余价值: {residual:.2f}")
        return residual
    
    elif node.type == 'semi':
        # 半成品剩余价值 = 拆解价值 - 拆解费用
        if disassemble_decisions.get(node.name, False):
            child_values = [calculate_residual_value(child, inspect_decisions, disassemble_decisions, depth+1) 
                          for child in node.children]
            total = sum(child_values) - node.disassemble_cost
            print(f"{indent}  - 拆解后总价值: {total:.2f}")
            return total
        else:
            print(f"{indent}  - 不拆解，剩余价值=0")
            return 0
    
    elif node.type == 'final':
        # 成品无剩余价值（已售出或调换）
        return 0

# ================== 4. 遗传算法完整实现 ==================
class GeneticOptimizer:
    def __init__(self, pop_size=50, elite_size=5, mutation_rate=0.1, generations=100):
        self.pop_size = pop_size
        self.elite_size = elite_size
        self.mutation_rate = mutation_rate
        self.generations = generations
        
    def initialize_population(self, decision_length=12):
        """生成随机初始种群（二进制编码）"""
        return [np.random.randint(0, 2, decision_length).tolist() for _ in range(self.pop_size)]
    
    def decode_decisions(self, individual):
        """解析基因编码为检测与拆解决策字典"""
        # 编码结构: [零配件1-8检测, 半成品1-3检测, 成品拆解]
        inspect_parts = {f"Part{i+1}": bool(individual[i]) for i in range(8)}
        inspect_semi = {f"Semi{i+1}": bool(individual[8+i]) for i in range(3)}
        disassemble_final = {"Final": bool(individual[11])}
        return {**inspect_parts, **inspect_semi, **disassemble_final}
    
    def calculate_total_cost(self, root_node, inspect_decisions):
        """计算总成本（包括检测成本）"""
        total_cost = 0
        for part in root_node.children[0].children:  # 遍历所有零配件
            if inspect_decisions.get(part.name, False):
                total_cost += part.inspect_cost
        return total_cost
    
    def monte_carlo_risk(self, root_node, inspect_decisions, n_sim=1000):
        """蒙特卡洛风险模拟"""
        profits = []
        for _ in range(n_sim):
            total_cost = self.calculate_total_cost(root_node, inspect_decisions)
            residual_value = calculate_residual_value(root_node, inspect_decisions, {})
            profit = (1000 * 200) - total_cost + residual_value  # 假设生产1000个成品
            profits.append(profit)
        return np.mean(profits), np.std(profits)
    
    def fitness(self, individual, root_node):
        """计算适应度：利润 + 0.3*剩余价值 - 0.2*风险标准差"""
        decisions = self.decode_decisions(individual)
        inspect_decisions = {k: v for k, v in decisions.items() if k.startswith('Part') or k.startswith('Semi')}
        disassemble_decisions = {k: v for k, v in decisions.items() if k == 'Final'}
        
        # 计算总成本与利润
        total_cost = self.calculate_total_cost(root_node, inspect_decisions)
        residual_value = calculate_residual_value(root_node, inspect_decisions, disassemble_decisions)
        profit = (1000 * 200) - total_cost + residual_value  # 假设生产1000个成品
        
        # 蒙特卡洛风险模拟
        _, profit_std = self.monte_carlo_risk(root_node, inspect_decisions)
        return profit + 0.3 * residual_value - 0.2 * profit_std
    
    def roulette_selection(self, population):
        """轮盘赌选择"""
        total_fitness = sum(ind[0] for ind in population)
        selection_probs = [ind[0] / total_fitness for ind in population]
        return np.random.choice([ind[1] for ind in population], size=self.pop_size - self.elite_size, p=selection_probs)
    
    def crossover(self, parent1, parent2):
        """单点交叉"""
        crossover_point = np.random.randint(1, len(parent1))
        return parent1[:crossover_point] + parent2[crossover_point:]
    
    def mutate(self, individual):
        """随机位变异"""
        for i in range(len(individual