Фабрика, склад, модель выделения продаж: основанный на решателе

В этом примере показано, как настроить и решить смешано-целочисленную задачу линейного программирования. Проблема состоит в том, чтобы найти оптимальные уровни производства и распределения среди набора фабрик, складов и торговых точек. Для подхода, основанного на проблеме смотрите Фабрику, Склад, Модель Выделения Продаж: основанный на проблеме.

Пример сначала генерирует случайные местоположения для фабрик, склады и торговые точки. Не стесняйтесь изменять масштабный коэффициент N, который масштабирует и размер сетки, в которой находятся производственные объекты и распределительные сети, но также и масштабирует количество этих средств так, чтобы плотность средств каждого типа на область сетки была независима от N.

Местоположения средства

Для данного значения масштабного коэффициента N, предположите, что там следующие:

  • fN2 фабрики

  • wN2 склады

  • sN2 торговые точки

Эти средства находятся на отдельных целочисленных узлах решетки между 1 и N в x и y направления. Чтобы средства имели отдельные местоположения, вы требуете этого f+w+s1. В этом примере взять N=20, f=0.05, w=0.05, и s=0.1.

Производство и распределение

Существуют P продукты сделаны фабриками. Взять P=20.

Спрос на каждый продукт p в торговой точке s d(s,p). Спрос является количеством, которое может быть продано во временном интервале. Одно ограничение на модель состоит в том, что требование соблюдается, означая, что система производит и распределяет точно количества в спросе.

Существуют полные ограничения на каждую фабрику и каждый склад.

  • Производство продукта p на фабрике f меньше pcap(f,p).

  • Способность склада w wcap(w).

  • Сумма продукта p это может быть транспортировано со склада w к торговой точке во временном интервале меньше turn(p)*wcap(w), где turn(p) текучесть кадров продукта p.

Предположим, что каждая торговая точка получает свои предоставления всего от одного склада. Часть проблемы должна определить самое дешевое отображение торговых точек в склады.

Затраты

Стоимость переноса продуктов от фабрики до склада, и от склада до торговой точки, зависит от расстояния между средствами, и на конкретном продукте. Если dist(a,b) расстояние между средствами a и b, затем стоимость поставки продукта p между этими средствами времена расстояния стоимость транспортировки tcost(p):

dist(a,b)*tcost(p).

Расстояние в этом примере является расстоянием сетки, также известным L1 расстояние. Это - сумма абсолютной разности в x координаты и y координаты.

Стоимость создания модуля продукта p на фабрике f pcost(f,p).

Задача оптимизации

Учитывая набор местоположений средства, и требования и полные ограничения, найдите:

  • Производственный уровень каждого продукта на каждой фабрике

  • Расписание распределения для продуктов от фабрик до складов

  • Расписание распределения для продуктов от складов до торговых точек

Эти количества должны гарантировать, что спросу удовлетворяют, и общая стоимость минимизирована. Кроме того, каждая торговая точка требуется, чтобы получать все свои продукты точно от одного склада.

Переменные и уравнения для задачи оптимизации

Контрольные переменные, означая тех, которых можно изменить в оптимизации,

  • x(p,f,w) = сумма продукта p это транспортируется от фабрики f в склад w

  • y(s,w) = бинарное переменное принимающее значение 1, когда торговая точка s сопоставлен со складом w

Целевая функция, чтобы минимизировать

fpwx(p,f,w)(pcost(f,p)+tcost(p)dist(f,w))

+swp(d(s,p)tcost(p)dist(s,w)y(s,w)).

Ограничения

wx(p,f,w)pcap(f,p) (способность фабрики).

fx(p,f,w)=s(d(s,p)y(s,w)) (требование соблюдается).

psd(s,p)turn(p)y(s,w)wcap(w) (способность склада).

wy(s,w)=1 (каждая торговая точка сопоставляет в один склад).

x(p,f,w)0 (неотрицательное производство).

y(s,w)ϵ{0,1} двоичный файл y).

Переменные x и y появитесь в цели, и ограничение функционирует линейно. Поскольку y ограничивается целочисленными значениями, проблемой является смешано-целочисленная линейная программа (MILP).

Сгенерируйте случайную проблему: местоположения средства

Установите значения N, f, w, и s параметры, и генерируют местоположения средства.

rng(1) % for reproducibility
N = 20; % N from 10 to 30 seems to work. Choose large values with caution.
N2 = N*N;
f = 0.05; % density of factories
w = 0.05; % density of warehouses
s = 0.1; % density of sales outlets

F = floor(f*N2); % number of factories
W = floor(w*N2); % number of warehouses
S = floor(s*N2); % number of sales outlets

xyloc = randperm(N2,F+W+S); % unique locations of facilities
[xloc,yloc] = ind2sub([N N],xyloc);

Конечно, не реалистично взять случайные местоположения для средств. Этот пример предназначается, чтобы показать методы решения, не то, как сгенерировать хорошие местоположения средства.

Постройте средства. Средства 1 через F являются фабриками, F+1 через F+W являются складами, и F+W+1 через F+W+S являются торговыми точками.

