параметризованный вектор с помощью PySCIPopt

Я пытаюсь использовать PySCIPopt для решения традиционной проблемы типа ограничений Ax-b +. У меня много значений b, и мне нужно запустить оптимизатор для каждого из них. Как я могу повторно использовать настройку? Второй вопрос: что эквивалентно norm в PySCIPopt? Или как правильно довести Ax-b до нуля, насколько это возможно? Смотрите ??? отметки ниже

import numpy as np
from pyscipopt import Model, quicksum

def make_program():
    A = ... load constant master matrix ...
    model = Model('Match_to_Master')
    x = []
    y = []
    for i in range(A.shape[1]):
        x.append(model.addVar(vtype='C', lb=0.0, ub=4.0, name='x(%s)' % i))
        y.append(model.addVar(vtype='B', name='y(%s)' % i))
        model.addCons(x[i] <= y[i]*4)
    for i in range(0, A.shape[1] - 20, 20):
       model.addCons(quicksum(y[i:i+20]) <= 1)

    #b = Parameter(A.shape[0], nonneg=True) ???
    model.setObjective(norm(A*x - b), sense='minimize') ???
    return b, x, model


def run_program(data, thresh=0.2):
    b, x, model = make_program()
    B = ... from data load matrix for analysis ...
    c = 0
    for column in B.T:
        b.value = column   ???
        model.optimize()  # reuse previous x values as starting point
        x.value[x.value < thresh] = 0.0
        for i in range(0, x.value.size - 20, 20):
            sum = np.sum(x.value[i:i+20])
            if sum > 0.2:
                print('  hit at ', str(i//20), ' for column ', str(c))
        c += 1

person Brannon    schedule 11.07.2019    source источник
comment
В общем, l1-норма потребует некоторой линеаризации функции abs (LP-совместимой), в то время как l2-норма - это проблема QP или (лучше) SOCP, и для этого потребуется решатель.   -  person sascha    schedule 13.07.2019


Ответы (1)


Что ж, я не думаю, что то, что вы решаете, является традиционной задачей типа Ax-b. Традиционный тип - min | Ax - b | ^ 2. Для этого вам не нужен SCIP.

Если вы хотите решить программу для нескольких значений b, вам следует соответствующим образом адаптировать функцию make_program. Кажется, вас интересуют интегральные решения. В то время как LP имеют некоторые возможности горячего запуска для разных правых сторон (с использованием двойного симплекса), целочисленные программы не имеют такой возможности.

Тем не менее, вы можете использовать решение из предыдущего запуска. Использовать

solution = model.createSol(None)

# for each nonzero variable
model.setSolVal(solution, var, val)

model.trySol(solution)

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

Что касается минимизации нормы, как указывает sascha, вы можете минимизировать только многогранные нормы, скорее всего, норму 1 / inf. В первом случае вы можете использовать

мин y_1 + ... + y_m + z_1 + ... + z_m

s.t. A x + y - z = b

y, z >= 0

В последнем вы можете использовать

разум

s.t. A x + y - z = b

y, z >= 0

d >= y_i, i=1,...,m

d >= z_i i=1,...,m

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

Изменить:

Что касается 2-нормы: во-первых, обратите внимание, что минимизация 2-нормы превращает проблему в смешанную целочисленную нелинейную программу (MINLP), в частности, в смешанную целочисленную квадратичную задачу (MIQP). Таким образом, проблему становится труднее решить :) Я не знаю, подходит ли SCIP в данном случае (я слышал хорошие отзывы о Pajarito). Тем не менее, SCIP может решать MINLP.

Чтобы смоделировать норму, вы должны отметить, что SCIP не поддерживает нелинейные цели, а только нелинейные ограничения. Таким образом, модель должна быть

мин y_1 + ... + y_m

s.t. (A x - b)^2 - y <= 0

Вы должны иметь возможность добавлять произвольные нелинейные ограничения с помощью model.addConst(...). В зависимости от степени выражения PySCIPOpt поступит правильно. В случае ограничения степени 2, т. Е.

model.addCons(x[i]*y[i] <= 0)

это означает добавление квадратичного ограничения.

Обратите внимание, что это минимизирует квадрат нормы вместо нормы. Вы получите правильное решение, но объективное значение будет неверным.

Кроме того, я думаю, что это поможет или, возможно, потребуется для компиляции SCIP с включенной поддержкой Ipopt. Ipopt - это решатель для нелинейных программ (NLP), которые являются релаксацией MINLP.

person hfhc2    schedule 14.07.2019
comment
Я отлично решаю min |Ax - b|^2. Как это сделать с помощью готового решения? - person Brannon; 15.07.2019