Как стать автором
Обновить

Комментарии 84

цитирую вас:

Задача о сумме подмножеств в общей формулировке звучит так:

Существует множество S чисел. Вопрос состоит в том, будет ли сумма некоторого подмножества от S равна заданному числу Т.

Известно, что данная задача NP-полная.

в вашей ссылке на вики https://ru.wikipedia.org/wiki/Задача_о_сумме_подмножеств

написано :

Вычислительная сложность задачи о сумме подмножеств зависит от двух параметров — числа N элементов множества и точности P (определяется как число двоичных разрядов в числах, составляющих множество). Замечание: здесь буквы N и P не имеют ничего общего с классом задач NP.

обратите внимание на Замечание: здесь буквы N и P не имеют ничего общего с классом задач NP.

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

Так как для чисел порядка миллиарда, число двоичных разрядов практически не имеет значения, то я несколько упростил повествование.

Эта задача действительно NP-полная.

То, что задача NP - полна коллега вроде и не оспаривал. Видимо он хотел указать, что сложность задачи зависит не только от числа элементов (N), но и от точности (P) и это не имеет отношение к NP-полноте. Я не верно понял?

Мне показалось, что SemenovVV увидел по ссылке N, P и подумал, что вы ошибочно интерпретировали их как NP. Другого прочтения этого комментария я придумать не могу.

Только в действительных числах, а такого в CS не бывает.
Для натуральных чисел (можно обобщить и для рациональных), сложность решения M*N, где N - число элементов, M - сумма элементов множества.
https://leetcode.com/discuss/study-guide/1152328/01-Knapsack-Problem-and-Dynamic-Programming

Но в общем случае M = O(2^N). Ведь в CS сложность считается от размера входных данных. И если у вас одно число на N бит, то вот у вас уже M=2^N. Если же вводить ограничения на числа, то это уже немного другая задача. Но все-равно, если там 64-битные числа, то этот метод не применим, ибо массива на 10^18 элементов вы не заведете. Да и если заведете, работать это будет не сильно быстро.

Да, точно. Я всё время забываю, что есть быдлокодерская сложность, где каждое сложение-умножение это O(1), а цикл for i=1 to N имеет сложность O(N), и есть сложность для адептов Машин Тьюринга, где тот же цикл - уже NP-полная задача, потому что число итераций растёт как степень двойки от длины входа.

Ну смотрите, задача разложения числа на множители тоже решается одним циклом. Перебрали все числа до корня - поделили. Даже не за линию, а за O(Sqrt(N)). Сублинейный алгоритм! Но вы же не будете спорить, что задача факторизации - NP-полная? Или вся криптограия идет лесом.

Вот тут тоже самое же.

Я и говорю, что есть две параллельные терминологии.
Часто, говоря о сложности задачи, считают её не от длины записи, а от значения входа.

  • Какая сложность поиска элмемента в массиве длины N?

  • O(N)

  • Какая сложность наивной факторизации числа N?

  • O(N^2).

Я и говорю, что есть две параллельные терминологии.

Нет. Терминология одна. У вас непонимание, что такое NP и его братья.

Вернемся к началу дисскуссии:

Я написал, что задача NP-полная, вы возразили, что есть решение за O(NM). Вообще, определения там запутанные, так что давайте посмотрим на определение P-задач:

множество задач, для которых существуют «быстрые» алгоритмы решения (время работы которых полиномиально зависит от размера входных данных)

Вот эти NP-полные - это те, для которых неизвестны "быстрые" алгоритмы, и их скорее всего их не существутет, потому что решение для любой такой задачи автоматически "быстро" решает их все. Обратите внимание на выделенные жирным слова в определении. Надо считать ассмтотику относительно размера входных данных. Эти N, M для ДП - не размер входных данных. Собственно, M будет экспоненциально от этого размера.

С практической точки зрения, вы можете на основе этой задачи сделать шифр, где публичный ключ - множество, а приватный - подмножество с заданной суммой. И жалкий 4096-битный ключ вам этим ДП придется взламывать тысячи лет. Хотя для действительно полиномиальных задач любой 4-килобайтный входной файлик будет пережеван весьма быстро.

Так что задача все еще остается "сложной" и она действительно является NP-полной.

Есть же решение через динамическое программирование за O(N sum). Надо завести массив размером с искоммую сумму, и там помечать те суммы, которые можно получить. Там же можно хранить, какой элемент множества можно взять.

Для вот таких вот чисел в задаче оно найдет решение за доли секунды. Сумма порядка 10 миллионов. С другой стороны, для 24 чисел можно полным перебором перебрать ~8 миллионов вариантов (2^24), что будет даже быстрее.

Имеет смысл применять что-то иное, вроде линейного программирования, когда у вас числа большие и их самих больше. Как долго ваше решение работает, если чисел 50, а сумма их - порядка миллиарда?

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

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

Если вы приложите набор тестовых данных, то могу протестировать.

Динамическое программирование требует очень много память на больших N и P.

Нет, от N оно линейно. ДП страдает, когда числа большие. А когда N маленькое - отлично работает полный перебор.

Вот и получается, что единственная ниша, где ваш метод может побороться - это относительно большое (>=30) количество больших (~миллиард) чисел. Есть ли там выигрыш?

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

n = 100
val = 77019193
source = [1885440, 1403958, 1794772, 1933488, 1441001, 1042450, 1271493, 1536110, 1509532, 1424604, 1962838, 1821872, 1870163, 1318046, 1499748, 1375441, 1611720, 1934973, 1952225, 1229053, 1529202, 1146039, 1295528, 1146534, 1792518, 1099437, 1648406, 1838234, 1262674, 1953938, 1558433, 1739426, 1849574, 1631140, 1945989, 1154100, 1325213, 1103560, 1765284, 1077324, 1942500, 1891786, 1717209, 1346236, 1495077, 1587007, 1105592, 1370977, 1455262, 1331556, 1640561, 1671532, 1957361, 1214410, 1579363, 1500181, 1464197, 1907343, 1546678, 1273145, 1065304, 1844132, 1963080, 1575352, 1960489, 1014723, 1097802, 1754665, 1880899, 1418196, 1744754, 1864912, 1823182, 1700609, 1655638, 1001198, 1641620, 1517553, 1868287, 1909747, 1349317, 1255759, 1765752, 1341001, 1737822, 1912755, 1066043, 1200348, 1961564, 1595078, 1232473, 1250206, 1842368, 1149416, 1842194, 1569366, 1469730, 1095646, 1084353, 1335601]

solution optimized = [1794772, 1933488, 1536110, 1509532, 1962838, 1821872, 1870163, 1611720, 1934973, 1952225, 1792518, 1648406, 1838234, 1953938, 1558433, 1631140, 1945989, 1942500, 1891786, 1717209, 1640561, 1671532, 1957361, 1579363, 1907343, 1546678, 1844132, 1963080, 1960489, 1880899, 1864912, 1823182, 1700609, 1641620, 1517553, 1868287, 1909747, 1765752, 1912755, 1961564, 1842368, 1842194, 1569366]
time = 4.474821

Видимо, линейная алгебра сильно продвинулась из-за развития ML.

Потому что, если решать в лоб, "по-школьному", дискриминант матрицы на 30 строк из таких больших чисел переполнит любые фиксированные типы. И ещё проблема, что решение должно быть точным, а любая микро-погрешность не найдёт решение Xi in {0,1}, посчитается какой-нибудь Xi=0.99999 и всё, вариант отсечётся. То есть, я думаю, тщательное тестирование на очень сложных случаях покажет, что этот метод не всегда даёт правильный результат.

Тут не линейная алгебра же. Это же целочисленное программирование.

Очень интересно. Я вообще не знал, что есть такое направление. Есть ли учебники, или это всё на уровне исследований?

Есть учебники. Это весьма зрелая область. Первое, что нагуглил. Вообще, начинают изучение не с целочисленного, а с простого линейного программирования. Симплекс методу из этой области уже лет 70. В целочисленных решателях часто используют итеративный подход со всякими отсеченями, что есть в общем-то надстройка над этим.

Симплекс-метод это старая тема, которую в ВУЗах дают даже не математикам. А вот в целочисленной области что придумали, это интересно.

Смутно помню метод, который запускает симплекс, получает нецелочисленное решение, добавляет ограничение, которое отсекает вот этот лишний оптимум. И так далее, пока не найдется целое решение. Вроде это оно. Ему лет 60. Что там нового в этой теме напридумывали - не в курсе.

