Используйте линейное программирование вручную

Я занимался проблемами линейного программирования в своем классе, рисуя их графики, но я хотел бы знать, как написать программу для конкретной задачи, чтобы решить ее для меня. Если есть слишком много переменных или ограничений, я бы никогда не смог сделать это с помощью графика.

Пример задачи, максимизируйте 5x + 3y с ограничениями:

  • 5x - 2y> = 0
  • х + у <= 7
  • х <= 5
  • х> = 0
  • y> = 0

Я взял это и получил видимую область с 3 углами. х = 5 у =2 - оптимальная точка.

Как мне превратить это в код? Я знаю о симплекс-методе. И что очень важно, все ли проблемы LP будут закодированы в одной структуре? Будет ли грубая сила работать?

Источник

2 ответа

Решение

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

В дополнение к тому, что упомянуто в комментарии (Численные рецепты на C), вы также можете найти:

  1. Симплекс-Солвер от Google
  2. Тогда есть МОНЕТА ИЛИ
  3. У GNU есть свой GLPK
  4. Если вы хотите реализацию C++, то эта в Google Code действительно доступна.
  5. В R есть много реализаций, включая загрузочный пакет. (В R вы можете увидеть реализацию функции, набрав ее без скобок.)

Чтобы ответить на два других вопроса:

  1. Будут ли кодироваться все LP? Да, универсальный решатель LP может быть написан для загрузки и решения любого LP. (Существуют стандартные форматы для чтения LP, как mps а также .lp

  2. Будет ли грубая сила работать? Имейте в виду, что многие компании и крупные организации тратят много времени на тонкую настройку решателей. Есть LP, которые имеют интересные свойства, которые многие решатели будут пытаться использовать. Кроме того, некоторые вычисления могут быть решены параллельно. Алгоритм экспоненциальный, поэтому при большом количестве переменных / ограничений грубая сила не будет работать.

Надеюсь, это поможет.

Источник

Я написал это Matlab вчера, который может быть легко переписан на C++, если вы используете библиотеку Eigen или напишите свой собственный матричный класс, используя std::vector из std::vector

function [x, fval] = mySimplex(fun, A, B, lb, up)

%Examples paramters to show that the function actually works 

% sample set 1 (works for this data set)

% fun = [8 10 7];
% A = [1 3 2; 1 5 1];
% B = [10; 8];
% lb = [0; 0; 0];
% ub = [inf; inf; inf];

% sample set 2 (works for this data set)

fun = [7 8 10];
A = [2 3 2; 1 1 2];
B = [1000; 800];
lb = [0; 0; 0];
ub = [inf; inf; inf];


% generate a new slack variable for every row of A 

numSlackVars = size(A,1); % need a new slack variables for every row of A 

% Set up tableau to store algorithm data 
tableau = [A; -fun];

tableau = [tableau, eye(numSlackVars + 1)];

lastCol = [B;0];

tableau = [tableau, lastCol];

% for convienience sake, assign the following: 

numRows = size(tableau,1);
numCols = size(tableau,2);

% do simplex algorithm 

% step 0: find num of negative entries in bottom row of tableau 

numNeg = 0; % the number of negative entries in bottom row

for i=1:numCols 
    if(tableau(numRows,i) < 0)
        numNeg = numNeg + 1;
    end
end

% Remark: the number of negatives is exactly the number of iterations needed in the
% simplex algorithm 

for iterations = 1:numNeg 
    % step 1: find minimum value in last row 
    minVal = 10000; % some big number 
    minCol = 1; % start by assuming min value is the first element 
    for i=1:numCols
        if(tableau(numRows, i) < minVal)
            minVal = tableau(size(tableau,1), i);
            minCol = i; % update the index corresponding to the min element 
        end
    end 

    % step 2: Find corresponding ratio vector in pivot column 
    vectorRatio = zeros(numRows -1, 1);
    for i=1:(numRows-1) % the size of ratio vector is numCols - 1
        vectorRatio(i, 1) = tableau(i, numCols) ./ tableau(i, minCol);
    end 

    % step 3: Determine pivot element by finding minimum element in vector
    % ratio

    minVal = 10000; % some big number 
    minRatio = 1; % holds the element with the minimum ratio 

    for i=1:numRows-1
        if(vectorRatio(i,1) < minVal)
            minVal = vectorRatio(i,1);
            minRatio = i;
        end 
    end 

    % step 4: assign pivot element 

    pivotElement = tableau(minRatio, minCol);

    % step 5: perform pivot operation on tableau around the pivot element 

    tableau(minRatio, :) = tableau(minRatio, :) * (1/pivotElement);

    % step 6: perform pivot operation on rows (not including last row)

    for i=1:size(vectorRatio,1)+1 % do last row last 
        if(i ~= minRatio) % we skip over the minRatio'th element of the tableau here 
            tableau(i, :) = -tableau(i,minCol)*tableau(minRatio, :) +  tableau(i,:);
        end
    end
end 

% Now we can interpret the algo tableau 

numVars = size(A,2); % the number of cols of A is the number of variables 

x = zeros(size(size(tableau,1), 1)); % for efficiency 

% Check for basicity 
for col=1:numVars
    count_zero = 0;
    count_one = 0;
    for row = 1:size(tableau,1)
        if(tableau(row,col) < 1e-2)
            count_zero = count_zero + 1;
        elseif(tableau(row,col) - 1 < 1e-2)
            count_one = count_one + 1;
            stored_row = row; % we store this (like in memory) column for later use 
        end
    end
    if(count_zero == (size(tableau,1) -1) && count_one == 1) % this is the case where it is basic 
        x(col,1) = tableau(stored_row, numCols);
    else 
        x(col,1) = 0; % this is the base where it is not basic 
    end
end

% find function optimal value at optimal solution 
fval = x(1,1) * fun(1,1); % just needed for logic to work here 
for i=2:numVars 
    fval = fval + x(i,1) * fun(1,i);
end


end
Источник
Другие вопросы по тегам