#python #sympy
#python #sympy
Вопрос:
Я пытаюсь написать реализацию многомерного метода Ньютона-Рафсона, в Python
. Для начала я пытаюсь решить систему:
10 * x 3 * y * y - 3
x * x - exp(y) - 2
Хотя я намерен распространить это на любую произвольную m x n
систему (при условии, что решение возможно / существует, конечно).
Мой код выглядит следующим образом:
from dataclasses import dataclass
from sympy import *
x, y, z = symbols('x y z')
@dataclass
class Solve:
@staticmethod
def newton_raphson():
F = Matrix([10 * x 3 * y * y - 3, x * x - exp(y) - 2])
v = Matrix([x, y])
print(J := F.jacobian(v))
print(J * -1)
xx, yy = 0, 0
for i in range(10):
A = J.subs({"x": xx, "y": yy})
b = F.subs({"x": xx, "y": yy})
update = linsolve((A, b), [x, y])
(xx, yy) = tuple(*update) # should be adding, not setting equal to
print(xx, yy)
mySolver = Solve
mySolver.newton_raphson()
Моя проблема в том, что, когда я выполняю .subs()
работу с матрицами F
и J
, вычисляется все, кроме экспоненциальной функции. То есть строка print(xx, yy)
выводит:
3*(573 - 50*exp(3))/(20*(27 - 25*exp(3))) 5*(13 20*exp(3))/(4*(-27 25*exp(3)))
(-28650000*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 9) - 7500000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 9) - 53550000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 6) - 105598890*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6))) 92498679*exp(5*(351 215*exp(3))/(4*(-1350*exp(3) 729 625*exp(6)))) 164313360*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 3) 15235500*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 6) 164358750*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 3))/(-81000000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 6) - 107010000*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 6) - 31492800*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6))) 72404280*exp(5*(351 215*exp(3))/(4*(-1350*exp(3) 729 625*exp(6)))) 38032200*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 3) 87480000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 3) 9000000*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 9) 25000000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 9)) (-16200000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 6) - 40086000*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 6) - 68751693*exp(5*(351 215*exp(3))/(4*(-1350*exp(3) 729 625*exp(6)))) - 6298560*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6))) 17496000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 3) 68874390*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 3) 10450000*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 9) 5000000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 9))/(8*(-2025000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 6) - 2675250*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 6) - 787320*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6))) 1810107*exp(5*(351 215*exp(3))/(4*(-1350*exp(3) 729 625*exp(6)))) 950805*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 3) 2187000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 3) 225000*exp(1755/(4*(-1350*exp(3) 729 625*exp(6))) 1075*exp(3)/(4*(-1350*exp(3) 729 625*exp(6))) 9) 625000*exp(625*exp(6)/(-1350*exp(3) 729 625*exp(6)) 9)))
...etc
Хотя мне нужно строго числовое значение, для обновления догадок xx
и yy
.
Комментарии:
1. В Python есть синтаксис для полномочий, поэтому вместо
3*y*y
вы можете написать3*y**2
.
Ответ №1:
Sympy специализируется на символьной оценке. exp(3) является иррациональным числом и, следовательно, не может быть представлено численно. Sympy знает это.
Чтобы принудительно выполнить численное приближение, вы можете использовать N
функцию из sympy, как в
A = N(J.subs({"x": xx, "y": yy}))
Комментарии:
1. Вы также можете использовать
J.evalf(subs={"x", ...})
, который потенциально более точен в числовом отношении. Обычно что-то подобное было бы сделаноlambdify
для ускорения (например, этоnsolve
то, что делает sympy).