Этот алгоритм это есть в книжке, которую вы привели выше.
Но там те проблемы, которые я упомянул - требуется работа в поле рациональных чисел с большими знаменателями.
Там дальше в книжке рассмотрен алгоритм, работающий только в целых числах, пока я ещё не докурил, свободен ли он от этих проблем.

И ещё проблема, что решение должно быть точным, а любая микро-погрешность не найдёт решение Xi in {0,1}, посчитается какой-нибудь Xi=0.99999 и всё, вариант отсечётся. То есть, я думаю, тщательное тестирование на очень сложных случаях покажет, что этот метод не всегда даёт правильный результат.

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

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

А динамическое программирование запускали? Сколько оно работает?

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

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

Да, вы какой-то странный алгоритм реализовали. Или numpy не использовали, без него массив 77М*100 очень много памяти занимает. Код ниже на вашем тесте отработал за 1.6 секунд. Запустите на вашей машине, чтобы сравнить.

import numpy as np
import time

def find_subset_sum(numbers, target):
    """
    Finds a subset of given numbers that sums up to the target using dynamic programming.
    :param numbers: List of integers
    :param target: Target sum
    :return: Subset of numbers that adds up to the target (or None if not possible)
    """
    n = len(numbers)
    dp = np.zeros(target + 1, dtype=bool)
    dp[0] = True

    for i in range(1, n + 1):
        dp |= np.roll(dp, numbers[i - 1])

    if not dp[target]:
        return None

    subset = []
    j = target
    for i in range(n, 0, -1):
        if j >= numbers[i - 1] and dp[j - numbers[i - 1]]:
            subset.append(numbers[i - 1])
            j -= numbers[i - 1]

    return subset[::-1]

target_sum = 77019193
source = [1885440, 1403958, 1794772, 1933488, 1441001, 1042450, 1271493, 1536110, 1509532, 1424604, 1962838, 1821872, 1870163, 1318046, 1499748, 1375441, 1611720, 1934973, 1952225, 1229053, 1529202, 1146039, 1295528, 1146534, 1792518, 1099437, 1648406, 1838234, 1262674, 1953938, 1558433, 1739426, 1849574, 1631140, 1945989, 1154100, 1325213, 1103560, 1765284, 1077324, 1942500, 1891786, 1717209, 1346236, 1495077, 1587007, 1105592, 1370977, 1455262, 1331556, 1640561, 1671532, 1957361, 1214410, 1579363, 1500181, 1464197, 1907343, 1546678, 1273145, 1065304, 1844132, 1963080, 1575352, 1960489, 1014723, 1097802, 1754665, 1880899, 1418196, 1744754, 1864912, 1823182, 1700609, 1655638, 1001198, 1641620, 1517553, 1868287, 1909747, 1349317, 1255759, 1765752, 1341001, 1737822, 1912755, 1066043, 1200348, 1961564, 1595078, 1232473, 1250206, 1842368, 1149416, 1842194, 1569366, 1469730, 1095646, 1084353, 1335601]
start = time.time()
result = find_subset_sum(source, target_sum)
end = time.time()
elapsed_time = end - start
print(f"Execution time: {elapsed_time:.6f} seconds")
print(f"Subset with sum {target_sum}: {result}")

Можно извратиться с битовыми оптимизациями и раз в 30 ускорить.

Алгоритм ДП, что я рассматривал, очень неэффективен по сравнению с тем, что предоставили Вы.

Беру свои слова назад, в такой реализации ДП на многих наборах обходит мою версию ЛП.

Любопытно, на некоторых наборах сильно лидирует ЛП, но если проигрывает, то часто значительно, тогда как ДП более стабилен. Как вариант интересно устраивать гонку алгоритмов, какой финиширует первым.

Однако меня огорчил тот факт, что на больших наборах входных данных, памяти ДП всё же не хватает, и он не выполняется с ошибкой:

MemoryError: Unable to allocate 14.2 GiB for an array with shape (15207959671,) and data type bool.

n = 200
t = 15207959670
source = [118034063, 176397250, 108470054, 134234785, 115826780, 166496171, 160329669, 163383683, 187455328, 150951092, 128179657, 112597620, 165479012, 103804733, 152319252, 158085012, 181528947, 100282669, 193393106, 159778857, 135746282, 196843463, 130703945, 179343270, 113720696, 142604684, 104105718, 102996023, 103415285, 187180606, 172667152, 101235465, 151164366, 192138303, 129071478, 156655527, 197422287, 103897788, 170817221, 129754951, 158772277, 166546792, 174203556, 131284065, 146399124, 130986382, 190845073, 129364293, 161686932, 138893829, 102884299, 155858725, 174686034, 186207290, 113421809, 124951916, 184470316, 197125183, 139780845, 116225575, 199743456, 144653591, 196835997, 195454543, 167216197, 156654242, 168144655, 189966890, 125481199, 140717432, 138139224, 178863733, 167023240, 167818046, 152795029, 179054544, 104633978, 164454973, 132580007, 199821838, 154262629, 155608283, 189220365, 123220660, 149274526, 173658522, 194360533, 190527955, 199081602, 150291788, 111605483, 158916432, 189088064, 168239848, 114486288, 121971213, 169919170, 152781805, 149730710, 165725551, 198350162, 103969484, 162991083, 105836765, 141410118, 194406345, 182518497, 179615772, 177601456, 152828055, 186859830, 122863882, 122628343, 167409318, 130459014, 101651090, 126778633, 172426227, 173596743, 131162152, 154285013, 168957265, 146147529, 177550306, 147415655, 161623617, 136142079, 188478314, 173550819, 181731190, 197898435, 100766266, 151497950, 199388685, 168786576, 117347564, 169615820, 175344177, 127579764, 157188922, 107532741, 164572392, 148954043, 176504015, 174410468, 126821992, 167742434, 155485614, 165085546, 147887538, 155623117, 146449792, 100212701, 172273400, 172492277, 183683337, 182201978, 144444516, 161491422, 180511200, 103754738, 130817065, 185278064, 123784892, 173921254, 178445010, 124264416, 112294532, 173957878, 134264986, 104356590, 190343768, 109456105, 111171496, 102240178, 160800468, 101954206, 137741579, 133495272, 136056483, 114695314, 183859516, 124777959, 146227654, 138961302, 109330196, 122477484, 121424575, 134254527, 170783798, 122568032, 188134944, 136629955, 187000307, 195507983, 139526151, 161029019, 194304805, 143218345, 166638250]

solution optimized = [176397250, 166496171, 160329669, 163383683, 187455328, 165479012, 181528947, 193393106, 159778857, 179343270, 187180606, 172667152, 192138303, 197422287, 170817221, 166546792, 174203556, 190845073, 161686932, 174686034, 186207290, 184470316, 197125183, 199743456, 196835997, 195454543, 167216197, 189966890, 178863733, 167023240, 179054544, 164454973, 199821838, 189220365, 173658522, 194360533, 190527955, 199081602, 158916432, 189088064, 169919170, 165725551, 198350162, 194406345, 182518497, 179615772, 177601456, 186859830, 167409318, 172426227, 173596743, 168957265, 177550306, 188478314, 173550819, 181731190, 197898435, 199388685, 168786576, 169615820, 175344177, 164572392, 176504015, 174410468, 155485614, 172273400, 172492277, 183683337, 182201978, 180511200, 185278064, 173921254, 178445010, 173957878, 190343768, 160800468, 183859516, 170783798, 188134944, 187000307, 195507983, 161029019, 194304805, 143218345, 166638250]

Да, я про это ограничение и писал. Если числа большие, то ДП вообще не работает. Но вот вопрос в том, хорошо ли в этом случае работает линейный решатель?

Не утверждаю, что решение с ЦЛП - серебряная пуля.

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

Немного доработал представленный алгоритм динамического программирования для Python.

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

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

Возможно, кому то пригодится.

import bitstring as bit
import random as rnd
import time

def find_subset_sum(numbers: list, target: int):
    dp = bit.BitArray(target + 1)
    dp[0] = 1

    for i in numbers:
        temp = dp.copy()
        dp.ror(i)
        dp |= temp
        if dp[target] == 1:
            break

    if not dp[target]:
        return []

    subset = []
    j = target
    for i in numbers[::-1]:
        if j >= i and dp[j - i]:
            subset.append(i)
            j -= i

    return subset[::-1]


