Система полиномиальных уравнений второй степени

[DMZ]

Требуется программно решать системы полиномиальных уравнений специального вида:

$$x_i \sum\limits_{j \in T_i} x_j = a_i \; i \in [1 ; N]$$

при следующих условиях:

$$ a_i > 0 \; \forall i \qquad (1) $$

$$ i \in T_j \to j \in T_i \; \forall i,j \qquad (2) $$

$$ i \in T_i \; \forall i \qquad (3) $$

Размер списков $T_i$ одинаков и примерно равен логарифму количества неизвестных $N$.
$N$ достаточно велико (от 10 до 1000, а лучше ещё больше). Требуется очень точное и быстро вычисляемое решение. Думаю, что если сгенерировать код, аналитически вычисляющий решение, то этого будет достаточно. Численно решать не хочу - к приемлемой точности сильно долго будет сходиться
Я попробовал скормить такие уравнения SymPy (библиотека символьных вычислений на Питоне) - повисла на 10 неизвестных. С другими символьниками я на Вы, а мне же кодогенерация потом понадобится.
В Википедии https://en.wikipedia.org/wiki/System_of ... _equations сильно общо написано и без конкретных рецептов - говорят к треугольному виду приводите. А как? Подстановкой 1000 переменных? К тому же у меня сильно специальный вид, может как-то проще можно сделать?

На маленьких системах видно, что всегда получается два решения - все $x$ положительные или все отрицательные. Мне оба варианта подходят, да и по уравнениям видно, что перемена знаков ничего не меняет.

-- 02.03.2016, 21:54 --

Уточнение относительно $T_i$.

Количество элементов в $T_i$ одинаково для любых $i$. Обозначим его $M$. При этом каждое $x_j$ встречается ровно в $M$ различных $T_i$.
Пример:

$N = 4$, $M = 3$

