Многомерный метод Ньютона для уравнений, полученных с помощью Sympy

Преобразование уравнений, производных Sympy, в совместимость с scipy.optimize.newton для многомерной задачи

Я сослался на Ошибка при использовании 'exp' в sympy -TypeError и Ошибка атрибута отображается. Я обнаружил, что использование lambdify сделает символьные уравнения совместимыми с другими пакетами, такими как numpy и scipy.

Это продолжение моей предыдущей публикации: Отображается ошибка использования" exp "в sympy -TypeError и Ошибка атрибута

import numpy as np
import sympy as sym
import scipy.optimize
from sympy import symbols, diff, exp, log, power
from sympy.utilities.lambdify import lambdify

data = [3, 33, 146, 227, 342, 351, 353, 444, 556, 571, 709, 759, 836, 860, 968, 1056, 1726, 1846, 1872, 1986, 2311, 2366, 2608, 2676, 3098, 3278, 3288, 4434, 5034, 5049, 5085, 5089, 5089, 5097, 5324, 5389,5565, 5623, 6080, 6380, 6477, 6740, 7192, 7447, 7644, 7837, 7843, 7922, 8738, 10089, 10237, 10258, 10491, 10625, 10982, 11175, 11411, 11442, 11811, 12559, 12559, 12791, 13121, 13486, 14708, 15251, 15261, 15277, 15806, 16185, 16229, 16358, 17168, 17458, 17758, 18287, 18568, 18728, 19556, 20567, 21012, 21308, 23063, 24127, 25910, 26770, 27753, 28460, 28493, 29361, 30085, 32408, 35338, 36799, 37642, 37654, 37915, 39715, 40580, 42015, 42045, 42188, 42296, 42296, 45406, 46653, 47596, 48296, 49171, 49416, 50145, 52042, 52489, 52875, 53321, 53443, 54433, 55381, 56463, 56485, 56560, 57042, 62551, 62651, 62661, 63732, 64103, 64893, 71043, 74364, 75409, 76057, 81542, 82702, 84566, 88682]
n = len(data)
tn = data[n-1]


b, c = sym.symbols('b c', real=True)

f = -(-n +sum([sym.log(b*c*(num**(c-1))*sym.exp(-b*(num**c))) for num in data]))

bh = lambdify((b,c),diff(f,b),"numpy")
ch = lambdify((b,c),diff(f,c),"numpy")

sol = scipy.optimize.newton([bh,ch],(0.00404,1.0))
print(sol)

Кажется, я не могу заставить работать метод Ньютона. Любая информация или ресурсы приветствуются.


python scipy sympy scipy-optimize newtons-method
person Lucky    schedule 20.05.2019    source источник
comment
newton работает для скалярных функций. Для многомерных функций см. здесь   -  person Mstaino    schedule 21.05.2019

Ответы (2)


arrow_upward
1
arrow_downward

Согласно scipy.optimize.newton ваш x0 должен быть скаляром, а не массивом или кортежем (это то, что вы передаете scipy.optimize.newton() в своем коде). К счастью, scipy.optimize.fsolve также может находить нулевую позицию и принимать массивы как x0.

Насколько я понимаю, вы усложняете жизнь, используя sympy, когда вы можете просто использовать numpy, что ускоряет ваш код (вы можете избежать цикла for). В следующем примере показано, как его можно использовать (надеюсь, это поможет):

import warnings 
import numpy as np
import scipy.optimize as opt
# filter warnings
warnings.filterwarnings("ignore")

data = np.array([3, 33, 146, 227, 342, 351, 353, 444, 556, 571, 709, 759, 836, 860, 968,
        1056, 1726, 1846, 1872, 1986, 2311, 2366, 2608, 2676, 3098, 3278, 3288,
        4434, 5034, 5049, 5085, 5089, 5089, 5097, 5324, 5389,5565, 5623, 6080, 
        6380, 6477, 6740, 7192, 7447, 7644, 7837, 7843, 7922, 8738, 10089,
        10237, 10258, 10491, 10625, 10982, 11175, 11411, 11442, 11811, 12559, 
        12559, 12791, 13121, 13486, 14708, 15251, 15261, 15277, 15806, 16185,
        16229, 16358, 17168, 17458, 17758, 18287, 18568, 18728, 19556, 20567,
        21012, 21308, 23063, 24127, 25910, 26770, 27753, 28460, 28493, 29361,
        30085, 32408, 35338, 36799, 37642, 37654, 37915, 39715, 40580, 42015,
        42045, 42188, 42296, 42296, 45406, 46653, 47596, 48296, 49171, 49416,
        50145, 2042,  52489, 52875, 53321, 53443, 54433, 55381, 56463, 56485,
        56560, 57042, 62551, 62651, 62661, 63732, 64103, 64893, 71043, 74364,
        75409, 76057, 81542, 82702, 84566, 88682])

n    = len(data)
tn   = data[n-1]

f        = lambda b, c: n - np.sum(np.log(b * c * data**(c-1) * np.exp(-b * data**c)))
obj_func = lambda x   : np.array([f(x[0], x[1]) - x[0], f(x[0], x[1]) - x[1]])


x0  = np.array([0.00404, 1.0])
sol = opt.fsolve(obj_func, x0 = x0)
print('b, c = ', sol)

Выход:

b, c =  [0.01737803 0.4300348 ]
person SuperKogito    schedule 21.05.2019
comment
Спасибо. Для целевой функции, которую я указал в своем коде, проще использовать fsolve или root. Но мне нужно использовать метод Ньютона для очень сложной целевой функции, поэтому я пытаюсь заставить эту простую функцию работать, прежде чем переходить к сложной, - person Lucky; 21.05.2019
comment
Проблема в том, что scipy.optimize.newton не поддерживает многомерные функции, поэтому вам нужно подумать о поиске другого подхода. Это или решение для каждой переменной за раз, используя цикл, но это было бы уродливым медленным решением. Удачи. - person SuperKogito; 21.05.2019

arrow_upward
0
arrow_downward

Вместо определения символов мы можем использовать DeferredVector в sympy для определения символов. Вместо того

b, c = sym.symbols('b c', real=True)
f = -(-n +sum([sym.log(b*c*(num**(c-1))*sym.exp(-b*(num**c))) for num in data]))
bh = lambdify((b,c),diff(f,b),"numpy")
ch = lambdify((b,c),diff(f,c),"numpy"

выполните следующие действия: x = DeferredVector ('x') f = - (- n + sum ([sym.log (x [0] x [1] (num ** (c-1) ) для числа в данных]))

Теперь используйте diff w.r.t x [0] и x [1] и запустите scipy.optimize.newton ().

person Lucky    schedule 21.05.2019