n = 100
# rnd.seed(0)
source = [rnd.randint(1000000, 9999999) for i in range(n)]
print('source:', source)
target_sum = sum(source) // 2
print('target_sum:', target_sum)

print()
start = time.time()
result = find_subset_sum(source, target_sum)
end = time.time()
elapsed_time = end - start
print(f"Execution time: {elapsed_time:.6f} seconds")
print(f"Subset with sum {target_sum}: {result}")

Вот только тут большая ошибка. Он правильно найдет, есть ли решение, или нет. Но ваше восстановление ответа не работатет. Потому что нельзя брать число i, только если dp[j-i]. Как пример, у вас есть набор чисел {4,2}, надо набрать 4. Если вы сначала посмотрите i=2, то вы его возьмете. И как бы вам пришлось бы брать 2 два раза, но это делать нельзя же. Понятно, что реальный пример должен быть сильно сложнее, ведь набирать надо половину всех чисел. Но подобные ситуации могут и там возникнуть.

Чтобы восстановить ответ, посмотрите мой комментарий с кодом на c++. Там пара опечаток, почитайте всю ветку, но для восстановления ответа вам нужен массив чисел, а не битов.

Если вам только проверить надо на существование ответа, а не сам список получить, то массив битов работает отлично.

Так понимаю достаточно сделать: numbers.sort() первым шагом функции.

Поскольку обход цикла идёт с конца списка, проблема полностью устраняется.

source: [2, 2, 2, 2, 4, 4, 6]
target_sum: 8

Subset with sum 8: [2, 6]

Нет, не достаточно. Смотрите, что у вас в dp[] хранится по итогу? У вас там пометка, можно ли набрать вот такую вут сумму всеми числами или нет. Только по этой информации никак не восстановить множество, ведь вы не знаете, какими числами каждую сумму можно набрать, и поэтому не знаете, можно ли брать заданное число.

И не важно, с какого конца их брать. Вот контр-пример к сортировке:

Вам надо набрать 8, есть числа {2,3,3,4}. Первым шагом вы возьмете 4, ведь можно же набрать 4 всеми этими числами. И дальше все. Оставшиеся 4 вы никак из {2,3,3} не наберете.

Вам надо набрать 8, есть числа {2,3,3,4}. Первым шагом вы возьмете 4, ведь можно же набрать 4 всеми этими числами. И дальше все. Оставшиеся 4 вы никак из {2,3,3} не наберете.

Это не так! Проиллюстрирую:

Так, т.е. берем от младших чисел к старшим? Ну вот пример:

Надо набрать 6. Есть {3,6}. Вы первым делом возьмете 3, потому что 3 набрать можно. И вот уже не правильно. Более сложный пример, когда надо половину набрать, как в исходной задаче из статьи: {3, 5, 6, 7, 11}. Надо набрать 16. Проверяем 3. 16-3=13=6+7 взять можно, вы беретe 3. Осталось набрать 13. 13-5=8 набрать можно? Можно: 3+5. Вы берете 5. Осталость набрать 8. А никак вы 8 оставшимеся {6,7,8} не наберете. Его можно набрать 3+5, но вы из оба уже использовали. А так-то ответ должен быть, например, 3+6+7 или 5+11.

Хоть минимальные сначала берите, хоть максимальные, хоть в случайном порядке числа проверяйте, все равно будет тест, который вам все сломает.

Более сложный пример, когда надо половину набрать, как в исходной задаче из статьи: {3, 5, 6, 7, 11}.

Обратите внимание на условие:

for i in numbers:
  temp = dp.copy()
  dp.ror(i)
  dp.set(0, (range(i)))    
  dp |= temp
  if dp[target]:
      break

Оно заранее завершает процесс заполнения dp, и исключает во многом те ситуации, что вы описывали в предыдущем ответе.

