Найдите правильный корень параметризованной функции, данное решение для одного набора параметровPython

Программы на Python
Ответить Пред. темаСлед. тема
Anonymous
 Найдите правильный корень параметризованной функции, данное решение для одного набора параметров

Сообщение Anonymous »

Допустим, у меня есть функция Foo (x, a, b) , и я хочу найти конкретный из его (потенциально множественных) корней x0 , то есть значение x0 , что foo (x0, a, b) == 0 . Я знаю, что для (a, b) == (0, 0) root i want is x0 == 0 и что функция непрерывно изменяется с и b , поэтому я могу «следить за« корнем из (0, 0) к желаемому (a, b) . class = "lang-py prettyprint-override">

Код: Выделить всё

def foo(x, a, b):
return (1 + a) * np.sin(a + b - x) - x
Изображение

For (a, b) == (0, 0) I want to the root at 0, for (2, 0) I want the one near 1.5 and Для (2, 1) я хочу, чтобы один около 2.2 .
Теперь эта проблема кажется такой, которая может быть достаточно общей, чтобы иметь подготовленный, быстрый решатель в Scipy или другой стандартный пакет (или инструменты для легко и эффективно построить один). Тем не менее, я не знаю, какие термины искать, чтобы найти его (или убедиться, что его не существует). есть готовый инструмент для этого? Как называется такая проблема? Когда это происходит, возвращение NAN является предпочтительным поведением, хотя это не очень важно. Я продолжу расчет корень для множества различных параметров, например, для построения и интеграции над ними. Конечно, это не очень быстро, что в некоторой степени ограничивает меня в моем реальном приложении. < /P>

Код: Выделить всё

import functools
import numpy as np
from scipy.optimize import root_scalar
from matplotlib import pyplot as plt

def foo(x, a, b):
return (1 + a) * np.sin(a + b - x) - x

fig, ax = plt.subplots()
ax.grid()
ax.set_xlabel("x")
ax.set_ylabel("foo")
x = np.linspace(-np.pi, np.pi, 201)
ax.plot(x, foo(x, 0, 0), label="(a, b) = (0, 0)")
ax.plot(x, foo(x, 2, 0), label="(a, b) = (2, 0)")
ax.plot(x, foo(x, 2, 1), label="(a, b) = (2, 1)")
ax.legend()
plt.show()

# Semi-bodged solver for reference:

def _dfoo(x, a, b):
return -(1 + a) * np.cos(a + b - x) - 1

def _solve_fooroot(guess, a, b):
if np.isnan(guess):
return np.nan
# Determine limits for finding the root.
# Allow for slightly larger limits to account for numerical imprecision.
maxlim = 1.1 * (1 + a)
y0 = foo(guess, a, b)
if y0 == 0:
return guess
dy0 = _dfoo(guess, a, b)
estimate = -y0 / dy0
# Too small estimates are no good.
if np.abs(estimate) < 1e-2 * maxlim:
estimate = np.sign(estimate) * 1e-2 * maxlim
for lim in np.arange(guess, guess + 10 * estimate, 1e-1 * estimate):
if np.sign(foo(lim, a, b)) != np.sign(y0):
bracket = sorted([guess, lim])
break
else:
return np.nan
sol = root_scalar(foo, (a, b), bracket=bracket)
return sol.root

@functools.cache
def _fooroot(an, astep, bn, bstep):
if an == 0:
if bn == 0:
return 0
guessan, guessbn = an, bn - 1
else:
guessan, guessbn = an - 1, bn
# Avoid reaching maximum recursion depth.
for thisbn in range(0, guessbn, 100):
_fooroot(0, astep, thisbn, bstep)
for thisan in range(0, guessan, 100):
_fooroot(thisan, astep, guessbn, bstep)
guess = _fooroot(guessan, astep, guessbn, bstep)
return _solve_fooroot(guess, an * astep, bn * bstep)

@np.vectorize
def fooroot(a, b):
astep = (-1 if a < 0 else 1) * 1e-2
bstep = (-1 if b < 0 else 1) * 1e-2
guess = _fooroot(int(a / astep), astep, int(b / bstep), bstep)
return _solve_fooroot(guess, a, b)

print(fooroot(0, 0))
print(fooroot(2, 0))
print(fooroot(2, 1))

fig, ax = plt.subplots()
ax.grid()
ax.set_xlabel("b")
ax.set_ylabel("fooroot(a, b)")
b = np.linspace(-3, 3, 201)
for a in [0, 0.2, 0.5]:
ax.plot(b, fooroot(a, b), label=f"a = {a}")
ax.legend()
plt.show()


Подробнее здесь: https://stackoverflow.com/questions/795 ... of-paramet
Реклама
Ответить Пред. темаСлед. тема

Быстрый ответ

Изменение регистра текста: 
Смайлики
:) :( :oops: :roll: :wink: :muza: :clever: :sorry: :angel: :read: *x)
Ещё смайлики…
   
К этому ответу прикреплено по крайней мере одно вложение.

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

Максимально разрешённый размер вложения: 15 МБ.

  • Похожие темы
    Ответы
    Просмотры
    Последнее сообщение

Вернуться в «Python»