$\left\{ \begin{array}{rcl} x_1(x_4 + x_1 + x_2) = a_1 \\ x_2(x_1 + x_2 + x_3) = a_2 \\ x_3(x_2 + x_3 + x_4) = a_3 \\ x_4(x_3 + x_4 + x_1) = a_4 \\ \end{array} \right.$

[Sonic86]

$T_i$ фиксированы или варьируются?

[Pphantom]

А чем обычный метод Ньютона плох? Матрица Якоби тут считается тривиально, все должно быть очень просто.

[DMZ]

Sonic86 в сообщении #1103792 писал(а):

$T_i$ фиксированы или варьируются?

Считаются один раз - при запуске программы. Здесь можно долгую подготовку сделать. Потом будут меняться только $a_i$.

[DeBill]

DMZ
В Вашем примере, все $T_i$ получаются друг из друга сдвигом. Это всегда так, или - случайность?

Можно из Ваших уравнений выразить $x_i = \frac{a_i}{(...)}$, так что ур-е примет вид $x=F(x)$ с довольно приличной $F$. Может, пойдет метод итераций?
Или: (если ТАК - верно) решать это последнее ур-е методом Ньютона: вроде бы , тогда в матрице $F'(x)$ строчки будут почти одинаковыми, так что есть надежда обратить ее - что-нить типа прогонки устроить....

[DMZ]

Pphantom в сообщении #1103805 писал(а):

А чем обычный метод Ньютона плох? Матрица Якоби тут считается тривиально, все должно быть очень просто.

Метод Ньютона же численный, то есть неизвестно сколько сходиться будет до нужной точности.
Или он для квадратных функций точен и сходится за гарантированное число итераций?

[Pphantom]

DMZ в сообщении #1103931 писал(а):

Метод Ньютона же численный, то есть неизвестно сколько сходиться будет до нужной точности.

В общем случае - да, но у Вас тут совсем не общий случай.

[DMZ]

DeBill в сообщении #1103860 писал(а):

DMZ
В Вашем примере, все $T_i$ получаются друг из друга сдвигом. Это всегда так, или - случайность?

Не совсем. Там в принципе может быть произвольный список, но можно ограничится такими $T_i$ который состоят из индексов, отличающиеся от $i$ только одним битом в двоичном представление - отсюда длина списков $T_i$ равна логарифму $N$ (то есть битности индексов).

Цитата:

Можно из Ваших уравнений выразить $x_i = \frac{a_i}{(...)}$, так что ур-е примет вид $x=F(x)$ с довольно приличной $F$. Может, пойдет метод итераций?
Или: (если ТАК - верно) решать это последнее ур-е методом Ньютона: вроде бы , тогда в матрице $F'(x)$ строчки будут почти одинаковыми, так что есть надежда обратить ее - что-нить типа прогонки устроить....

Это численный метод. А требуется большая точность и скорость. Плюс в ряде важных случаев будет сильный разброс по масштабу значений $x$ - кажется это называется плохой обусловленностью (пожалуйста, поправьте, если я не правильно термин применяю).
Я скорее возьмусь за реализацию какого-нибудь rational univariate representation. Хотя терминология в соотвествующей литературе меня в депрессию вгоняет.

[DeBill]

DMZ
1.Во первых, я хотел бы развеять Ваше предубеждение касательно численных методов. Реально, точно комп вычисляет лишь суммы-произведения, все остальное он делает приближенно - но быстро, и с высочайшей точностью. Просто он работает по быстрым и хорошим алгоритмам - и в результате у нас создается иллюзия о "точных" решениях. Даже когда мы нашли решение типа $x = 1+ \sqrt{5}$ , или $x =\sin(\frac{\pi}{7})$ - вычисляется ответ по численным алгоритмам. Далее

DMZ в сообщении #1103931 писал(а):

Или он для квадратных функций точен и сходится за гарантированное число итераций?

Нет, конечно. Но если Вы скажете, какая точность Вас устраивает, можно что-то и гарантировать...
2. Для Вашей задачи: по моим прикидкам, точное решение (я под этим имею в виду - решение в радикалах) у Вас не получится (даже для Вашей простейшей системы из черырех ур-й, получается ур-е высокой степени).
3. Метод Ньютона, похоже, для Вашей системы подходит просто идеально. Только надо ее немножко препарировать. Именно, перепишем ее в виде

$x_i +\sum\limits_{k\approx i}^{} x_k = \frac{a_i}{x_i}$, где сумма берется по всем "соседям" $i$ (т.е., по индексам, отличающимся от $i$ ровно одним битом)

Решение будем искать последовательными приближениями $x_i^{n+1} = x_i^n + h_i$; в качестве начального мне нравится $x_i^0 = \sqrt{a_i}$. В соответствии с методом Ньютона, получим для $h_i$ систему

$h_i \cdot (1 + \frac{a_i}{(x_i^n)^2})= - \sum\limits_{k\approx i}^{} h_k -\Delta _i$

где $\Delta _i = x_i ^n+\sum\limits_{k\approx i}^{} x_k^n - \frac{a_i}{x_i ^n}$ - невязка (ошибка) на предыдущем шаге. Эту систему я тоже рекомендую решать приближенно: подставлять в правую часть текущее приближение (взяв начальное - нулевым), и - поделив, находить следующее.
4. Прикидки: для $M=10, a_i \approx 100$, точность $10^{-8}$ будет достигаться: за 3-4 итерации для поправок $h_i$, и за 5-10 операций - для $x_i$.
Попробуйте!

-- 04.03.2016, 21:39 --

DMZ
Критерий истины - практика. Так что попробовать стоит. И - если все будет хорошо - метод чрезвычайно легко может быть адаптирован и на случай других "окружений" $T_i$ - просто во всех суммах надо будет поменять области суммирования - и все!

[Евгений Машеров]

По-моему, тут ситуация, когда итеративный метод быстрее "прямого" и приближённый точнее "точного".
Кстати, даже есть получится решение через треугольную систему - приведение к ней будет численным.