В прошлом возражении вы начинали брать числа с меньших (для теста {2,3,3,4}, а тут берете с больших. Как так-то?

Ну хорошо, если вы рассуждения не принимаете и с завидным упорством цепляетесь за частности, то запустите ваш код на ([2, 3, 3, 4, 52, 52, 100], 108), ([110,70,60,50,30,1,1], 161) и ([2,2,1,1], 3). И в конце убедитесь, что сумма ответа равна искомой target_sum. Ни в одном из этих трех тестов ответ не сойдется ([3,4,100], [1,1], [60,50,30,1,1]). Конкретнее некуда.

например так
def test(numbers: list):
   total_sum = sum(numbers)
   if total_sum % 2 == 1: 
     print ("!!! Odd sum " + str(total_sum))
   ans = find_subset_sum(numbers, total_sum // 2)
   if sum(ans) != total_sum // 2:
     print( "Wrong answer for " + str(numbers) + " : " + str(ans))

Оно заранее завершает процесс заполнения dp, и исключает во многом те ситуации, что вы описывали в предыдущем ответе.

Но не исключает их в принципе, все еще оставляя алгоритм неправильным. Всегда можно все числа умножить на 10 и добавить в конец две единицы. И тогда, пока вы одну из этих единиц не обработаете, вы нечетную нужную сумму никак не наберете, а значит все интересные числа обработаются. Если же числа по возрастанию сортируются, то можно в конец добавить несколько больших чисел, так что все интересные для ошибки числа все-равно обработаются, как в первом примере.

Я привел тесты, где числа и по возрастанию и по убыванию упорядочены. Похожие тесты можно сделать, даже если вы в восстановлении ответа будете числа не с конца брать, а с начала списка (без [::-1] в последнем цикле). Например, на тесте ([2, 3, 3, 4, 52, 52, 100], 108) исправленный так код вернет [52,52,2], что не дает нужную сумму.

Какую бы вы эвристику не придумали, пока вы ее не докажете - она по умолчанию не работает. Вы должны четко и с математичесской строгостью показать, почему эта эвристика делает так, что dp[i-j] можно будет собрать не используя j и уже ранее взятые числа.

Благодарю Вас в вашем упорстве находить истину и желании помочь другим в этом же.

Действительно в моих первоначальных рассуждениях закралась ошибка: Если мы усекаем список просматриваемых чисел на первом этапе, то его и на втором этапе нужно усекать в том же объёме.

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

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

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

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

import bitstring as bit

def find_subset_sum(numbers: list, target: int):    
    numbers.sort()
     
    dp = bit.BitArray(target + 1)
    dp[0] = 1
    count = -1
    
    for i in numbers:
        temp = dp.copy()
        dp.ror(i)
        dp.set(0, range(i)) 
        dp |= temp
        print('i:', i, dp.bin)
        count += 1
        if dp[target] == 1:
            break

    if not dp[target]:
        return []
    print()

    subset = []
    j = target
    for i in numbers[count::-1]:   
        if j >= i and dp[j - i]:
            print('i:', i, 'j:', j, '- берём')
            subset.append(i)
            j -= i
        else:
            print('i:', i, 'j:', j, '- не берём')
    return subset[::-1]


source = [2, 3, 3, 4, 52, 52, 100]
target_sum = 108

print('source:', source)
print('target_sum:', target_sum)
print()

result = find_subset_sum(source, target_sum)
print('result:', sum(result), result)

Ну вот новый контрпример: [2, 3, 3, 4, 52, 52, 100, 1000, 1000]

Вы опять придумали заплатку, которая никак логически не может влиять на основную проблему. Еще раз, проблема в том, что можно получить ситуацию, что вам надо набрать 8, а у вас остались числа [2, 3, 3, 4]. И тут вы никак не можете понять, что 4 брать нельзя. Или что-то подобное. Вы сделали так, что нельзя брать 100 в этом конкретном примере, чтобы попасть в 108-100=8? Ну вот я добавил пару 1000 в массив и ваша эвристика сломалась.

Пока вы логически не покажете, почему вы исключаете ситуацию, что dp[i-j] требует использовать j в сумме - вы не правы. В общем случае, а не затыкая конкретные контр-примеры.

И я не вижу никакого способа по лишь одной строке dp в конце этого избежать.

Ну вот новый контрпример: ([2, 3, 3, 4, 52, 52, 100, 1000, 1000], 108) ...

Ну вот я добавил пару 1000 в массив и ваша эвристика сломалась.

Где же сломалась? Работает, проверьте исходные данные.

Еще раз, проблема в том, что можно получить ситуацию, что вам надо набрать 8, а у вас остались числа [2, 3, 3, 4]

Я понимаю резон в ваших словах, но и вы поймите, что смысл всех танцев как раз в том, чтобы исключить возможность выбора тех величин, что приведут к тупику в дальнейшем.

И да вы правы нужно доказать логически, а это не простая задача.

Где же сломалась? Работает, проверьте исходные данные.

Вы сумму до 1108 исправили? В этом случае ваш код набирает 1106.

И сортировка тут вообще никак не поможет. В [2,3,3,4]->8 вы можете случайно взять 4, что делать нельзя. В [2,1,1] -> 3 вы можете взять 1 и потом опять 1. Тут все еще хуже, ибо тут надо использовать не ту 1, которую вы берете сейчас, а что-то, что вы рассмотрели еще раньше. При этом сами числа могут быть как больше, так и меньше текущего.

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

Мне не остаётся ничего кроме, как проиндексировать биты dp, тогда восстановление идет всегда правильно. Памяти требует, как в оригинальном решении, при длине входного массива до 256 элементов, и в двое больше в ином случае.

import numpy as np

def find_subset_sum(numbers: list, target: int):
    dp = np.zeros(target + 1, dtype='uint8' if len(numbers) < 256 else 'uint16')
    dp[0] = len(numbers)

    for i, val in enumerate(numbers):
        temp = np.roll(dp, val)
        temp[:val] = 0
        temp[temp > 0] = len(numbers) - i
        dp = np.maximum(dp, temp)
        if dp[target]:
            break        

    if not dp[target]:
        return []

    subset = []
    j = target   
    while j > 0:
        i = numbers[len(numbers) - dp[j]]
        subset.append(i)
        j -= i

    return subset[::-1]


source = [2, 3, 3, 4, 52, 52, 100, 1000, 1000]
target_sum = 1108
print('source:', source)
print('target_sum:', target_sum)
print()
result = find_subset_sum(source, target_sum)
print('result:', sum(result), result)

Да. Это, похоже, лучшее по памяти решение. И если вы будете его сравнивать с линейным решателем, то надо делить время на 10. Ведь решатель-то написан на каком-нибудь с++, а из питона вы лишь вызваете библиотеку. А циклы в питоне, да и вот это перекладывание из одного numpy массива в другой - это медленно. Там еще массив temp будет выделятся на каждй итерации, да и всякие конструкции вроде temp[temp>0] похоже выделяют память под список индексов. А решение на С++ будет вообще без лишних выделений памяти и сильно быстрее.

Я придумал, кажется, весьма интуитивное его объяснение этому решению. Очевидно, что можно хранить прямоугольный битовый dp[][] - не одну строчку, а всю матрицу. По ней восстановление ответа легко сделать - можно взять число j, если dp[i-j][j-1]. Потом можно заметить, что каждый столбец там начинается с нулей, но где-то станет равен 1 и дальше будут только 1. поэтому можно вместо прямоугольной матрицы хранить строку чисел, где помечать, где именно начинаются 1 в каждом столбце. У вас в решении примерно это и хранится, только перевернутое - вы храните, сколько там единиц в столбце. Мой код на С++ где-то в комментах тут хранит само число с номером равным этой строке.

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

И если вы будете его сравнивать с линейным решателем, то надо делить время на 10. Ведь решатель-то написан на каком-нибудь с++, а из питона вы лишь вызываете библиотеку

Ну это такой аргумент, numpy если я не ошибаюсь написан на том же си и так же использует SIMD векторизацию. Так что разница будет не столь значительна.

У вас в решении примерно это и хранится, только перевернутое - вы храните, сколько там единиц в столбце.

Не совсем так, я храню индекс числа в оригинальном массиве, в перевёрнутом виде. И оставляю больший индекс проходя функцией max (оставляем более старый) чтоб не затереть уже найденный маршрут.  Мне такой подход напоминает алгоритм Дейкстры для поиска пути в графе.

Хотелось бы сравнить алгоритм с вашей реализацией по скорости.

Ну это такой аргумент, numpy если я не ошибаюсь написан на том же си и так же использует SIMD векторизацию.

Ну если бы в numpy была операция для всех вычислений, то да. Но все, что вы из встроенного используете - это np.roll. А дальше там пересылка индексов через питон. Всякие максимумы, которые не факт что векторизованы.

Не совсем так, я храню индекс числа в оригинальном массиве, в перевёрнутом виде.

Ну так это же и есть количество единиц. Индекс в перевернутом виде - номер строки с конца. Вот это взятое число - это та строка, с которой начинаются единицы. Вот и получается, что номер числа с конца это и есть количество единиц. Ну, может там +1 откуда-то получается везде, но это не принципиально.

Хотелось бы сравнить алгоритм с вашей реализацией по скорости.

Попозже сравню у себя питоновскую реализацию с си++.

Я был не прав. Мое изначальное решение на C++, которое шло по массиву задом-наперед, было гораздо медленнее. 50с вместо 10с у numpy. Развернув цикл, оно стало работать чуть быстрее - за 35с.

Действительно, векторизация внутри numpy перевешивает накладные расходы из-за питона.

Сравнивалась реализация с bitstring? Тогда может, дело в нагрузке на память?

Да, не, там не в памяти дело. Вообще, 1000 раз пройтись по массиву размером 80м, это 80 миллиардов операций. В секунду влезает 2-5 миллиардов. Ну вот и получается. что где-то 16-40 секунд оно и должно работать. Операции очень непросытые, ибо там условие, несколько проверок и присваение. Так что ближе к 40с. С векторизацией в 8 раз быстрее из-за векторизации, но в пару раз медленнее из-за дополнительных проходов (надо там отдельно roll сделать, отдельно занулить что-то и отдельно максимумы взять). Так что все более менее сходится.

Попробуйте изменить на это:

for i, val in enumerate(numbers):
  k = np.uint8(len(numbers) - i) if len(numbers) < 256 else np.uint16(len(numbers) - i)
  temp = np.roll(dp, val)
  temp[:val] = 0
  np.maximum(dp, temp.astype(bool) * k, dp)
  if dp[target]:
      break

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

Всякие максимумы, которые не факт что векторизованы.

np.maximum 100% векторизированна, там внутри должна быть интерпретация к PMAXUB/PMAXUW. А в зависимости от процессора там прирост может быть значительный.

12.5 секунд. На 25% медленнее, чем изначальное решение. Выделение этого k все портит. Ну и я не очень понимаю, где вы тут один проход убрали. Раньше у вас было temp[temp > 0] = len(numbers) - i , вместо него стало temp.astype(bool)*k

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

С последним вариантом неувязка вышла. Зато новая попытка в среднем 10% выигрыш по времени.

Основная идея в том, что не нужно трогать ту часть массива, которая остаётся неизменной на каждом этапе.

for i, val in enumerate(numbers):
    temp = np.roll(dp, val).view()[val:]
    temp[temp > 0] = len(numbers) - i
    dp_slice = dp[val:]
    np.maximum(dp_slice, temp, dp_slice)
    if dp[target]:
        break

А ларчик просто открывался!

Прирост скорости от 50% - 100% от базового решения.

Мне только так и не нравится, что по факту три прохода по массиву.

1)     Копирование части массива

2)     Присваивание значимых элементов

3)     Максимизация

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

Как понимаю на C можно в один проход сделать, но мне кажется всё равно медленнее будет чем текущая версия с векторизацией.

import random as rnd
import numpy as np
import time

def find_subset_sum(numbers: list, target: int):
    dp = np.zeros(target + 1, dtype='uint8' if len(numbers) < 256 else 'uint16')    
    dp[0] = 1
    offset = 1

    for i, val in enumerate(numbers):
        offset += val
        if offset > target + 1:
            offset = target + 1
        temp = dp[:offset-val].copy()
        temp[temp > 0] = len(numbers) - i
        np.maximum(dp[val:offset], temp, dp[val:offset])
        if dp[target]:
            break

    if not dp[target]:
        return []

    subset = []
    j = target   
    while j > 0:
        i = numbers[len(numbers) - dp[j]]
        subset.append(i)
        j -= i

    return subset[::-1]


n = 50
rnd.seed(1)
source = [rnd.randint(1000000, 9999999) for i in range(n)]

target_sum = sum(source) // 2
print('source:', source)
print('target_sum:', target_sum)
print()
start = time.time()
result = find_subset_sum(source, target_sum)
end = time.time()
elapsed_time = end - start
print(f"Execution time: {elapsed_time:.6f} seconds")
print('result:', sum(result), result)

Оставлю тут для истории.

Если заменить сравнение:

temp[temp > 0] = len(numbers) – i

на векторный аналог:

