容量约束设施选址问题(CFLP)¶

类别： 选址

使用 OptAgent 的 Python 接口描述变量、约束与目标。

问题与原始示例来源：Hexaly Code Templates。

问题描述¶

**在容量约束设施选址问题(Capacitated Facility Location Problem, CFLP)**中,若干具有已知需求的客户点必须被分配给可用的设施。分配给每个设施的客户点总需求不得超过其容量。每个设施都有一个固定的开启费用,只要该设施至少服务一个客户点就必须支付。目标是最小化总成本,该成本为所有已开启设施的开启费用与所有客户点的分配费用之和。

学习要点¶

添加集合决策变量以建模分配给每个设施的客户点
使用 lambda 表达式计算每个设施的总需求与总费用

数据¶

所提供的数据文件来自 OR-LIB。数据格式如下:

第一行:候选设施数量与客户点数量
每个设施的容量与开启费用
每个客户点的需求
每个设施与每个客户点之间的分配费用

建模思路¶

容量约束设施选址问题(CFLP)的 OptAgent 模型使用集合决策变量来表示分配给每个设施的客户点。借助 partition 算子,我们确保每个客户点恰好被分配到一个设施。

我们可以使用需求数组上的 at 算子来访问序列中每个客户点的需求。每个设施所服务的总需求通过一个 lambda 函数计算,该函数将 sum 算子应用于所有关联客户点的需求。该总需求必须不超过该设施的容量。

类似地,我们将每个设施的分配费用计算为其所服务的所有客户点的分配价格之和。使用 count 算子,我们检查每个设施是否至少服务一个客户点。如果是,则还需支付该设施的开启费用。

目标函数为所有设施的开启费用与分配费用之和。

Python 实现¶

from pathlib import Path

from optagent import OptModel, solve




def read_data(filename):
    """Read a CFLP instance in the OR-LIB format."""
    with open(filename, encoding="utf-8") as f:
        tokens = f.read().split()
    it = iter(tokens)
    nb_max_facilities = int(next(it))
    nb_sites = int(next(it))
    capacity_data = []
    opening_price_data = []
    for _ in range(nb_max_facilities):
        capacity_data.append(float(next(it)))
        opening_price_data.append(float(next(it)))
    demand_data = [float(next(it)) for _ in range(nb_sites)]
    allocation_price_data = [
        [float(next(it)) for _ in range(nb_sites)] for _ in range(nb_max_facilities)
    ]
    return (
        nb_max_facilities,
        nb_sites,
        capacity_data,
        opening_price_data,
        demand_data,
        allocation_price_data,
    )


def main(instance_file, output_file=None, time_limit=10):
    (
        nb_max_facilities,
        nb_sites,
        capacity_data,
        opening_price_data,
        demand_data,
        allocation_price_data,
    ) = read_data(instance_file)

    model = OptModel()
    # Each facility assignment is the set of sites served by that facility.
    facility_assignments = [model.set(nb_sites) for _ in range(nb_max_facilities)]
    # Every site is served by exactly one facility.
    model.constraint(model.partition(facility_assignments))

    demand = model.array(demand_data)
    allocation_price = model.array(allocation_price_data)

    cost = []
    for f in range(nb_max_facilities):
        facility = facility_assignments[f]
        size = model.count(facility)

        # Capacity constraint: total demand served by this facility.
        demand_lambda = model.lambda_function(lambda i: demand[i])
        model.constraint(model.sum(facility, demand_lambda) <= capacity_data[f])

        # Allocation cost + opening cost when the facility serves at least one site.
        cost_selector = model.lambda_function(
            lambda i: model.at(allocation_price, f, i)
        )
        cost.append(
            model.sum(facility, cost_selector) + opening_price_data[f] * (size > 0)
        )

    total_cost = model.sum(*cost)
    model.minimize(total_cost)

    solution = solve(model, time_limit_s=float(time_limit))
    result_values = {'total_cost': total_cost.value, **{f'facility_{facility}': assignment.value for facility, assignment in enumerate(facility_assignments)}}

    lines = [
        f"Facilities = {nb_max_facilities}; Sites = {nb_sites}; "
        f"Total cost = {result_values['total_cost']}; Status = {solution.feasible}"
    ]
    open_facilities = 0
    for f in range(nb_max_facilities):
        sites = result_values[f"facility_{f}"]
        if sites:
            open_facilities += 1
            lines.append(f"Facility {f} -> sites {sorted(sites)}")
    lines.append(f"Open facilities: {open_facilities}/{nb_max_facilities}")
    result_text = "\n".join(lines)
    print(result_text)
    if output_file is not None:
        Path(output_file).write_text(result_text + "\n", encoding="utf-8")
    return solution

本地运行¶

Notebook 直接调用 main 并显式传入实例路径。以下代码格相互独立，可以按需要单独运行；调整 time_limit 可以控制每个实例的求解时间。

from pathlib import Path

INSTANCE_DIR = Path.cwd() / "instances"
print("Instances:", INSTANCE_DIR)

solution_cap61 = main(INSTANCE_DIR / "cap61", time_limit=1)

solution_cap62 = main(INSTANCE_DIR / "cap62", time_limit=1)

solution_cap63 = main(INSTANCE_DIR / "cap63", time_limit=1)