Корректная последовательность команд для уравнения symboli c с использованием sympy

Question

Примечание: я новичок в симпати и пытаюсь понять, как это работает.

Что у меня сейчас: Я получаю правильные решения, но мне нужно 35 - 50 секунд.

Цель: Чтобы ускорить вычисления, определив уравнение symboli c один раз, а затем повторно используя его с различными переменными.

Настройка: Мне нужно вычислить полином G (t) (t = 6 корней) для каждой итерации l oop. (Всего 220 итераций) G (t) имеет 6 других переменных, которые вычисляются и известны на каждой итерации. Эти переменные отличаются на каждой итерации.

Первая попытка (медленно): Я просто помещаю каждую в одну python функцию, где я определил Gt символически и решил для т . Он работал около 35-40 секунд. function_G вызывается на каждой итерации.

def function_G(f1, f2, a, b, c, d):
    t = sp.symbols('t')
    left = t * ((a * t + b)**2 + f2**2 * (c*t+d)**2)**2 
    right = (a*d-b*c) * (1+ f1**2 * t**2)**2 * (a*t+b) * (c*t+d)
    eq = sp.expand(left - right)
    roots = sp.solveset(Gt, t)
    return roots

• Тогда человек дал мне подсказку, что:

Вам нужно только (символически ) определите коэффициенты многочлена один раз, как шаг предварительной обработки. После этого, при обработке каждой итерации, вы просто вычисляете коэффициенты полинома, а затем решаете для корней.

• Я попросил уточнить, кто добавил:

Итак, я определил функцию g (t), а затем использовал sympy.expand для определения всех скобок / показателей степени, а затем sympy.collect для сбора терминов по степеням t. Наконец, я использовал .coeff на выходе команды collect, чтобы получить коэффициенты для подачи в numpy. root.

Вторая попытка: Чтобы следовать совету, я определил G (t) сначала символически и передал его функции, которая запускает l oop вместе с его параметрами symboli c. Таким образом, функция constructGt () вызывается только один раз.

def constructGt():
    t, a, b, c, d, f1, f2 = sp.symbols('t a b c d f1 f2')
    left = t * ((a * t + b)**2 + f2**2 * (c*t+d)**2)**2 
    right = (a*d-b*c) * (1+ f1**2 * t**2)**2 * (a*t+b) * (c*t+d)
    gt = sp.Eq(left - right, 0)
    expanded = sp.expand(gt)
    expanded = sp.collect(expanded, t)

    g_vars = {
      "a": a,
      "b": b,
      "c": c,
      "d": d,
      "f1": f1,
      "f2": f2
    }

    return expanded, g_vars

, затем на каждой итерации я передавал функцию и ее параметры для получения корней:

#Variables values:
#a = 0.00011713490404073987 
#b = 0.00020253296124588926 
#c = 4.235688216068313e-07 
#d = 0.012262546040805029 
#f1= -0.012553203944721956 
#f2 = 0.018529776776949003

def function_G(f1_, f2_, a_, b_, c_, d_, Gt, v):
    Gt = Gt.subs([(v['a'], a_), (v['b'], b_), 
                  (v['c'], c_), (v['d'], d_),
                  (v['f1'], f1_), (v['f2'], f2_)])

    roots = sp.solveset(Gt, t)
    return roots 

Но она становилась еще медленнее около 56 секунд.

Вопрос: Я не понимаю, что я делаю не так? Я также не понимаю, как этот человек использует .coeff (), а затем np.roots на результаты.

Answer 1

Даже если ваши f1 и f2 являются линейными по переменной, вы работаете с полтином четверти c, а корни его очень длинные. Если это одномерное, то лучше просто использовать выражение, решить его при некотором значении констант, где решение известно, а затем использовать это значение и новые константы, которые относительно близки к старым, и использовать nsolve , чтобы получить следующий root. Если вы заинтересованы в более чем одном решении, вам, возможно, придется «следовать» за каждым root отдельно с nsolve ... но я думаю, что вы будете намного счастливее с общей производительностью. Использование real_roots - это еще один вариант, особенно если выражение является просто полиномом от некоторой переменной.

Учитывая, что вы работаете с квартирой c, вы должны помнить об этом: общее решение так долго и сложный (за исключением очень особых случаев), что неэффективно работать с общим решением и подставлять в значения, как они известны. Однако это очень легко найти для числовых значений и намного быстрее:

Сначала создайте выражение symboli c, в которое будут подставляться значения; Безошибочные символы «ванили» используются:

t, a, b, c, d, f1, f2 = symbols('t a b c d f1 f2')
left = t * ((a * t + b)**2 + f2**2 * (c*t+d)**2)**2 
right = (a*d-b*c) * (1+ f1**2 * t**2)**2 * (a*t+b) * (c*t+d)
eq = left - right

Затем определите словарь замен для замены в выражении, отметив, что dict(x=1) создает {'x': 1} и когда это используется с subs Символ ванили будет создан для "x":

reps = dict(
a = 0.00011713490404073987 ,
b = 0.00020253296124588926 ,
c = 4.235688216068313e-07 ,
d = 0.012262546040805029 ,
f1= -0.012553203944721956 ,
f2 = 0.018529776776949003)

Оцените действительные корни выражения:

from time import time
t=time();[i.n(3) for i in real_roots(eq.subs(reps))];'%s sec' % round(time()-t)
[-11.5, -1.73, 8.86, 1.06e+8]
'3 sec'

Найдите все 6 корней выражения, но возьмите только действительные части:

>>> roots(eq.subs(reps))
{-11.4594523988215: 1, -1.73129179415963: 1, 8.85927293271708: 1, 106354884.4365
42: 1, -1.29328524826433 - 10.3034942999005*I: 1, -1.29328524826433 + 10.3034942
999005*I: 1}
>>> [re(i).n(3) for i in _]
[-11.5, -1.73, 8.86, 1.06e+8, -1.29, -1.29]

Измените одно или несколько значений и сделайте это снова

reps.update(dict(a=2))
[i.n(3) for i in real_roots(eq.subs(reps))]
[-0.0784, -0.000101, 0.0782, 3.10e+16]

Обновите значения в al oop:

>>> a = 1
>>> for i in range(3):
...     a += 1
...     reps.update(dict(a=a))
...     a, real_roots(eq.subs(reps))[0].n(3)
...
(2, -0.0784)
(3, -0.0640)
(4, -0.0554)

Примечание: при использовании roots, реальные корни будут появляться сначала в отсортированном порядке, а затем воображаемые корни будут приходить в сопряженных парах (но не в любом порядке).