np.clip(temp, None, len(numbers) - i , temp)

то можно получить ещё 60% прироста производительности.

import random as rnd
import numpy as np
import time

def find_subset_sum(numbers: list[int], target: int) -> list[int]:
    """
    Задача о сумме подмножеств натуральных чисел (Subset sum problem - SSP)
    Решение с использованием алгоритма динамического программирования
    """
    dp = np.zeros(target + 1, dtype='uint8' if len(numbers) < 256 else 'uint16')
    dp[0] = len(numbers)
    offset = 1

    for i, val in enumerate(numbers):
        offset += val
        if offset > target + 1:
            offset = target + 1
        temp = dp[:offset-val].copy()   
        np.clip(temp, None, len(numbers) - i , temp)
        np.maximum(dp[val:offset], temp, dp[val:offset])
        if dp[target]:
            break
        
        if dp[target]:
            break        
            
    if not dp[target]:
        return []

    subset = []
    j = target   
    while j > 0:
        i = numbers[len(numbers) - dp[j]]
        subset.append(i)
        j -= i

    return subset[::-1]


n = 50
rnd.seed(1)
source = [rnd.randint(1000000, 9999999) for i in range(n)]

target_sum = sum(source) // 2
print('source:', source)
print('target_sum:', target_sum)
print()
start = time.time()
result = find_subset_sum(source, target_sum)
end = time.time()
elapsed_time = end - start
print(f"Execution time: {elapsed_time:.6f} seconds")
print('result:', sum(result), result)

Цикл 70млн*100 итераций, который имеет примерно такой же футпринт по памяти, работает 4516 ms на моей машине.

Динамическое программирование требует очень много память на больших N и P

Нет, в этой задаче таблица заполняется сверху вниз, и нужно хранить только последнюю заполненную строку для расчёта следующей.

Там все сложно с восстановлением ответа. Если вам надо знать только есть ли разбиение, то такая оптимизация отлично работает. А если надо само подмножество вывести, то все совсем плохо. Ибо вы не знаете, какие числа вам понадобятся потом и какие были использованы ранее. Так что во время восстановления можно прийти в ситуацию, что вычтя какое-то число вы оказываетесь в сумме, которую набрать нельзя. Ну, или вы одно и то же число несколько раз возьмете.

Вообще, придумалось как хранить только одну строку, но там есть хитрость:

vector<int> SubsetSum(const vector<int> &numbers, int target)
  vector<int> dp(target+1);
  dp[0] = 1;
  for (auto number in numbers) {
    for (int j = target; j >= number; --j) {
      if (dp[j] == 0 && dp[j-number] != 0) dp[j] = number;
    }
  }
  vector<int> answer;
  if (dp[target] == 0) return answer;
  while (target > 0) {
    answer.push_back(target);
    cur_sum -= dp[target];
  }
  return answer;
}

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

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

Видимо, где-то ошибка. Последний цикл бесконечный, потому что условие выхода target не обновляется. А cur_sum вообще не объявлен, наверное он и есть target. Вообще, интересная идея идти от конца к началу, чтобы не хранить "новую" и "предыдущую" строку, меняя их местами, а считать, что "новая" - справа, предыдущая - слева от итератора.

Опечатка да. Должно быть target -= dp[target];

Похоже на правду. Должно быть
answer.push_back(dp[target]);

И тут тоже правда. Решил убрать лишнюю переменную, и накосячил.

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

Можно ускорить где-то в 2 раза, если заполнять только "диагональ", а числа будут отсортированы