h = figure;
plot(xloc(1:F),yloc(1:F),'rs',xloc(F+1:F+W),yloc(F+1:F+W),'k*',...
    xloc(F+W+1:F+W+S),yloc(F+W+1:F+W+S),'bo');
lgnd = legend('Factory','Warehouse','Sales outlet','Location','EastOutside');
lgnd.AutoUpdate = 'off';
xlim([0 N+1]);ylim([0 N+1])

Сгенерируйте случайные мощности, затраты и требования

Сгенерируйте случайные производственные затраты, мощности, текучесть кадров и требования.

P = 20; % 20 products

% Production costs between 20 and 100
pcost = 80*rand(F,P) + 20;

% Production capacity between 500 and 1500 for each product/factory
pcap = 1000*rand(F,P) + 500;

% Warehouse capacity between P*400 and P*800 for each product/warehouse
wcap = P*400*rand(W,1) + P*400;

% Product turnover rate between 1 and 3 for each product
turn = 2*rand(1,P) + 1;

% Product transport cost per distance between 5 and 10 for each product
tcost = 5*rand(1,P) + 5;

% Product demand by sales outlet between 200 and 500 for each
% product/outlet
d = 300*rand(S,P) + 200;

Эти случайные требования и мощности могут привести к неосуществимым проблемам. Другими словами, иногда спрос превышает производство и складские полные ограничения. Если вы измените некоторые параметры и получите неосуществимую проблему, во время решения вы получите exitflag-2.

Сгенерируйте матрицы цели и ограничения и векторы

Вектор целевой функции obj в intlincon состоит из коэффициентов переменных x(p,f,w) и y(s,w). Таким образом, существует естественно P*F*W + S*W коэффициенты в obj.

Один способ сгенерировать коэффициенты состоит в том, чтобы начаться с P-by-F-by-W массив obj1 для x коэффициенты и S-by-W массив obj2 для y(s,w) коэффициенты. Затем преобразуйте эти массивы в два вектора и объедините их в obj путем вызова

obj = [obj1(:);obj2(:)];

obj1 = zeros(P,F,W); % Allocate arrays
obj2 = zeros(S,W);

В течение генерации цели и ограничительных векторов и матриц, мы генерируем (p,f,w) массив или (s,w) массив, и затем преобразует результат в вектор.

Чтобы начать генерировать входные параметры, сгенерируйте массивы расстояния distfw(i,j) и distsw(i,j).

distfw = zeros(F,W); % Allocate matrix for factory-warehouse distances
for ii = 1:F
    for jj = 1:W
        distfw(ii,jj) = abs(xloc(ii) - xloc(F + jj)) + abs(yloc(ii) ...
            - yloc(F + jj));
    end
end

distsw = zeros(S,W); % Allocate matrix for sales outlet-warehouse distances
for ii = 1:S
    for jj = 1:W
        distsw(ii,jj) = abs(xloc(F + W + ii) - xloc(F + jj)) ...
            + abs(yloc(F + W + ii) - yloc(F + jj));
    end
end

Сгенерируйте записи obj1 и obj2.

for ii = 1:P
    for jj = 1:F
        for kk = 1:W
            obj1(ii,jj,kk) = pcost(jj,ii) + tcost(ii)*distfw(jj,kk);
        end
    end
end

for ii = 1:S
    for jj = 1:W
        obj2(ii,jj) = distsw(ii,jj)*sum(d(ii,:).*tcost);
    end
end

Объедините записи в один вектор.

obj = [obj1(:);obj2(:)]; % obj is the objective function vector

Теперь создайте ограничительные матрицы.

Ширина каждой линейной матрицы ограничений является длиной obj вектор.

matwid = length(obj);

Существует два типа линейных неравенств: ограничения производственной мощности и складские полные ограничения.

Существует P*F ограничения производственной мощности и W складские полные ограничения. Ограничительные матрицы довольно разреженны на порядке ненулевого 1%, поэтому сохраните память при помощи разреженных матриц.

Aineq = spalloc(P*F + W,matwid,P*F*W + S*W); % Allocate sparse Aeq
bineq = zeros(P*F + W,1); % Allocate bineq as full

% Zero matrices of convenient sizes:
clearer1 = zeros(size(obj1));
clearer12 = clearer1(:);
clearer2 = zeros(size(obj2));
clearer22 = clearer2(:);

% First the production capacity constraints
counter = 1;
for ii = 1:F
    for jj = 1:P
        xtemp = clearer1;
        xtemp(jj,ii,:) = 1; % Sum over warehouses for each product and factory
        xtemp = sparse([xtemp(:);clearer22]); % Convert to sparse
        Aineq(counter,:) = xtemp'; % Fill in the row
        bineq(counter) = pcap(ii,jj);
        counter = counter + 1;
    end
end

% Now the warehouse capacity constraints
vj = zeros(S,1); % The multipliers 
for jj = 1:S
    vj(jj) = sum(d(jj,:)./turn); % A sum of P elements
end

