Решите полиномиальное уравнение с отрицательной степенью или нецелой степенью с питоном

Question

Я пытался решить уравнение

1600 = 0.41 + 6.31*d**-1.54 + 2.42*d**-3 

Mathematica может дать мне d=6.4673 менее чем за 1 секунду. Но я не смог получить ответ с символическим вычислением python.

Использование «решить» Sympy взять навсегда. Есть ли способ решить это уравнение, используя символьные вычисления python? Кажется, что проблемы возникают в основном из нецелой отрицательной силы.

Answer 1

SymPy предоставляет численные расчеты с помощью библиотеки mpmath; это включает нахождение числового корня через nsolve. В этом случае, поскольку в знаменателе есть d, мы делаем, как подсказывает строка документа nsolve, и работаем с числителем выражения и даем начальное предположение. Тот же самый корень, который уже был процитирован, быстро найден:

>>> f
-6.31*d**(-1.54) + 1599.59 - 2.42/d**3
>>> nsolve(f.as_numer_denom()[0], 1)
0.119100503944930

Answer 2

Первый вопрос, который нужно задать: хотите ли вы символьное или числовое решение?

Для символического решения: его нет. После подстановки x = d**(-1/50) уравнение становится A*x**150 + B*x**77 + C == 0. Не существует символической формулы для решения таких полиномиальных уравнений высокой степени.

Для числового решения: вам не нужен SymPy, потому что SymPy предназначен для символьных вычислений. Найти рут с помощью SciPy. В качестве отправной точки:

from scipy.optimize import root
root(lambda d: 0.41 + 6.31*d**(-1.54) + 2.42*d**(-3) - 1600, 0.1)

Это дает 0.1191005 как решение. Начальная точка должна быть небольшим положительным числом, иначе решатель не сможет сходиться. Как сказал WIP, Mathematica потерпела неудачу таким образом, и ее ответ является поддельным.

Но для скалярных уравнений лучше использовать специализированный решатель, такой как brentq, тем более что у вас здесь монотонная функция. Для этого решателя для начала требуется интервал брекетинга: одна точка, где функция положительная, и другая, где она отрицательная. Без калькулятора можно заметить, что 0,1 дает положительное значение (один из терминов 2.42*1000), а 1 дает отрицательное (три небольших числа минус 1600). Итак,

from scipy.optimize import brentq
brentq(lambda d: 0.41 + 6.31*d**(-1.54) + 2.42*d**(-3) - 1600, 0.1, 1)

, который быстро и надежно возвращается с 0.11910050394499523.