max_filled = 0;
for (auto number in numbers) {
    for (int j = min(target, max_filled+number); j >= number; --j) {
      if (dp[j] == 0 && dp[j-number] != 0) dp[j] = number, max_filled = max(max_filled, j);

Да, действительно. Можно не перебирать слишком большие j, ведь больше суммы текущих чисел все-равно собрать мы никак не сможем.

Делала два раза подобную задачу в сверхбольшой размерности: при отгрузке листопроката (но там кроме набора заданного веса было ограничение на количество перекладываний в штабелях и дополнительно нужно было освобождать рабочее пространство склада). А вот второй раз при госзакупках нужно было выбрать выделенную на подразделение сумму с помощью списка: товары, их стоимость, требуемое количество. Суммы были огромные, но программа набирала с точностью нескольких рублей (программу писала работающая в бухгалтерии дипломница, алгоритм описан в статье, если кому интересно, дам ссылку). Идея гибридного алгоритма. Основная часть суммы набирается довольно хитрым не совсем жадным алгоритмом, т.к. производятся некоторые возвраты. Добор небольшой суммы не более чем 20 слагаемыми производится точным алгоритмом. Работает быстро и довольно точно, верхняя оценка погрешности величина последнего добавленного слагаемого, числа не целые (рубли, копейки).

Было бы замечательно, если бы вы описали свой опыт в отдельной статье. С удовольствием читаю реальные истории.

И да, как вы описываете ограничения, очень похоже на задачу о рюкзаке. Она родственна текущей.

Нет, там не рюкзак, там чистый набор заданной суммы заданными слагаемыми по недостатку (недобрать можно, перебрать нельзя). Программы у меня нет, т.к. писали под какую-то конкретную бухгалтерскую программу. Прям отдельную статью про это писать не хочется (тем более, она опубликована в "Информационных технологиях"), просто выложу сюда описание первой, приближенной части алгоритма. И на примере покажу, как он работает.

У вас там чистый рюкзак, только словесами окружён.

На скорую руку набросал модель.
import pulp as pl
import random as rnd

summa = 5000000.0    
print('Сумма к которой стремимся, руб:', summa)

n = 20 # Количество наименований
# Кортеж из максимального числа товара для заказа и цены
goods = [(rnd.randint(20, 100), rnd.randint(100000, 1000000)/100) for i in range(n)]
print('Товарная номенклатура:', goods)
print()
# Объявляем модель
model = pl.LpProblem(name="Задача о госзакупках", sense=pl.LpMaximize)
# Подключаем решатель CBC
solver = pl.PULP_CBC_CMD(msg=True)
# Объявляем переменны
x = [pl.LpVariable(name=f'x{i:03}', cat='Integer', lowBound=0, upBound=goods[i][0]) for i in range(len(goods))]
# Целевая функция
model += pl.lpSum([x[i] * goods[i][1] for i in range(n)])
# Ограничение на сумму закупки
model += pl.lpSum([x[i] * goods[i][1] for i in range(n)]) <= summa
# Выводим модель
print(model)
# Находим решение
status = model.solve(solver)
if status == 1:
    print('Фактически набрали, руб:', round(model.objective.value(), 2))
    print('Заказать шт:')
    for i, val in enumerate(model.variables()):
        print(f'Товар{i} {int(val.value())} из {goods[i][0]} шт, на сумму: {round(int(val.value()) * goods[i][1], 2)}')
Сумма к которой стремимся, руб: 5000000.0
Товарная номенклатура: [(84, 8737.57), (79, 5993.63), (70, 1959.58), (75, 4109.07), (92, 4111.34), (75, 5417.81), (74, 9230.46), (50, 4267.64), (41, 9351.66), (98, 1643.8), (82, 4038.91), (64, 9046.13), (50, 5948.99), (69, 1295.87), (80, 8750.16), (97, 4335.1), (49, 6789.67), (22, 8726.63), (35, 4979.46), (50, 1953.29)]

Задача_о_госзакупках:
MAXIMIZE
8737.57*x000 + 5993.63*x001 + 1959.58*x002 + 4109.07*x003 + 4111.34*x004 + 5417.81*x005 + 9230.46*x006 + 4267.64*x007 + 9351.66*x008 + 1643.8*x009 + 4038.91*x010 + 9046.13*x011 + 5948.99*x012 + 1295.87*x013 + 8750.16*x014 + 4335.1*x015 + 6789.67*x016 + 8726.63*x017 + 4979.46*x018 + 1953.29*x019 + 0.0
SUBJECT TO
_C1: 8737.57 x000 + 5993.63 x001 + 1959.58 x002 + 4109.07 x003 + 4111.34 x004
 + 5417.81 x005 + 9230.46 x006 + 4267.64 x007 + 9351.66 x008 + 1643.8 x009
 + 4038.91 x010 + 9046.13 x011 + 5948.99 x012 + 1295.87 x013 + 8750.16 x014
 + 4335.1 x015 + 6789.67 x016 + 8726.63 x017 + 4979.46 x018 + 1953.29 x019
 <= 5000000

VARIABLES
0 <= x000 <= 84 Integer
0 <= x001 <= 79 Integer
0 <= x002 <= 70 Integer
0 <= x003 <= 75 Integer
0 <= x004 <= 92 Integer
0 <= x005 <= 75 Integer
0 <= x006 <= 74 Integer
0 <= x007 <= 50 Integer
0 <= x008 <= 41 Integer
0 <= x009 <= 98 Integer
0 <= x010 <= 82 Integer
0 <= x011 <= 64 Integer
0 <= x012 <= 50 Integer
0 <= x013 <= 69 Integer
0 <= x014 <= 80 Integer
0 <= x015 <= 97 Integer
0 <= x016 <= 49 Integer
0 <= x017 <= 22 Integer
0 <= x018 <= 35 Integer
0 <= x019 <= 50 Integer

Фактически набрали, руб: 5000000.0
Заказать шт:
Товар0 84 из 84 шт, на сумму: 733955.88
Товар1 0 из 79 шт, на сумму: 0.0
Товар2 56 из 70 шт, на сумму: 109736.48
Товар3 60 из 75 шт, на сумму: 246544.2
Товар4 87 из 92 шт, на сумму: 357686.58
Товар5 1 из 75 шт, на сумму: 5417.81
Товар6 2 из 74 шт, на сумму: 18460.92
Товар7 50 из 50 шт, на сумму: 213382.0
Товар8 41 из 41 шт, на сумму: 383418.06
Товар9 98 из 98 шт, на сумму: 161092.4
Товар10 82 из 82 шт, на сумму: 331190.62
Товар11 63 из 64 шт, на сумму: 569906.19
Товар12 50 из 50 шт, на сумму: 297449.5
Товар13 69 из 69 шт, на сумму: 89415.03
Товар14 80 из 80 шт, на сумму: 700012.8
Товар15 97 из 97 шт, на сумму: 420504.7
Товар16 12 из 49 шт, на сумму: 81476.04
Товар17 22 из 22 шт, на сумму: 191985.86
Товар18 15 из 35 шт, на сумму: 74691.9
Товар19 7 из 50 шт, на сумму: 13673.03

Тестирование: несколько сотен номенклатур и миллионы рублей справляется за десяток секунд – минуту.

У рюкзака каждый элемент имеет вес и стоимость, здесь - только стоимость (просто я называю ее весом, так короче). Если взять целые веса и положить D=0, то получим задачу о наборе заданной суммы. Специально (для вас) взяла модельный целочисленный пример с D=0. Вы до него дочитали?

А вообще в этом нет ничего удивительного, NP-полные задачи друг к другу полиномиально сводятся. Сейчас пишу с телефона, неудобно, доберусь до десктопа, докажу полиномиальную сводимость для Рюкзака и Суммы размеров.

Стоимость это и есть вес, то, что тут нулевая масса у товара только облегчает расчёт. То, что NP-полные задачи сводятся друг к другу известный факт.

Модель выше я привёл только для того чтобы показать, что не нужно изобретать гибридные алгоритмы, когда задача сводится к каноническому рюкзаку. D не всегда будет равно 0, но модель гарантирует, что оно будет максимально малым из всех возможных вариантов. Вам только остаётся принять результат или отбросить, если D больше заданной величины.

Мой модельный пример. Нужно набрать 100; рабочий массив A = [11(10), 8 (20), 3, 2(10)], в круглых скобках указана кратность слагаемого. Результат работы алгоритма - [11(8), 8, 2(2)], набираем ровно 100.

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

Теперь второе. Умный любит учиться, как Антон Павлович Чехов однажды заметил. Покажите, как набрать нецелочисленную сумму порядка 10^7 с помощью нескольких сотен нецелочисленных слагаемых, хоть рюкзаком, хоть ДП, хоть ЦЛП, хоть чем.

Ну как же в вашем примере это и есть рюкзак, только с единичной стоимостью за каждый товар :), те каждый товар имеет одинаковый вес (стоимость).

Если 100 это не единицы, а например общий вес в килограммах, тогда ещё проще. Просто поменяем тип переменных на нецелочисленный.

Покажите, как набрать нецелочисленную сумму порядка 10^7 с помощью нескольких сотен нецелочисленных слагаемых

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

Код
import pulp as pl
import random as rnd

summa = 100000000
print('Сумма к которой стремимся:', summa)

n = 200 # Количество наименований
# Кортеж из максимального числа товара для заказа и цены
goods = [rnd.uniform(100000, 1000000) for i in range(n)]
print('Товарная номенклатура:', goods)
print()
# Объявляем модель
model = pl.LpProblem(name="test", sense=pl.LpMaximize)
# Подключаем решатель CBC
solver = pl.PULP_CBC_CMD(msg=True)
# Объявляем переменны
x = [pl.LpVariable(name=f'x{i:03}', cat='Continuous', lowBound=1000, upBound=goods[i]) for i in range(len(goods))]
# Целевая функция
model += pl.lpSum([x[i] for i in range(n)])
# Ограничение
model += pl.lpSum([x[i] for i in range(n)]) <= summa
# Выводим модель
# print(model)
# Находим решение
status = model.solve(solver)
if status == 1:
    print('Фактически набрали:', sum(val.value() for val in model.variables()))
    print('Заказать:')
    for i, val in enumerate(model.variables()):
        print(f'Товар{i} {val.value()} из {goods[i]}')
Сумма к которой стремимся: 100000000
Товарная номенклатура: [645691.2931885093, 317858.3130222785, 744057.5532077567, 133639.75382448838, 562529.939680595, 220121.72357789087, 868511.0402679854, 667797.8031388721, 388828.67949273606, 465363.4296413606, 403835.1071355331, 491265.3752803172, 862382.3549521483, 626751.8561946657, 752075.9914224785, 250814.25645500244, 694641.3500011705, 785298.1534501243, 756109.8824451261, 709645.6754570607, 917972.7722415725, 955554.7921899806, 844965.3743659591, 343735.3504507794, 460236.4328581181, 233724.30587792073, 893685.5343474532, 787940.3364260227, 719747.4481560593, 839276.7076206218, 171350.49046382695, 774521.6771186346, 769371.4382454185, 243399.19725691518, 818371.7224263975, 263151.0702352981, 676299.8690587819, 687833.2955874939, 112600.19384111586, 499453.031874757, 661994.231672027, 225596.36411810524, 958092.6936240986, 681207.4839389888, 221200.14850950055, 200157.84548828035, 286583.0234027279, 577355.6256299126, 333399.2612021096, 345060.8183317391, 240376.03111233027, 680518.4663040837, 463773.4974714662, 603426.9349440299, 585451.9677209563, 818634.3504291275, 171432.37501094883, 495861.8017951375, 103335.55560669722, 545972.0778357927, 937260.2539286873, 799451.2343092249, 417260.2715133255, 120591.66216110905, 358925.2599059303, 855285.4256332519, 347277.713113735, 372653.72638729645, 613341.4875585965, 834008.5598722058, 113592.21629697594, 194557.76160520766, 800146.2055155531, 156028.23439463944, 234383.65903390953, 621669.0491833334, 130080.20308397856, 544095.0768277377, 965060.2690712278, 199402.61458110923, 461476.2426870921, 271535.88928093534, 748673.3531951199, 291829.2366716161, 975848.7145694873, 485652.46757367696, 737538.760422223, 189907.22583928623, 109988.40031527661, 984319.965890239, 618869.2020713241, 362838.3468138749, 376860.5702819198, 322127.9759303717, 979302.096856456, 820451.3462925756, 189339.12006266386, 193907.17287269933, 748611.581533799, 655311.699386001, 796814.15879224, 423842.34665262274, 952255.412758336, 581889.2083885402, 722644.4024255745, 451917.5543063265, 785003.6822073464, 357255.87479735096, 951802.8429228318, 374079.6320391303, 424580.48065054807, 281567.7474588673, 298954.75890279224, 432199.85272635537, 581294.8771369768, 173356.4211978094, 706161.0307623808, 605862.960038322, 541192.3019834808, 329730.3633221672, 687708.8338998506, 330420.67869774863, 879502.9486966063, 808008.6649083216, 301853.2235122584, 203009.2664988107, 364144.00981451105, 627419.9685893131, 176633.14050900837, 340103.1123771363, 802280.9580293891, 444502.25690602284, 863546.4821482579, 248756.34119952025, 825107.119392471, 395423.0745651665, 708798.0964840548, 831190.0211546295, 661983.2835724461, 290013.1429606739, 262265.4727330602, 892882.5141660307, 155820.2331731233, 701627.5126059169, 805011.8389548, 687742.7342606192, 951532.8624331568, 573488.7620339128, 703435.5461830014, 825949.8310937249, 474429.0545429463, 622619.9695408016, 687086.0058438272, 934332.8887548604, 256323.9250997613, 593558.238855277, 949542.4556368465, 771376.9966740404, 786506.1936954544, 628402.9217995974, 231923.77522107676, 831762.7684432977, 414479.4428540272, 523589.8938237171, 258845.04426821976, 100903.96307715664, 163034.36959932314, 653401.590558145, 721290.1491665429, 186619.42698570632, 147015.38818919464, 141790.47086389558, 636990.7669476826, 736812.1654561599, 664545.766602985, 102627.77708906103, 786403.1857887985, 358518.36412963097, 635916.4400603629, 485805.0121598375, 468858.9812794371, 556002.1127324968, 973357.8578341331, 974833.2465184958, 242180.86936195672, 322935.57640718075, 748520.5761595789, 496686.76150578726, 325659.26421606675, 996816.1515029438, 765943.5165366514, 464941.89187614346, 784662.0608521587, 599685.5985937037, 555830.582672107, 273451.6046468464, 827005.4104140169, 772435.8386630347, 301427.2057494812, 211326.18839592487]

Фактически набрали: 99999999.99000001
Заказать:
Товар0 1000.0 из 645691.2931885093
Товар1 1000.0 из 317858.3130222785
Товар2 1000.0 из 744057.5532077567
Товар3 1000.0 из 133639.75382448838
Товар4 1000.0 из 562529.939680595
Товар5 1000.0 из 220121.72357789087
Товар6 1000.0 из 868511.0402679854
Товар7 1000.0 из 667797.8031388721
Товар8 1000.0 из 388828.67949273606
Товар9 1000.0 из 465363.4296413606
Товар10 1000.0 из 403835.1071355331
Товар11 1000.0 из 491265.3752803172
Товар12 1000.0 из 862382.3549521483
Товар13 1000.0 из 626751.8561946657
Товар14 165825.29 из 752075.9914224785
Товар15 250814.26 из 250814.25645500244
Товар16 694641.35 из 694641.3500011705
Товар17 785298.15 из 785298.1534501243
Товар18 756109.88 из 756109.8824451261
Товар19 709645.68 из 709645.6754570607
Товар20 917972.77 из 917972.7722415725
Товар21 955554.79 из 955554.7921899806
Товар22 844965.37 из 844965.3743659591
Товар23 343735.35 из 343735.3504507794
Товар24 460236.43 из 460236.4328581181
Товар25 233724.31 из 233724.30587792073
Товар26 893685.53 из 893685.5343474532
Товар27 787940.34 из 787940.3364260227
Товар28 719747.45 из 719747.4481560593
Товар29 839276.71 из 839276.7076206218
Товар30 171350.49 из 171350.49046382695
Товар31 774521.68 из 774521.6771186346
Товар32 769371.44 из 769371.4382454185
Товар33 243399.2 из 243399.19725691518
Товар34 818371.72 из 818371.7224263975
Товар35 263151.07 из 263151.0702352981
Товар36 676299.87 из 676299.8690587819
Товар37 687833.3 из 687833.2955874939
Товар38 112600.19 из 112600.19384111586
Товар39 499453.03 из 499453.031874757
Товар40 661994.23 из 661994.231672027
Товар41 225596.36 из 225596.36411810524
Товар42 958092.69 из 958092.6936240986
Товар43 681207.48 из 681207.4839389888
Товар44 221200.15 из 221200.14850950055
Товар45 200157.85 из 200157.84548828035
Товар46 286583.02 из 286583.0234027279
Товар47 577355.63 из 577355.6256299126
Товар48 333399.26 из 333399.2612021096
Товар49 345060.82 из 345060.8183317391
Товар50 240376.03 из 240376.03111233027
Товар51 680518.47 из 680518.4663040837
Товар52 463773.5 из 463773.4974714662
Товар53 603426.93 из 603426.9349440299
Товар54 585451.97 из 585451.9677209563
Товар55 818634.35 из 818634.3504291275
Товар56 171432.38 из 171432.37501094883
Товар57 495861.8 из 495861.8017951375
Товар58 103335.56 из 103335.55560669722
Товар59 545972.08 из 545972.0778357927
Товар60 937260.25 из 937260.2539286873
Товар61 799451.23 из 799451.2343092249
Товар62 417260.27 из 417260.2715133255
Товар63 120591.66 из 120591.66216110905
Товар64 358925.26 из 358925.2599059303
Товар65 855285.43 из 855285.4256332519
Товар66 347277.71 из 347277.713113735
Товар67 372653.73 из 372653.72638729645
Товар68 613341.49 из 613341.4875585965
Товар69 834008.56 из 834008.5598722058
Товар70 113592.22 из 113592.21629697594
Товар71 194557.76 из 194557.76160520766
Товар72 800146.21 из 800146.2055155531
Товар73 156028.23 из 156028.23439463944
Товар74 234383.66 из 234383.65903390953
Товар75 621669.05 из 621669.0491833334
Товар76 130080.2 из 130080.20308397856
Товар77 544095.08 из 544095.0768277377
Товар78 965060.27 из 965060.2690712278
Товар79 199402.61 из 199402.61458110923
Товар80 461476.24 из 461476.2426870921
Товар81 271535.89 из 271535.88928093534
Товар82 748673.35 из 748673.3531951199
Товар83 291829.24 из 291829.2366716161
Товар84 975848.71 из 975848.7145694873
Товар85 485652.47 из 485652.46757367696
Товар86 737538.76 из 737538.760422223
Товар87 189907.23 из 189907.22583928623
Товар88 109988.4 из 109988.40031527661
Товар89 984319.97 из 984319.965890239
Товар90 618869.2 из 618869.2020713241
Товар91 362838.35 из 362838.3468138749
Товар92 376860.57 из 376860.5702819198
Товар93 322127.98 из 322127.9759303717
Товар94 979302.1 из 979302.096856456
Товар95 820451.35 из 820451.3462925756
Товар96 189339.12 из 189339.12006266386
Товар97 193907.17 из 193907.17287269933
Товар98 748611.58 из 748611.581533799
Товар99 655311.7 из 655311.699386001
Товар100 796814.16 из 796814.15879224
Товар101 423842.35 из 423842.34665262274
Товар102 952255.41 из 952255.412758336
Товар103 581889.21 из 581889.2083885402
Товар104 722644.4 из 722644.4024255745
Товар105 451917.55 из 451917.5543063265
Товар106 785003.68 из 785003.6822073464
Товар107 357255.87 из 357255.87479735096
Товар108 951802.84 из 951802.8429228318
Товар109 374079.63 из 374079.6320391303
Товар110 424580.48 из 424580.48065054807
Товар111 281567.75 из 281567.7474588673
Товар112 298954.76 из 298954.75890279224
Товар113 432199.85 из 432199.85272635537
Товар114 581294.88 из 581294.8771369768
Товар115 173356.42 из 173356.4211978094
Товар116 706161.03 из 706161.0307623808
Товар117 605862.96 из 605862.960038322
Товар118 541192.3 из 541192.3019834808
Товар119 329730.36 из 329730.3633221672
Товар120 687708.83 из 687708.8338998506
Товар121 330420.68 из 330420.67869774863
Товар122 879502.95 из 879502.9486966063
Товар123 808008.66 из 808008.6649083216
Товар124 301853.22 из 301853.2235122584
Товар125 203009.27 из 203009.2664988107
Товар126 364144.01 из 364144.00981451105
Товар127 627419.97 из 627419.9685893131
Товар128 176633.14 из 176633.14050900837
Товар129 340103.11 из 340103.1123771363
Товар130 802280.96 из 802280.9580293891
Товар131 444502.26 из 444502.25690602284
Товар132 863546.48 из 863546.4821482579
Товар133 248756.34 из 248756.34119952025
Товар134 825107.12 из 825107.119392471
Товар135 395423.07 из 395423.0745651665
Товар136 708798.1 из 708798.0964840548
Товар137 831190.02 из 831190.0211546295
Товар138 661983.28 из 661983.2835724461
Товар139 290013.14 из 290013.1429606739
Товар140 262265.47 из 262265.4727330602
Товар141 892882.51 из 892882.5141660307
Товар142 155820.23 из 155820.2331731233
Товар143 701627.51 из 701627.5126059169
Товар144 805011.84 из 805011.8389548
Товар145 687742.73 из 687742.7342606192
Товар146 951532.86 из 951532.8624331568
Товар147 573488.76 из 573488.7620339128
Товар148 703435.55 из 703435.5461830014
Товар149 825949.83 из 825949.8310937249
Товар150 474429.05 из 474429.0545429463
Товар151 622619.97 из 622619.9695408016
Товар152 687086.01 из 687086.0058438272
Товар153 934332.89 из 934332.8887548604
Товар154 256323.93 из 256323.9250997613
Товар155 593558.24 из 593558.238855277
Товар156 949542.46 из 949542.4556368465
Товар157 771377.0 из 771376.9966740404
Товар158 786506.19 из 786506.1936954544
Товар159 628402.92 из 628402.9217995974
Товар160 231923.78 из 231923.77522107676
Товар161 831762.77 из 831762.7684432977
Товар162 414479.44 из 414479.4428540272
Товар163 523589.89 из 523589.8938237171
Товар164 258845.04 из 258845.04426821976
Товар165 100903.96 из 100903.96307715664
Товар166 163034.37 из 163034.36959932314
Товар167 653401.59 из 653401.590558145
Товар168 721290.15 из 721290.1491665429
Товар169 186619.43 из 186619.42698570632
Товар170 147015.39 из 147015.38818919464
Товар171 141790.47 из 141790.47086389558
Товар172 636990.77 из 636990.7669476826
Товар173 736812.17 из 736812.1654561599
Товар174 664545.77 из 664545.766602985
Товар175 102627.78 из 102627.77708906103
Товар176 786403.19 из 786403.1857887985
Товар177 358518.36 из 358518.36412963097
Товар178 635916.44 из 635916.4400603629
Товар179 485805.01 из 485805.0121598375
Товар180 468858.98 из 468858.9812794371
Товар181 556002.11 из 556002.1127324968
Товар182 973357.86 из 973357.8578341331
Товар183 974833.25 из 974833.2465184958
Товар184 242180.87 из 242180.86936195672
Товар185 322935.58 из 322935.57640718075
Товар186 748520.58 из 748520.5761595789
Товар187 496686.76 из 496686.76150578726
Товар188 325659.26 из 325659.26421606675
Товар189 996816.15 из 996816.1515029438
Товар190 765943.52 из 765943.5165366514
Товар191 464941.89 из 464941.89187614346
Товар192 784662.06 из 784662.0608521587
Товар193 599685.6 из 599685.5985937037
Товар194 555830.58 из 555830.582672107
Товар195 273451.6 из 273451.6046468464
Товар196 827005.41 из 827005.4104140169
Товар197 772435.84 из 772435.8386630347
Товар198 301427.21 из 301427.2057494812
Товар199 211326.19 из 211326.18839592487

Рюкзак, задача размена монетками, одномерной упаковки - это все очень похожие задачи. Например, набрать заданную сумму или как можно к ней ближе - частный случай рюкзака, где стоимость равна весу. Это, похоже, ваш случай. И ДП там практически идентичное. Правда потом возникают ньюансы, в задаче размена можно сжать память до линии с восстановлением ответа, в задаче о рюкзаке - уже нет.

Покажите, как набрать нецелочисленную сумму порядка 10^7 с помощью нескольких сотен нецелочисленных слагаемых, хоть рюкзаком, хоть ДП, хоть ЦЛП, хоть чем.

Зависит от точности. Можно все числа домножить на общий знаменатель. Но будет очень медленно, потому что сумма и слагаемые могут быль очень большими. И быстрого точного решения, по-моему, нет. Ваш метод -аппроксимация. Аппроксимационных алгоритмов можно придумать кучу.

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

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

Но все же это другой класс алгоритмов. В статье тут рассматривается точный алгоритм.

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

Я подумаю, понимаю я или еще нет))

