2013-06-25 1 views
8

У меня есть 12 наборов векторов (около 10-20 векторов каждый), и я хочу выбрать один вектор каждого набора, чтобы максимизировать функцию f, которая принимает сумму этих векторов в качестве аргумента. Кроме того, у меня есть ограничения для некоторых компонентов этой суммы.Оптимизация с дискретными параметрами в Matlab

Пример:

a_1 = [3 2 0 5], a_2 = [3 0 0 2], a_3 = [6 0 1 1], ... , a_20 = [2 12 4 3] 
b_1 = [4 0 4 -2], b_2 = [0 0 1 0], b_3 = [2 0 0 4], ... , b_16 = [0 9 2 3] 
... 
l_1 = [4 0 2 0], l_2 = [0 1 -2 0], l_3 = [4 4 0 1], ... , l_19 = [3 0 9 0] 

s = [s_1 s_2 s_3 s_4] = a_x + b_y + ... + l_z 

Ограничения:

s_1 > 40 
s_2 < 100 
s_4 > -20 

Цель: Выбрали х, у, ..., г, чтобы максимизировать F (S):

f(s) -> max 

где Р является нелинейной функцией, которая принимает вектор s и возвращает скаляр.

Bruteforcing занимает слишком много времени, потому что существует около 5,9 трлн комбинаций, и поскольку мне нужно максимальное (или даже лучшее из 10 лучших комбинаций), я не могу использовать какие-либо из жадных алгоритмов, которые пришли мне на ум.

Векторы довольно разреженные, около 70-90% - это нули. Если это как-то помогает ...?

Инструментарий оптимизации Matlab не помог ни потому, что он не поддерживает отдельную оптимизацию.

+1

Вы можете что-то рассказать о нелинейной функции 'f (s)'? что ты знаешь об этом? что вы можете предположить? – Shai

+0

Если вы можете предоставить полностью воспроизводимый пример (с ограничениями и подробной функцией obj), мы можем предложить способы решения проблемы. –

+1

С ограничениями пространство поиска не так велико. – Bytemain

ответ

5

В основном это проблема блокировки, когда контакты замка имеют 20 различных положений, а также 12 контактов. Также:

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

... интересно!

на основе подхода RASMAN и комментарии Phpdna в, и предположении, что вы используете int8 как тип данных, при заданных ограничениях есть

>> d = double(intmax('int8')); 
>> (d-40) * (d+100) * (d+20) * 2*d 
ans = 
    737388162 

возможные векторы s (дать или взять несколько, убежища думал о + 1 и т. д.). ~ 740   миллионов оценок вашего относительно простого f(s) не должны занимать более 2 секунд, и, найдя всеs, которые максимизируют f(s), у вас остается проблема поиска линейных комбинаций в вашем векторном наборе, которые составляют одну из эти решения s.

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

int16: ans = 2.311325368800510e+018 
int32: ans = 4.253529737045237e+037 
int64: ans = 1.447401115466452e+076 

Итак, я рассмотрю более прямой и более общий подход здесь.

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

К сожалению, в инструменте оптимизации (или в Файловом обмене, насколько я мог найти) ничего подобного нет. fmincon - это не-go, так как он использует информацию о градиенте и гессе (которая обычно будет равна нулю для целых чисел), а fminsearch - это не-go, так как вам понадобится действительно хорошая начальная оценка, а ставка конвергенции (примерно) O(N), то есть для этой 20-мерной проблемы вам придется подождать довольно долго до конвергенции, без гарантией поиска глобального решения.

interval method может быть, однако, у меня лично мало опыта с этим. В MATLAB или любом из его наборов инструментов нет родственных материалов, связанных с интервалом, но есть свободно доступный INTLAB.

Так что, если вы не чувствуете, как реализации собственного нелинейного алгоритма двоичного целочисленного программирования, или не в настроении для приключений с INTLAB, есть на самом деле только одна вещь: эвристические методы. В this link существует аналогичная ситуация с контуром решения: используйте генетический алгоритм (ga) из инструментария Global Optimization.

Я бы реализовать эту проблему примерно так:

function [sol, fval, exitflag] = bintprog_nonlinear() 

    %// insert your data here 
    %// Any sparsity you may have here will only make this more 
    %// *memory* efficient, not *computationally* 
    data = [... 
     ... %// this will be an array with size 4-by-20-by-12 
     ... %// (or some permutation of that you find more intuitive) 
     ]; 

    %// offsets into the 3D array to facilitate indexing a bit 
    offsets = bsxfun(@plus, ... 
     repmat(1:size(data,1), size(data,3),1), ... 
     (0:size(data,3)-1)' * size(data,1)*size(data,2)); %//' 

    %// your objective function 
    function val = obj(X) 

     %// limit "X" to integers in [1 20] 
     X = min(max(round(X),1),size(data,3)); 

     %// "X" will be a collection of 12 integers between 0 and 20, which are 
     %// indices into the data matrix 

     %// form "s" from "X"   
     s = sum(bsxfun(@plus, offsets, X*size(data,1) - size(data,1))); 


     %// XxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxX   
     %// Compute the NEGATIVE VALUE of your function here 
     %// XxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxXxX 


    end 

    %// your "non-linear" constraint function 
    function [C, Ceq] = nonlcon(X) 

     %// limit "X" to integers in [1 20] 
     X = min(max(round(X),1),size(data,3)); 

     %// form "s" from "X"   
     s = sum(bsxfun(@plus, offsets, X(:)*size(data,1) - size(data,1))); 

     %// we have no equality constraints 
     Ceq = []; 

     %// Compute inequality constraints 
     %// NOTE: solver is trying to solve C <= 0, so: 
     C = [... 
      40 - s(1) 
      s(2) - 100 
      -20 - s(4) 
     ]; 

    end 

    %// useful GA options 
    options = gaoptimset(... 
     'UseParallel', 'always'... 
     ... 
    ); 

    %// The rest really depends on the specifics of the problem. 
    %// Useful to look at will be at least 'TolCon', 'Vectorized', and of course, 
    %// 'PopulationType', 'Generations', etc. 

    %// THE OPTIMZIATION 
    [sol, fval, exitflag] = ga(... 
     @obj, size(data,3), ... %// objective function, taking a vector of 20 values 
     [],[], [],[], ...  %// no linear (in)equality constraints 
     1,size(data,2), ...  %// lower and upper limits 
     @nonlcon, options);  %// your "nonlinear" constraints 


end 

Обратите внимание, что даже если ваши ограничения являются по существу линейного, путь, по которому вы должны вычислить значение для вашего s требует использования обычая функция ограничения (nonlcon).

Особо следует отметить, что это в настоящее время (возможно) неоптимальный способ использования ga - Я не знаю особенностей вашей целевой функции, поэтому может быть намного больше. Например, в настоящее время я использую простой round() для преобразования ввода X в целые числа, но с использованием 'PopulationType', 'custom' (с обычным 'CreationFcn', и т. Д.) Может привести к лучшим результатам. Кроме того, 'Vectorized', скорее всего, ускорит работу, но я не знаю, легко ли ваша функция векторизации.

И да, я использую вложенные функции (мне просто нравятся эти вещи!); он предотвращает эти огромные, обычно идентичные списки входных аргументов, если вы используете подфункции или автономные функции, и они действительно могут быть повышением производительности, поскольку копирование данных происходит немного. Но я понимаю, что их правила определения области зрения делают их несколько схожими с конструкциями goto, и поэтому они - «не все чашки чая» ... вы можете захотеть преобразовать их в подфункции, чтобы предотвратить длительные и бесполезные дискуссии с ваши сотрудники :)

В любом случае, это должно быть хорошим местом для начала. Дайте мне знать, если это вообще полезно.

+0

+1 для подробного ответа – bla

0

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

Скажите, что вы нашли s s.t. f (s) - максимальные заданные ограничения s, вам все равно нужно выяснить, как построить s с двенадцатью 4-элементными векторами (переопределенная система, если она когда-либо была), где каждый вектор имеет 20 возможных значений.Редкость может помочь, хотя я не уверен, как можно иметь вектор с четырьмя элементами: 70-90% ноль, а разреженность будет полезна только в том случае, если была какая-то еще не описанная методология в том, как вектор организован

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

+0

Тот факт, что '' 'должен быть линейной комбинацией заданных векторов (со всеми коэффициентами единства), на самом деле является самым жестким ограничением в этой задаче. Можно даже сказать, что это * проблема; минимум 'f (s)' только вторичное рассмотрение. –

+0

О разреженности и размере векторов: Мои векторы имеют 24 элемента, я просто использовал более мелкие, чтобы сделать пример. Вероятно, я должен был упомянуть об этом. – Johannes

0

Я знаю, этот ответ доходит до вас действительно late.

К сожалению, проблема, как показывает не так много моделей, которые будут эксплуатироваться, кроме грубой силы ветвпа & Bound, Master & невольника, etc.- Попытки подхода -i.e. рабовладельца решая сначала функцию непрерывной нелинейной задачи как ведущего и решающую дискретную селекцию в качестве подчиненного, могли бы помочь, но с таким количеством комбинаций и без дополнительной информации по векторам не так много места для работы.

Но на основе данных непрерывных функций почти всюду, основанных на комбинациях сумм и операторов умножения и их обратных, разреженность - это явная точка, которая будет использоваться здесь. Если 70-90% векторов равны нулю, почти значительная часть пространства решений будет близка к нулю или близка к infinite. Следовательно, псевдо-решение 80-20 легко отбрасывает «нулевые» комбинации и использует только «бесконечные».

Таким образом, грубая сила может направляться.