建模实战|第七期：电信公司人员调度优化问题-使用ortools求解器（附python代码)

1.问题背景

技术人员路由和排程(Technician Routing and Scheduling-TRS)是电信企业面临的一个普遍问题，电信企业必须能够为大量的客户提供服务。由于资源有限，这些公司必须尽力提供及时、负担得起和可靠的服务，以最大限度地提高客户满意度。TRS问题涉及多技能技术人员的分配、调度和路由，以独特的优先级、服务时间窗口、处理时间、技能要求和地理位置为客户服务。在目前的实践中，这些决定往往是由操作员做出的，有时是手动的——这种方法不仅耗时，而且可能大大偏离最优值。

这个示例我们将用数学优化的方法帮助电信公司解决TRS问题，使电信运营商自动化和改进其技术人员分配、调度和路由决策，以确保最高水平的客户满意度。

2.问题描述

电信公司经营多个服务中心为其客户服务。每个服务中心都有自己的技术人员，这些技术人员从服务中心派遣出去去完成指定的工作，待所有工作完成后再返回服务中心。技术人员有多种技能和可用的工作能力，在调度范围内不能超过他们的这些技能和可用的工作能力。服务订单/作业具有已知的处理时间、客户指定的启动服务的时间窗口、完成服务的截止日期及其技能要求。根据工作的性质(日常维护或紧急情况)以及与客户的关系(现有的或新的)，公司评估工作的重要性并为其分配优先级分数。一项工作最多分配给一名拥有所需技能的技术人员，但该技术人员在调度范围内的可用能力不能被超过。

为了解决TRS问题，电信公司必须能够同时做出三种类型的决策:

• (1)将工作分配给所有服务中心的技术人员;
• (2)每个技术人员的路线，即技术人员拜访客户的顺序;
• (3)工作调度，即技术人员到达客户地点并完成相应工作的时间安排。

目标是最小化所有工作的总加权延迟，其优先级为权重。

必须满足以下约束条件:

• (1)技术人员在使用的情况下，离开他/她所在的服务中心，并在他/她被分配的工作完成后返回同一个服务中心。
• (2)在调度范围内，不能超过技术人员的可用容量。
• (3)一项工作，如果被选中，最多分配给一名拥有所需技能的技术人员。
• (4)技术人员必须在客户指定的时间间隔(时间窗口)内到达客户地点，并在客户要求的最后期限前完成工作。这是保证顾客满意的重要约束。

从上述描述可以看出，基本的TRS问题是带时间窗的多车场车辆路径问题的一种变体，在文献[1]中称为MDVRPTW。

3.人员调度与排程问题建模

 目标函数说明：我们希望将满足作业需求的约束和在时间窗口内启动作业的约束视为软约束。软约束可以被违反，但违反将招致巨大的惩罚。考虑到每个作业具有不同的优先级，我们使用优先级惩罚来计算与违反软约束相关的惩罚。我们假设不满足工作需求会极大地降低客户满意度，因此应该招致与软约束相关的最大惩罚。回想一下，参数M的值是一个很大的数字，那么不满足工作需求的相关惩罚系数确定为:  。违反时间窗口约束也会使客户满意度降低，但其程度要小一些，所以违反时间窗约束的相关惩罚系数确定为: 

 

4.人员调度与排程问题的一个实例

在这个场景中，我们考虑下一个工作日的技术人员的路由和调度问题，以使客户安排的延迟最小化。这家电信公司有7名技术人员:Albert, Bob, Carlos, Doris, Ed, Flor, and Gina.。有两个服务中心:Heidelberg and Freiburg im Breisgau。技术人员只在其中一个服务中心工作。每个技术人员的工作小时数和服务中心如下表所示：

	Albert	Bob	Carlos	Doris	Ed	Flor	Gina
Minutes（分钟）	480	480	480	480	480	360	360
Depot（服务中心）	Heidelberg	Heidelberg	Freiburg im Breisgau	Freiburg im Breisgau	Heidelberg	Freiburg im Breisgau	Heidelberg

电信公司有不同类型的工作。下表显示了作业类型的优先级(最高优先级为4，最不重要的为1)和持续时间(以小时为单位)：

	Priority	Duration (min)