Это приближенный алгоритм, у которого абсолютная погрешность фиксирована и не зависит от количества слагаемых и величины набираемой суммы, соответственно, асимптотическая погрешность равна 1. По-моему, для NP-hard лучше не бывает.

Если всё-таки превратить эту задачу в рюкзак. Тогда она будет звучать так. Дано мн-во A предметов, у к-х нецелый вес=стоимость, и дано нецелое число B. Вопрос. Найти рюкзак с суммарным весом <= B и суммарной стоимостью >= B. ДП эту задачу решать нельзя. Для ДП нужна целочисленность хотя бы по одному параметру.

Хорошо, для ДП есть т.н. метод масштабирования. Если попытаться применить его к этой задаче. Вначале придется умножить на 100 (рубли-копейки). Затем всё нужно уменьшить в scale раз, где scale - к-т масштабирования, и отбросить дробную часть. Величина погрешности n*scale, где n - число слагаемых. Для этой задачи scale = 10000 и n = несколько сотен.

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

Но все же это другой класс алгоритмов. В статье тут рассматривается точный алгоритм.

Спасибо вам на добром слове! Что он плох, мне уже сказал ТС, для к-го я два дня писала статью))

Да, это другой класс алгоритмов, обычная бухгалтерская программа для оптимизационной постановки задачи о наборе суммы. Наверно, было бы странно, если бы бухгалтер запустил программу, она бы прострекотала и сообщила, что решения нет. За сим откланиваюсь. Рюкзак так рюкзак.

