Как реализовать CVRPTW в IBM ILOG CPLEX с различными данными и размерами

Question

Не удается реализовать CVRPTW (Задача маршрутизации транспортных средств с ограниченной мощностью и временными окнами) в IBM ILOG CPLEX с адекватным решением.

Мой код .mod:

int n = ...; // Множество вершин.
int nVehicles = ...; // Количество транспортных средств.
range N = 0..n-1; // Диапазон всех точек (включая склад).
range C = 1..n-1; // Диапазон всех клиентов (исключая склад).
range V = 0..nVehicles-1; // Диапазон всех транспортных средств.

int c[N][N] = ...;    // Транспортные расходы.
dvar boolean x[N][N][V]; // 1, если k-е транспортное средство движется от i к j, в противном случае 0.
dvar int+ Z[V];
int Q = ...; // Вместимость транспортного средства.
int d[C] = ...; // Спрос клиентов.
dvar int+ s[N][V]; // Время прибытия.
int t[N][N] = ...;    // Время транспортировки.
int a[N] = ...;       // Начало временных окон.
int b[N] = ...;       // Конец временных окон.

// 1. Целевая функция: минимизировать общие транспортные расходы.
minimize sum(k in V) Z[k];

subject to {
    // Все транспортные средства выезжают со склада в 0 времени.
    //forall(k in V)
        //s[0][k] == 0; 
        
    // 2. Каждый клиент посещается ровно один раз.
    forall(i in C)
        sum(k in V, j in N) x[i][j][k] == 1;

    // 3. Ограничение по мощностям.
    forall(k in V)
        sum(i in C, j in N) (d[i] * x[i][j][k]) <= Q;
        
    // 4. Каждое транспортное средство начинает свой маршрут со склада.
    forall(k in V)
        sum(j in N) x[0][j][k] == 1;
        
    // 5. Передвижение между клиентами.
    forall(k in V, h in C)
        sum(i in N) x[i][h][k] == sum(j in N) x[h][j][k];
        
    // 6. Каждое транспортное средство возвращается на склад.
    forall(k in V)
        sum(i in N) x[i][0][k] == 1; // Возможно 0 вместо n-1.
        
    // 7. Время выезда из одной точки и время приезда в другую точку.
    forall(i in N, j in N, k in V)
        s[i][k] + t[i][j] <= s[j][k] + (1 - x[i][j][k]) * 10000; // Большое число для отключения ограничения, когда x[i][j][k] = 0.
        
    // 8. Ограничения временных окон.
    forall(i in N, k in V)
        a[i] <= s[i][k] <= b[i];
    
    // 10. Ограничение на использование указанного количества транспортных средств.
    sum(k in V, j in N) x[0][j][k] <= nVehicles;
    
    // Ограничение, связывающее переменную Z[k] с общими транспортными расходами
    forall(k in V)
        Z[k] == sum(i in N, j in N) (c[i][j] * x[i][j][k]);
}

Код .dat:

n = 6; // Множество вершин.
nVehicles = 4; // Количество транспортных средств.

// Транспортные расходы.
c = [
    [0, 10, 15, 20, 30, 25],
    [10, 0, 9, 14, 23, 18],
    [15, 9, 0, 10, 12, 25],
    [20, 14, 10, 0, 11, 15],
    [30, 23, 12, 11, 0, 12],
    [25, 18, 25, 15, 12, 0]
];

// Время транспортировки.
t = [
    [0, 2, 3, 4, 5, 6],
    [2, 0, 1, 3, 5, 7],
    [3, 1, 0, 2, 4, 6],
    [4, 3, 2, 0, 3, 5],
    [5, 5, 4, 3, 0, 2],
    [6, 7, 6, 5, 2, 0]
];

Q = 100; // Вместимость транспортного средства.
d = [2, 3, 4, 2, 1]; // Спрос.

// Временные окна
a = [0, 3, 5, 7, 9, 12];
b = [100, 10, 12, 14, 16, 18];

В результате я получаю нули в Z. Когда я меняю матрицу c диагонально на большие числа (например, 100000), я получаю 600000 100000 100000 100000.
Результаты, которые я хотел бы получить, должны выглядеть так:
Пример вывода

Предположим, у нас есть 4 транспортных средства и 5 клиентов. Результаты могут выглядеть так:

Маршрут для Транспортного Средства 1: Склад → Клиент 3 → Клиент 1 → Склад
Маршрут для Транспортного Средства 2: Склад → Клиент 2 → Клиент 4 → Клиент 5 → Склад
Маршрут для Транспортного Средства 3: Склад (не использовано)
Маршрут для Транспортного Средства 4: Склад (не использовано)
Общие расходы: 150 (например, общая дистанция или время маршрута)
Времена прибытия:
- Транспортное средство 1: Клиент 3 — 5, Клиент 1 — 9
- Транспортное средство 2: Клиент 2 — 4, Клиент 4 — 8, Клиент 5 — 12

Как реализовать CVRPTW в IBM ILOG CPLEX с различными данными и размерами

Вопрос или проблема

Ответ или решение

Ошибки и недостатки в модели

Корректировки в Коде

Оптимизация решения