均值优化的超集子集划分方法与实现

1. 问题定义与挑战

我们面临的核心问题是：给定一个包含 M 个元素的超集 S，以及 N 个预设的子集大小 x0, x1, ..., xn-1（其中 sum(x0, ..., xn-1) == M），如何将超集 S 中的所有元素无重复地分配到这 N 个子集中，使得每个子集的均值与超集 S 的均值尽可能接近。我们的目标是最小化所有子集均值与超集均值之间绝对偏差的总和。超集中的元素通常是实数（浮点数），且多为正数。

例如，如果超集 S = {100, 100, 100, 100, 100, 101, ..., 101 (10次), 102, ..., 102 (5次)}，其均值为 101。我们需要创建 3 个子集，大小分别为 2, 4, 14。一个“完美”的分配方案是使得每个子集的均值都为 101。

这个问题的挑战在于其组合爆炸性。随着超集元素数量和子集数量的增加，可能的分配方案呈指数级增长，暴力枚举变得不可行。此外，我们还需要在合理的时间内（例如1秒内）找到一个解决方案，尤其是在子集数量为10-25，超集元素数量可能高达1000-10000个唯一值的情况下。

2. 数学建模：集合划分问题与混合整数线性规划 (MILP)

这个特定的划分问题可以被建模为一个集合划分问题（Set Partitioning Problem），并通过混合整数线性规划（Mixed-Integer Linear Programming, MILP）来求解。MILP是一种优化技术，它允许目标函数和约束条件是线性的，并且部分或全部变量必须是整数。

2.1 变量定义

我们引入二进制决策变量 y_{ij}：

• y_{ij} = 1：如果超集 S 中的第 j 个元素被分配到第 i 个子集。
• y_{ij} = 0：否则。

其中，i 遍历 0 到 N-1（子集索引），j 遍历 0 到 M-1（超集元素索引）。

2.2 目标函数

首先计算超集 S 的总和 Sum_S = sum(S) 和均值 Mean_S = Sum_S / M。 对于每个子集 i，其目标总和应为 TargetSum_i = x_i * Mean_S。 我们的目标是最小化所有子集实际总和与其目标总和之间的绝对偏差之和。 即：Minimizing Sum_{i=0}^{N-1} | (sum_{j=0}^{M-1} y_{ij} * S[j]) - TargetSum_i |

为了将绝对值项线性化，我们引入辅助变量 e_i（表示第 i 个子集的绝对误差），并将目标函数改为： Minimizing Sum_{i=0}^{N-1} e_i

并添加以下约束：

• e_i >= (sum_{j=0}^{M-1} y_{ij} * S[j]) - TargetSum_i
• e_i >= -((sum_{j=0}^{M-1} y_{ij} * S[j]) - TargetSum_i)

2.3 约束条件

1. 子集大小约束： 每个子集 i 必须包含预设的 x_i 个元素。 sum_{j=0}^{M-1} y_{ij} = x_i，对于每个 i = 0, ..., N-1。

2. 元素唯一性约束： 超集 S 中的每个元素 j 必须且只能被分配到一个子集。 sum_{i=0}^{N-1} y_{ij} = 1，对于每个 j = 0, ..., M-1。

3. 使用 PuLP 进行 Python 实现

PuLP 是一个强大的 Python 库，用于建模和解决线性规划问题。它支持多种求解器（如 CBC、GLPK、Gurobi 等）。

下面是使用 PuLP 解决上述问题的示例代码：

from statistics import mean
import pulp