Найти рюкзак с суммарным весом <= B и суммарной стоимостью >= B.

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

Из-за нецелых чисел - да ДП в лоб не работает, да.

Затем всё нужно уменьшить в scale раз, где scale - к-т масштабирования, и отбросить дробную часть.

Зачем? Просто считайте в копейках. Вот у вас целые числа и ДП работает.

Средняя госзакупочная сумма - 10 миллионов рублей. Переводим в копейки - миллиард. Табличка с миллиардом столбиков. Предположим, что точное решение существует. Предположим, что бухгалтер его когда-нибудь дождется. А на самом деле его устроило бы любое решение с отклонением не более, чем 100 руб.

Вот и возникает вопрос - Зачем?

А зачем вам использовать ДП для данной задачи, когда переменные не целые? Используйте симплекс метод, это и проще и быстрее. Конечно, у него есть проблема с точностью при float, но если использовать натуральные дроби в расчётах, то всё решаемо.

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

Например, отгрузка готовой продукции со склада. Стоит под погрузкой состав, 10 вагонов, за простой - огромный штраф. Как только его отгрузят, подъедет следующий состав. И при этом еще параллельно на складе идет размещение готовой продукции, теми же кранами. И там тоже нужен план. Заказы разбросаны по складу. В условиях договора прописано: заказ отгружается с точностью до 5 тонн. Никому не нужен точный набор или лучшее приближение - как только сгенерировано допустимое решение - формируется план. И тоже самое с госзакупками, правда, здесь задача однокритериальная, нужно не оптимальное решение, а первое решение с допустимым отклонением.

Если мы сейчас обсуждаем мои скромные проекты - спасибо за советы, возможно, когда-нибудь это пригодится. Пока у меня всё работает, в точном соответствии с ТЗ.

Сверхбольшая размерность это сколько? Cplex при размере матрицы более ста миллионов ненулевых значений начинает захлёбываться и не может решить систему из за постоянно накапливающихся ошибок вычислений, которые он уже не может подавить, т.к. ошибки накапливаются быстрее их подавления. При этом, насколько я понимаю, для решения СЛАУ лучше сиплекса ничего нет.

Cplex при размере матрицы более ста миллионов ненулевых значений начинает захлёбываться.

Cplex, захлёбываться. А вот CBC вытягивает. Не знаю, как они этого добиваются, сверхдлинной арифметикой или натуральными дробями или ещё как-то. Главное чтобы памяти хватало, и конечно большие матрицы это долго.

При этом, насколько я понимаю, для решения СЛАУ лучше сиплекса ничего нет.

А как же эллиптические кривые и метод внутренней точки? Гораздо быстрее симплекса, но они не подходят для целочисленного программирования.

  1. Про решение СЛАУ я вам ничего сказать не могу.

  2. И спорить в вами не буду, вы правы. Для сверхбольшой размерности делают особые алгоритмы. Пусть будет не сверхбольшая, а большая производственная размерность.

Зарегистрируйтесь на Хабре, чтобы оставить комментарий

Публикации

Истории