Частично факторинг выражения в Sympy

Question

Предположим, у меня есть выражение вида . Я знаю, что могу упростить выражение следующим образом: . Однако sympy.simplify и sympy.factor оба возвращают исходное выражение.

Чтобы обойти это, я оперировал выражением на низком уровне:

factor_map = defaultdict(set)
additive_terms = expr.as_coeff_add()[-1]
for term1, term2 in combinations(additive_terms, 2):
    common_terms = (
        set(term1.as_coeff_mul()[-1])
        & set(term2.as_coeff_mul()[-1])
    )
    if common_terms:
        common_factor = sympy.Mul(*common_terms)
        factor_map[common_factor] |= {term1, term2}

factor_map теперь выглядит так:

{
    a: {a⋅x, -a⋅y}, 
    b: {b⋅x, -b⋅y}, 
    c: {-c⋅x, c⋅y}, 
    x: {a⋅x, b⋅x, -c⋅x}, 
    y: {-a⋅y, -b⋅y, c⋅y}
}

Я сортирую это по количеству операций, представленных терминами:

factor_list = sorted(
    factor_map.items(),
    key = lambda i: (i[0].count_ops() + 1) * len(i[1])
)[::-1]

Затем я просто перестроил выражение:

used = set()
new_expr = 0
for item in factor_list:
    factor = item[0]
    appearances = item[-1]
    terms = 0
    for instance in appearances:
        if instance not in used:
            terms += instance.as_coefficient(factor)
            used.add(instance)
    new_expr += factor * terms
for term in set(additive_terms) - used:
    new_expr += term

Это дает new_expr = d + x*(a + b - c) + y*(-a - b + c). Не отлично, но лучше.

Я также могу улучшить это, разделив каждую комбинацию аддитивных терминов друг на друга, проверив, является ли результат числом, и используя эту информацию для дальнейшего уменьшения вывода до new_expr = d + (x - y)*(a + b - c).

Я также пытался применить sympy.factor ко всем возможным комбинациям аддитивных терминов, но, очевидно, это очень быстро взрывается для любого достаточно большого выражения.

---

Редактировать: Вот реализация, которая использует sympy.factor на всех разделах набора дополнительных терминов (функция разделов заимствована у этого ответа ):

def partition(collection):
    if len(collection) == 1:
        yield [ collection ]
        return

    first = collection[0]
    for smaller in partition(collection[1:]):
        # insert `first` in each of the subpartition's subsets
        for n, subset in enumerate(smaller):
            yield smaller[:n] + [[ first ] + subset]  + smaller[n+1:]
        # put `first` in its own subset 
        yield [ [ first ] ] + smaller


def partial_factor(expr):
    args = list(expr.as_coeff_add()[-1])
    # Groupings is just a cache to avoid duplicating factor operations
    groupings = {}

    unique = set()

    for p in partition(args):
        new_expr = 0
        for group in p:
            add_group = sympy.Add(*group)
            new_expr += groupings.setdefault(add_group, sympy.factor(add_group))
        unique.add(new_expr)
    return sorted(unique, key=sympy.count_ops)[0]

Для такого выражения, как a*x + b*y + c*z + d + e*x + f*y + h*z, на моем компьютере требуется 7,8 секунды, тогда как другой метод занимает 378 микросекунд и дает тот же результат. Похоже, должен быть способ быть более строгим, чем первый метод, не требуя в 20 000 раз больше времени для его решения.

---

Я чувствую, что не должно быть так сложно получить то, что я хочу. Есть ли более простой способ сделать это?

Answer 1

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

Дано выражение: попытаться его учесть. Если неудачно, посмотрите на коэффициент каждого символа, который он содержит; метод Expr.coeff(Symbol) делает это. Символ с наименьшим коэффициентом (измеряемый количеством содержащихся символов) считается препятствием для факторинга и удаляется из выражения. Повторение.

Эта логика закодирована ниже, а partial_factor(a*x + b*x - c*x - a*y - b*y + c*y + d) возвращает d + (x - y)*(a + b - c).

def partial_factor(expr):
    to_factor = expr
    non_factorable = 0
    while to_factor != 0:
        if factor(to_factor) != to_factor:
            return factor(to_factor) + non_factorable
        coeffs = {v: to_factor.coeff(v) for v in to_factor.free_symbols}
        min_size = min([len(coeffs[v].free_symbols) for v in coeffs])
        for v in coeffs:
            if len(coeffs[v].free_symbols) == min_size:
                non_factorable += v*coeffs[v]
                to_factor -= v*coeffs[v]
                break
    return expr