Equipment Installation（设备安装）	2	60
Equipment Setup（设备设置）	3	30
Inspect/Service Equipment（设备检查）	1	60
Repair - Regular（设备修复-常规）	1	60
Repair - Important（设备修复-重要）	2	120
Repair - Urgent（设备修复-紧急）	3	90
Repair - Critical（设备修复-致命）	4	60

各技术人员所能胜任的作业类型如下表所示：

	Albert	Bob	Carlos	Doris	Ed	Flor	Gina
Equipment Installation（设备安装）				1		1	1
Equipment Setup（设备设置）	1	1	1	1			1
Inspect/Service Equipment（设备修复）	1		1		1
Repair - Regular（设备修复-常规）		1	1		1	1	1
Repair - Important（设备修复-重要）				1		1	1
Repair - Urgent（设备修复-紧急）		1	1		1	1	1
Repair - Critical（设备修复-致命）				1			1

电信公司从客户那里接收不同的作业类型。这些作业有不同的作业类型、预约(到期)时间和服务时间窗口（技术人员可以满足的）。客户的订单及要求如下表所示。对于每个客户，有不同的位置。

	C1:Mannheim	C2: Karlsruhe	C3: Baden-Baden	C4: Bühl	C5: Offenburg	C6: Lahr/Schwarzwald	C7: Lörrach
Job type（作业类型）	Equipment Setup	Equipment Setup	Repair - Regular	Equipment Installation	Equipment Installation	Repair - Critical	Inspect/Service Equipment
Due time（到期时间）	8:00	10:00	11:00	12:00	14:00	15:00	16:00
Time Window（服务时间窗口）	7:00-7:30	7:30-9:30	8:00-10:00	9:00-11:00	11:00-13:00	12:00-14:00	13:00-15:00

作业的时间范围从7:00到17:00，即10个小时。时间段以分钟为单位，那么到期时间和时间窗口将转换为分钟，从0分钟开始到600分钟结束。例如，客户C2的到期时间为10:00(180分钟)，时间窗口为7:30到9:30(30分钟到150分钟)。

下表显示了从任何服务中心或客户地点到任何服务中心或客户地点的行驶时间(以分钟为单位)：

5.python调用Ortools求解器实现问题的求解

class Technician():
    def __init__(self, name, cap, depot):
        self.name = name
        self.cap = cap
        self.depot = depot

    def __str__(self):
        return f"Technician: {self.name}\n  Capacity: {self.cap}\n  Depot: {self.depot}"


class Job():
    def __init__(self, name, priority, duration, coveredBy):
        self.name = name
        self.priority = priority
        self.duration = duration
        self.coveredBy = coveredBy

    def __str__(self):
        about = f"Job: {self.name}\n  Priority: {self.priority}\n  Duration: {self.duration}\n  Covered by: "
        about += ", ".join([t.name for t in self.coveredBy])
        return about


class Customer():
    def __init__(self, name, loc, job, tStart, tEnd, tDue):
        self.name = name
        self.loc = loc
        self.job = job
        self.tStart = tStart
        self.tEnd = tEnd
        self.tDue = tDue

    def __str__(self):
        coveredBy = ", ".join([t.name for t in self.job.coveredBy])
        return f"Customer: {self.name}\n  Location: {self.loc}\n  Job: {self.job.name}\n  Priority: {self.job.priority}\n  Duration: {self.job.duration}\n  Covered by: {coveredBy}\n  Start time: {self.tStart}\n  End time: {self.tEnd}\n  Due time: {self.tDue}"


import sys
import pandas as pd
from ortools.linear_solver import pywraplp