for ii = 1:W
    xtemp = clearer2;
    xtemp(:,ii) = vj;
    xtemp = sparse([clearer12;xtemp(:)]); % Convert to sparse
    Aineq(counter,:) = xtemp'; % Fill in the row
    bineq(counter) = wcap(ii);
    counter = counter + 1;
end

Существует два типа линейных ограничений равенства: ограничение, что требование соблюдается, и ограничение, что каждая торговая точка соответствует одному складу.

Aeq = spalloc(P*W + S,matwid,P*W*(F+S) + S*W); % Allocate as sparse
beq = zeros(P*W + S,1); % Allocate vectors as full

counter = 1;
% Demand is satisfied:
for ii = 1:P
    for jj = 1:W
        xtemp = clearer1;
        xtemp(ii,:,jj) = 1;
        xtemp2 = clearer2;
        xtemp2(:,jj) = -d(:,ii);
        xtemp = sparse([xtemp(:);xtemp2(:)]'); % Change to sparse row
        Aeq(counter,:) = xtemp; % Fill in row
        counter = counter + 1;
    end
end

% Only one warehouse for each sales outlet:
for ii = 1:S
    xtemp = clearer2;
    xtemp(ii,:) = 1;
    xtemp = sparse([clearer12;xtemp(:)]'); % Change to sparse row
    Aeq(counter,:) = xtemp; % Fill in row
    beq(counter) = 1;
    counter = counter + 1;
end

Связанные ограничения и целочисленные переменные

Целочисленные переменные - те от length(obj1) + 1 в конец.

intcon = P*F*W+1:length(obj);

Верхние границы от length(obj1) + 1 в конец также.

lb = zeros(length(obj),1);
ub = Inf(length(obj),1);
ub(P*F*W+1:end) = 1;

Выключите итеративное отображение так, чтобы вы не получали сотни линий выхода. Включайте функцию построения графика, чтобы контролировать прогресс решения.

opts = optimoptions('intlinprog','Display','off','PlotFcn',@optimplotmilp);

Решите задачу

Вы сгенерировали все входные параметры решателя. Вызовите решатель, чтобы найти решение.

[solution,fval,exitflag,output] = intlinprog(obj,intcon,...
                                     Aineq,bineq,Aeq,beq,lb,ub,opts);

if isempty(solution) % If the problem is infeasible or you stopped early with no solution
    disp('intlinprog did not return a solution.')
    return % Stop the script because there is nothing to examine
end

Исследуйте решение

Решение выполнимо, к в данных допусках.

exitflag
exitflag = 1
infeas1 = max(Aineq*solution - bineq)
infeas1 = 9.0949e-13
infeas2 = norm(Aeq*solution - beq,Inf)
infeas2 = 1.0232e-12

Проверяйте, что целочисленные компоненты являются действительно целыми числами или достаточно близки, что это разумно к раунду их. Чтобы изучить, почему эти переменные не могут быть точно целыми числами, см. документацию.

diffint = norm(solution(intcon) - round(solution(intcon)),Inf)
diffint = 1.1857e-13

Некоторые целочисленные переменные не являются точно целыми числами, но все очень близки. Так вокруг целочисленных переменных.

solution(intcon) = round(solution(intcon));

Проверяйте выполнимость округленного решения и изменение в значении целевой функции.

infeas1 = max(Aineq*solution - bineq)
infeas1 = 9.0949e-13
infeas2 = norm(Aeq*solution - beq,Inf)
infeas2 = 3.0667e-11
diffrounding = norm(fval - obj(:)'*solution,Inf)
diffrounding = 1.1176e-08

Округление решения не заметно изменило свою выполнимость.

Можно исследовать решение наиболее легко путем изменения его назад к его исходным размерностям.

solution1 = solution(1:P*F*W); % The continuous variables
solution2 = solution(intcon); % The integer variables
solution1 = reshape(solution1,P,F,W);
solution2 = reshape(solution2,S,W);

Например, сколько торговых точек сопоставлено с каждым складом? Заметьте, что в этом случае некоторые склады имеют 0 связанных выходов, означая, что склады не используются в оптимальном решении.

outlets = sum(solution2,1) % Sum over the sales outlets
outlets = 1×20

     3     0     3     2     2     2     3     2     3     1     1     0     0     3     4     3     2     3     2     1

Постройте связь между каждой торговой точкой и ее складом.

figure(h);
hold on
for ii = 1:S
    jj = find(solution2(ii,:)); % Index of warehouse associated with ii
    xsales = xloc(F+W+ii); ysales = yloc(F+W+ii);
    xwarehouse = xloc(F+jj); ywarehouse = yloc(F+jj);
    if rand(1) < .5 % Draw y direction first half the time
        plot([xsales,xsales,xwarehouse],[ysales,ywarehouse,ywarehouse],'g--')
    else % Draw x direction first the rest of the time
        plot([xsales,xwarehouse,xwarehouse],[ysales,ysales,ywarehouse],'g--')
    end
end
hold off

title('Mapping of sales outlets to warehouses')

Черный цвет * без зеленых линий представляет неиспользованные склады.

Похожие темы

Для просмотра документации необходимо авторизоваться на сайте