def partition_superset_by_mean(superset_data, set_sizes):
    """
    将超集划分为指定大小的子集，使每个子集的均值尽可能接近超集均值。

    Args:
        superset_data (list): 包含所有元素的超集列表。
        set_sizes (list): 包含每个子集所需元素数量的列表。

    Returns:
        list: 包含划分后子集列表的列表。
    """

    target_sum = sum(superset_data)
    N = len(set_sizes)
    M = len(superset_data)

    # 验证子集大小总和是否等于超集元素总数
    assert sum(set_sizes) == M, "子集大小总和必须等于超集元素总数"

    # 创建 PuLP 问题实例
    set_partitioning_model = pulp.LpProblem("Set_Partitioning_Model", pulp.LpMinimize)

    # 定义决策变量 y_ij
    # covering[s][i] 表示超集中的第 i 个元素是否分配给第 s 个子集
    covering = {}
    for s in range(N):
        vals = []
        for i, v in enumerate(superset_data):
            vals.append(
                pulp.LpVariable(
                    f"assign_set_{s}_element_idx_{i:>02}_val_{v}",
                    lowBound=0,
                    upBound=1,
                    cat=pulp.LpInteger,  # 二进制变量
                )
            )
        covering[s] = vals

    # 定义绝对误差变量 e_s
    abs_sum_errs = []
    for s_i in range(N):
        set_sum_err_abs = pulp.LpVariable(f"set_{s_i}_sum_error_abs")
        abs_sum_errs.append(set_sum_err_abs)

    # OBJECTIVE: 最小化所有子集绝对误差之和
    set_partitioning_model += pulp.lpSum(abs_sum_errs), "Total_Absolute_Error"

    # 添加绝对值线性化约束
    # set_sum_err = (当前子集总和 - 目标子集总和)
    # e_s >= set_sum_err
    # e_s >= -set_sum_err
    superset_mean = target_sum / M
    for s_i, st_vars in covering.items():
        current_set_sum = pulp.lpSum([p * superset_data[i] for i, p in enumerate(st_vars)])
        target_set_sum = set_sizes[s_i] * superset_mean

        # 定义一个中间变量来表示偏差
        set_sum_deviation = pulp.LpVariable(f"set_{s_i}_sum_deviation")
        set_partitioning_model += set_sum_deviation == current_set_sum - target_set_sum, \
                                  f"Deviation_Constraint_Set_{s_i}"

        # 绝对值约束
        set_partitioning_model += abs_sum_errs[s_i] >= set_sum_deviation, \
                                  f"Abs_Error_Positive_Set_{s_i}"
        set_partitioning_model += abs_sum_errs[s_i] >= -set_sum_deviation, \
                                  f"Abs_Error_Negative_Set_{s_i}"

    # 约束1: 子集大小是预设的
    for n, st_vars in zip(set_sizes, covering.values()):
        set_partitioning_model += pulp.lpSum(st_vars) == n, \
                                  f"Set_Size_Constraint_{n}"

    # 约束2: 每个超集元素只能被使用一次
    # zip(*covering.values()) 将所有子集的变量列表转置，以便按元素索引迭代
    for element_idx, element_assignment_vars in enumerate(zip(*covering.values())):
        set_partitioning_model += (
            pulp.lpSum(element_assignment_vars) == 1,
            f"Element_{element_idx}_Used_Once",
        )

    # 求解模型
    set_partitioning_model.solve()

    # 提取结果
    if pulp.LpStatus[set_partitioning_model.status] == 'Optimal':
        result_subsets = []
        print(f"超集均值: {superset_mean}")
        for k, v in covering.items():
            subset_elements = [superset_data[idx] for idx, var in enumerate(v) if var.value() == 1]
            result_subsets.append(subset_elements)
            print(f"子集 {k} ({len(subset_elements)}个元素): {subset_elements}, 均值 = {mean(subset_elements)}")
        return result_subsets
    else:
        print(f"未能找到最优解。状态: {pulp.LpStatus[set_partitioning_model.status]}")
        return None

# 示例 1: 完美分配
print("--- 示例 1: 完美分配 ---")
superset_ex1 = [100]*5 + [101]*10 + [102]*5
set_sizes_ex1 = [2, 4, 14]
partition_superset_by_mean(superset_ex1, set_sizes_ex1)

# 示例 2: 最佳拟合 (完美分配不可能)
print("\n--- 示例 2: 最佳拟合 ---")
superset_ex2 = [100]*5 + [103]*10 + [104]*5
set_sizes_ex2 = [2, 4, 14]
partition_superset_by_mean(superset_ex2, set_sizes_ex2)

示例 1 输出：

--- 示例 1: 完美分配 ---
超集均值: 101.0
子集 0 (2个元素): [101, 101], 均值 = 101
子集 1 (4个元素): [100, 100, 102, 102], 均值 = 101
子集 2 (14个元素): [100, 100, 100, 101, 101, 101, 101, 101, 101, 101, 101, 102, 102, 102], 均值 = 101

示例 2 输出：

--- 示例 2: 最佳拟合 ---
超集均值: 102.5
子集 0 (2个元素): [103, 103], 均值 = 103
子集 1 (4个元素): [100, 100, 104, 104], 均值 = 102
子集 2 (14个元素): [100, 100, 100, 103, 103, 103, 103, 103, 103, 103, 103, 104, 104, 104], 均值 = 102.57142857142857

从输出可以看出，PuLP 成功地为我们找到了最优（或接近最优）的划分方案，使得子集均值尽可能地接近超集均值。

4. 启发式方法与性能考量

虽然 MILP 能够找到最优解，但其计算复杂度较高。对于大规模问题，求解时间可能会很长，尤其当超集元素数量和子集数量都很大时。在实际应用中，如果对求解速度有严格要求，或者问题规模超出 MILP 的有效处理范围，可以考虑使用启发式方法。

4.1 贪婪分配策略

这是一种简单且快速的启发式方法，但不保证全局最优。

• 从小子集开始分配： 优先为最小的子集分配元素，尽量使其均值接近超集均值。然后处理下一个最小的子集，依此类推。
• 预分配与调整： 可以先将超集元素均匀地随机分配到各个子集，以使它们的初始均值接近超集均值。