class TSR():
    # 初始化函数
    def __init__(self):
        self.technicians = []  # 技术人员个数
        self.customers = []  # 客户地点个数
        self.dist = {}  # 各个位置之间的距离矩阵

    def read_data(self):
        # Read Excel workbook
        excel_file = 'data-Sce0.xlsx'
        df = pd.read_excel(excel_file, sheet_name='Technicians')
        df = df.rename(columns={df.columns[0]: "name", df.columns[1]: "cap", df.columns[2]: "depot"})

        df1 = df.drop(df.columns[3:], axis=1).drop(df.index[[0, 1]])
        # Create Technician objects
        self.technicians = [Technician(*row) for row in df1.itertuples(index=False, name=None)]
        # Read job data
        jobs = []
        for j in range(3, len(df.columns)):
            coveredBy = [t for i, t in enumerate(self.technicians) if df.iloc[2 + i, j] == 1]
            thisJob = Job(df.iloc[2:, j].name, df.iloc[0, j], df.iloc[1, j], coveredBy)
            jobs.append(thisJob)

        # Read location data
        df_locations = pd.read_excel(excel_file, sheet_name='Locations', index_col=0)  # , skiprows=1, index_col=0)

        # Extract locations and initialize distance dictionary
        locations = df_locations.index
        self.dist = {(l, l): 0 for l in locations}
        # Populate distance dictionary
        for i, l1 in enumerate(locations):
            for j, l2 in enumerate(locations):
                if i < j:
                    self.dist[l1, l2] = df_locations.iloc[i, j]
                    self.dist[l2, l1] = self.dist[l1, l2]

        # Read customer data
        df_customers = pd.read_excel(excel_file, sheet_name='Customers')
        for i, c in enumerate(df_customers.iloc[:, 0]):
            job_name = df_customers.iloc[i, 2]

            # Find the corresponding Job object
            matching_job = next((job for job in jobs if job.name == job_name), None)

            if matching_job is not None:
                # Create Customer object using corresponding Job object
                this_customer = Customer(c, df_customers.iloc[i, 1], matching_job, *df_customers.iloc[i, 3:])
                self.customers.append(this_customer)

    def solve_trs0(self):
        # Build useful data structures
        K = [k.name for k in self.technicians]
        J = [j.loc for j in self.customers]
        L = list(set([l[0] for l in self.dist.keys()]))
        D = list(set([d.depot for d in self.technicians]))
        W = {k.name: k.cap for k in self.technicians}
        # loc = {j.name: j.loc for j in self.customers}
        depot = {k.name: k.depot for k in self.technicians}
        canCover = {j.loc: [k.name for k in j.job.coveredBy] for j in self.customers}  # 1.这个客户地点内可达的技术人员的name集合
        p = {j.loc: j.job.duration for j in self.customers}
        tStart = {j.loc: j.tStart for j in self.customers}
        tEnd = {j.loc: j.tEnd for j in self.customers}
        tDue = {j.loc: j.tDue for j in self.customers}
        priority = {j.loc: j.job.priority for j in self.customers}

        ### Create solver
        solver = pywraplp.Solver.CreateSolver('SCIP')

        ### Decision variables
        # Customer-technician assignment
        x = {}
        u = {}
        y = {}
        for j in J:  # 2.range(len(J)) 改成了 J
            for k in K:
                x[(j, k)] = solver.BoolVar(f'x[{j},{k}]')  # 3.solver.NewBoolVar(f'x{j}{k}')  # 变量x[j, k]改成了这种格式的x[(j, k)]

        # Technician assignment
        for k in K:
            u[k] = solver.BoolVar(f'u[{k}]')  #

        # Edge-route assignment to technician
        for i in L:
            for j in L:
                for k in K:
                    y[(i, j, k)] = solver.BoolVar(f'y[{i},{j},{k}]')  #

        # Start time of service
        s = {}
        for i in L:
            s[i] = solver.IntVar(lb=float('-inf'), ub=600, name=f's[{i}]')  # 4.声明整数变量时，需要设定lb,ub
        # Lateness of service
        z = {}
        for j in J:
            z[j] = solver.IntVar(lb=0, ub=float('inf'), name=f'z[{j}]')

        # Artificial variables to correct time window upper and lower limits
        xa = {}
        xb = {}
        for j in J:
            xa[j] = solver.IntVar(lb=0, ub=float('inf'), name=f'xa[{j}]')
            xb[j] = solver.IntVar(lb=0, ub=float('inf'), name=f'xb[{j}]')

        # Unfilled jobs
        g = {}
        for j in J:
            g[j] = solver.BoolVar(f'g_{j}')  #

        # Technician cannot leave or return to a depot that is not its base
        for k in self.technicians:
            for d in D:
                if k.depot != d:
                    for i in L:
                        solver.Add(y[i, d, k.name] == 0)
                        solver.Add(y[d, i, k.name] == 0)

        ### Constraints

        # (1) A technician must be assigned to a job, or a gap is declared
        for j in J:
            solver.Add(sum(x[(j, k)] for k in canCover[j]) + g[j] == 1, name = 'assignToJob')

        # (2) At most one technician can be assigned to a job
        for j in J:
            solver.Add(sum(x[(j, k)] for k in canCover[j]) <= 1, name = 'assignOne')

        # (3) Technician capacity constraints
        for k in K:
            capLHS = sum(p[j] * x[(j, k)] for j in J) + sum(self.dist[(i, j)] * y[(i, j, k)] for i in L for j in L)
            solver.Add(capLHS <= W[k] * u[k], name = 'techCapacity')

        #####部分约束省略######

        ### Objective function
        M = 6100
        active = sum(priority[j] * z[j] for j in J)
        per_mw = sum(0.01 * M * priority[j] * (xa[j] + xb[j]) for j in J)
        startup_obj = sum(M * priority[j] * g[j] for j in J)
        solver.Minimize(active + per_mw + startup_obj)
        # Solve
        # Solve problem
        status = solver.Solve()

        if status != pywraplp.Solver.OPTIMAL:
            print("Model is either infeasible or unbounded.")
            sys.exit(0)

        ### Print results
        # Assignments
        print("")
        for c in self.customers:
            if g[c.loc].solution_value() > 0.5:
                jobStr = "Nobody assigned to {} ({}) in {}".format(c.name, c.job.name, c.loc)
            else:
                for k in K:
                    if x[(c.loc, k)].solution_value() > 0.5:
                        jobStr = "{} assigned to {} ({}) in {}. Start at s={:.2f}.".format(k, c.name, c.job.name, c.loc,
                                                                                           s[c.loc].solution_value())
                        if z[c.loc].solution_value() > 1e-6:
                            jobStr += " {:.2f} minutes late.".format(z[c.loc].solution_value())
                        if xa[c.loc].solution_value() > 1e-6:
                            jobStr += " Start time corrected by {:.2f} minutes.".format(xa[c.loc].solution_value())
                        if xb[c.loc].solution_value() > 1e-6:
                            jobStr += " End time corrected by {:.2f} minutes.".format(xb[c.loc].solution_value())
            print(jobStr)

        # Technicians
        print("")
        for k in self.technicians:
            if u[k.name].solution_value() > 0.5:
                cur = k.depot
                route = k.depot
                while True:
                    for c in self.customers:
                        if y[(cur, c.loc, k.name)].solution_value() > 0.5:
                            route += " -> {} (dist={}, s={:.2f}, proc={})".format(c.loc, self.dist[(cur, c.loc)],
                                                                                  s[c.loc].solution_value(),
                                                                                  c.job.duration)
                            cur = c.loc
                    for i in D:  # D是服务中心的集合
                        if y[(cur, i, k.name)].solution_value() > 0.5:
                            route += " -> {} (dist={})".format(i, self.dist[(cur, i)])
                            cur = i
                            break
                    if cur == k.depot:
                        break
                print("{}'s route: {}".format(k.name, route))
            else:
                print("{} is not used".format(k.name))

    def printScen(self, scenStr):
        sLen = len(scenStr)
        print("\n" + "*" * sLen + "\n" + scenStr + "\n" + "*" * sLen + "\n")


if __name__ == "__main__":
    # Base solver
    tsr = TSR()
    tsr.printScen("Solving base scenario solver")
    tsr.read_data()
    tsr.solve_trs0()

6.参考文献

[1] S. Salhi, A. Imran, N. A. Wassan.The multi-depot vehicle routing problem with heterogeneous vehicle fleet: Formulation and a variable neighborhood search implementation. Computers & Operations Research 52 (2014) 315-325.

7.相关阅读

Gurobi-modeling-examples

 完整代码：在这里