полином n-ой степени, численный поиск корней

[wizardofnothing]

Уважаемые знатоки. У меня есть задача, решить которую самостоятельно, не хватает опыта. Поэтому искрене прошу советов

Итак. У меня есть некоторая компьютерная программа, в которой производятся некоторые статистические вычисления. На одном из шагов мне нужно найти все корни полинома n-степени (причем может быть даже что-то вроде такого k*x^50+x^49+...=0).
Собственно, нужно выбрать наиболее отпимальный (в реализации\быстродейтсвии) способ нахождения корней.

1. Какой способ можете посоветовать?
2. Может быть кто-нибудь владеет \ знает где скачать программный код (на c++ или delphi) для решения поставленной задачи?

[Klad33]

Чем Вам такое не нравится?
http://www.wolframalpha.com/input/?i=2x ... 9%2Bx%2B10

[Joker_vD]

Лично я бы посоветовал последовательно использовать модифицированный метод Ньютона (который не давится кратными корнями).

[VAL]

Klad33 в сообщении #496055 писал(а):

Чем Вам такое не нравится?
http://www.wolframalpha.com/input/?i=2x ... 9%2Bx%2B10

Мне не нравится тем, что корней маловато. Их, все же 50. И кратных нет:

> f:=2*x^50-5*x^49+3*x^48+x+10:df:=diff(f,x);

                           49        48        47
                df := 100 x   - 245 x   + 144 x   + 1

> gcd(f,df);
> evalf(solve(2*x^50-5*x^49+3*x^48+x+10),20);

  1.1070783976502150137, 1.4999999594194430383,

        1.0887526302702221950 + 0.12832512894031300781 I,

        1.0477544511244468966 + 0.26112099742817438460 I,

        0.99286177117267352192 + 0.38964662292861444653 I,

        0.92491753395247690926 + 0.50999630544203193082 I,

        0.84440846523346705649 + 0.62038090459553681292 I,

        0.75218958623981078991 + 0.71938932194096141145 I,

        0.64944355849241159127 + 0.80574170757065662541 I,

        0.53761167686690813376 + 0.87829403300288033756 I,

        0.41833722779476571614 + 0.93607274113308065674 I,

        0.29341812777692649500 + 0.97830543984969870585 I,

        0.16476477545918509721 + 1.0044436424052003874 I,
        
        0.034360770458306630867 + 1.0141781219314906225 I,

        -0.095774606469229708680 + 1.0074478071694568603 I,

        -0.22362252981448731435 + 0.98444293380485576513 I,

        -0.34719935415217188307 + 0.94560298425854776394 I,

        -0.46458968217603803295 + 0.89160986053610185589 I,

        -0.57397742798925992156 + 0.82337671514230618447 I,

        -0.67367449928204851731 + 0.74203288444702871222 I,

        -0.76214674521744164696 + 0.64890540378621331089 I,

        -0.83803687662356082262 + 0.54549761604446718058 I,

        -0.90018414845240605465 + 0.43346540192908269239 I,

        -0.94764068764374588853 + 0.31459155165145178832 I,

        -0.97968443752576770068 + 0.19075876113745373988 I,

        -0.99582875803027256802 + 0.063921675067892839985 I,

        -0.99582875803027256802 - 0.063921675067892839985 I,

        -0.97968443752576770068 - 0.19075876113745373988 I,

        -0.94764068764374588853 - 0.31459155165145178832 I,

        -0.90018414845240605465 - 0.43346540192908269239 I,

        -0.83803687662356082262 - 0.54549761604446718058 I,

        -0.76214674521744164696 - 0.64890540378621331089 I,

        -0.67367449928204851731 - 0.74203288444702871222 I,

        -0.57397742798925992156 - 0.82337671514230618447 I,

        -0.46458968217603803295 - 0.89160986053610185589 I,

        -0.34719935415217188307 - 0.94560298425854776394 I,

        -0.22362252981448731435 - 0.98444293380485576513 I,

        -0.095774606469229708680 - 1.0074478071694568603 I,

        0.034360770458306630867 - 1.0141781219314906225 I,

        0.16476477545918509721 - 1.0044436424052003874 I,

        0.29341812777692649500 - 0.97830543984969870585 I,

        0.41833722779476571614 - 0.93607274113308065674 I,

        0.53761167686690813376 - 0.87829403300288033756 I,

        0.64944355849241159127 - 0.80574170757065662541 I,

        0.75218958623981078991 - 0.71938932194096141145 I,

        0.84440846523346705649 - 0.62038090459553681292 I,

        0.92491753395247690926 - 0.50999630544203193082 I,

        0.99286177117267352192 - 0.38964662292861444653 I,

        1.0477544511244468966 - 0.26112099742817438460 I,

        1.0887526302702221950 - 0.12832512894031300781 I

[ewert]

wizardofnothing в сообщении #496052 писал(а):

На одном из шагов мне нужно найти все корни полинома n-степени (причем может быть даже что-то вроде такого k*x^50+x^49+...=0).
Собственно, нужно выбрать наиболее отпимальный (в реализации\быстродейтсвии) способ нахождения корней.

Ничего -- на таких степенях любой алгоритм будет численно крайне неустойчив.

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

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

Присоединяюсь к совету. Если в задаче что-то требует нахождения корней полинома 50-й степени, то что-то в постановке ложно.

(Оффтоп)

Был у меня один знакомый - пол-аспирантуры убил на разработку метода нахождения корней характеристического полинома корреляционной матрицы, с построением как-то справился, методом Фадеева-Леверье, хоть и медленно, а корни оказались крайне неустойчивы, библиотеку повышенной точности сам написал, диссертация горела.
Даже не знаю, сделал ли я доброе дело, рассказав о стандартной программе нахождения С.З. методом Якоби. С одной стороны, диссер он в срок вполне успешно сдал, с другой - обидел я его сильно, жалко смотреть было на последовавшую истерику...

[svv]

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

В книге С.Л.Марпл, "Цифровой спектральный анализ и его приложения" для этой цели предлагается использовать алгоритм Дженкинса-Трауба, там приводится и программа (правда, на Фортране). Я её когда-то давно переводил на C++. Программа довольно сложная, занимает несколько страниц, но хорошие алгоритмы на дороге не валяются. Находит все корни полинома -- вещественные, комплексные. Коэффициенты полинома тоже в общем случае комплексные.

[Shadow]

Евгений Машеров в сообщении #496188 писал(а):

Присоединяюсь к совету. Если в задаче что-то требует нахождения корней полинома 50-й степени, то что-то в постановке ложно.

Игра: Генерирутся последовательно 100 случайных чисел в интервале $(0;1)$. После каждого наибольшего числа можем или выбрать, или подождать. Цель-выбрать наибольшее число. При каких "критических" значениях будем выбирать.

[VAL]

Shadow в сообщении #496363 писал(а):

Евгений Машеров в сообщении #496188 писал(а):

Присоединяюсь к совету. Если в задаче что-то требует нахождения корней полинома 50-й степени, то что-то в постановке ложно.

Игра: Генерирутся последовательно 100 случайных чисел в интервале $(0;1)$. После каждого наибольшего числа можем или выбрать, или подождать. Цель-выбрать наибольшее число. При каких "критических" значениях будем выбирать.

Известная стратегия - пропустить примерно $\frac{100}e$ чисел, запомнив наибольшее из них, и выбрать первое большее того, что запомнили.
Только причем здесь исходная задача?

[Shadow]

Цитата:

Известная стратегия - пропустить примерно $\frac{100}e$ чисел, запомнив наибольшее из них, и выбрать первое большее того, что запомнили.
Только причем здесь исходная задача?

Это когда нет ограничений. А когда числа ограничены в интервале $(0;1)$ задача совершенно другая.
Т.е понимаете, что если с первого раза попадется 0.99999999991, ждать $\frac{100}{e}$ попыток не буду.

[Cash]

Так ведь она и легче гораздо.
Если появился на $k$-м числе максимум $a$, то вероятность появления еще хотя бы одного максимума $1-a^{100-k}$. Ну и сравниваем с $1/2$. Каким образом полиномы выскакивают?

[VAL]

Shadow в сообщении #496384 писал(а):

Это когда нет ограничений. А когда числа ограничены в интервале $(0;1)$ задача совершенно другая.
Т.е понимаете, что если с первого раза попадется 0.99999999991, ждать $\frac{100}{e}$ попыток не буду.

Угу. Невнимательно прочел.
Сколько раз зарекался по ночам не отвечать...

[Shadow]

Cash в сообщении #496389 писал(а):

Так ведь она и легче гораздо.
Если появился на $k$-м числе максимум $a$, то вероятность появления еще хотя бы одного максимума $1-a^{100-k}$. Ну и сравниваем с $1/2$. Каким образом полиномы выскакивают?

Результат у Вас получится не плохой, но не оптимальный. А зачем сравнивать с $\frac 1 2$. С другим надо сравнивать.....и все не так просто

[Cash]

Чувствую себя злостным оффтопщиком...
Можно попросить модераторов выделить задачу в отдельную тему?

Цитата:

Результат у Вас получится не плохой, но не оптимальный

В задаче надо ставить четкие критерии, а не расплывчатые цели.
Я исходил из следующего:если игрок угадал максимум, то он получает рубль, если нет - рубль отдает. Цель - максимизировать прибыль.
У вас есть сомнения, что предложенная мной стратегия не оптимальна?
Любая другая стратегия будет ошибочна.
Допустим появляется на $k$-м шаге максимум $a$ и вы принимаете решение брать или не брать.
Пусть $p=1-a^{100-k}$. Если $p>1/2$ и вы отказываетесь или $p<1/2$ и вы берете, то вы играете в минус, при следовании стратегии - в плюс. Какие еще могут быть вопросы?
А вот величина этой прибыли - это отдельная история...

[Shadow]

Cash в сообщении #496406 писал(а):

Чувствую себя злостным оффтопщиком...
Можно попросить модераторов выделить задачу в отдельную тему?

Цитата:

Результат у Вас получится не плохой, но не оптимальный

В задаче надо ставить четкие критерии, а не расплывчатые цели.
Я исходил из следующего:если игрок угадал максимум, то он получает рубль, если нет - рубль отдает. Цель - максимизировать прибыль.

если игрок угадал максимум, то он получает рубль, если нет - ничего не получает и ничего не отдает. Правильная постановка.
Вый, кстати, получили оптимальную стратегию задачи: "Можете выбрать одно число. И ето ваша прибыль" (Без условие, чтобы было наибольшее). Но задача то другая. Интересный вариант - как изменится стратегия, если угадать максимум получаем не рубль, а стоимость числа