Comments 19
Немного царапают глаз опечатки в математических терминах. Мы все-таки пишем "коэффициент" и "Теорема Штурма является важным теоретическим результатом алгебры.", наверное не только алгебры, но и всей науки математики в целом. Поправьте, пожалуйста: "Таким образом, в первом случае индекс точки убудет ровно на 
единиц, в то время как во втором на 
или на 
единиц, тоесть на одно из близжайших 
четных чисел." "Тем самым, корни 
, чередуясь с корнями 
, как-бы отделяют корни полинома 
друг от друга." "Несмотря на простое определение сепаратора, увы, мне до сих пор не удалось найти разумный способ его построения. Под разумным я понимаю такой, где каждый коеффициент сепаратора является многочленом (возможно, однородным) от коеффициентов исходного полинома." - спорное утверждение, особенно на вещественной прямой с конечной точностью вычисления. Скорее всего, здесь придется использовать какую-либо быстро сходящуюся формулу нахождения корней по схеме Ньютона с численными производными.
Не совсем понял, что именно требуется уточнить.
Скорее всего, здесь придется использовать какую-либо быстро сходящуюся формулу нахождения корней по схеме Ньютона с численными производными.
Как раз нет, метод Ньютона совершенно неприменим здесь. По многим причинам. Как минимум, он может вовсе не сойтись ни к одному корню полинома.
Насчет метода Ньютона ошибся, извиняюсь. Но сам принцип "чередуем корни производной с корнями самого уравнения" для ряда плохо вычисляемых левых частей опять же ведет к модификации теоремы Штурма уже в работах года 1996-1998, где и сходимость и скорость сходимости поправили изрядно.
Опять же не совсем понял, что имеется ввиду под левыми частями. Может остатки от деления полиномов ?
модификации теоремы Штурма уже в работах года 1996-1998
О каких модификацях идет речь ? Если еще остались ссылки, или еще помните названия работ, можете их назвать ? Весьма интересует все что связано с этой темой.
По окончании процесса рекурсивного деления считается, что все корни найдены.
Противоречит требованию
КРИТЕРИЙ 1: Алгоритм должен быть абсолютно надежен, тоесть для любых коеффициентов, являющихся вещественными числами, все корни полинома гарантированно должны быть найдены. Исключение могут составлять случаи, когда корни настолько велики, что значение полинома в них не может быть вычислено по тем или иным причинам.
А именно, корни кратности два всегда будут пропущены.
Нет, корни полинома P(x) кратности 2 и любой другой никогда не будут пропущены. Другое дело, что невозможно, найдя корень, точно утверждать, какой он кратности. На мой взгляд, практически для любых практических применений это совершенно не критично и вполне удовлетворяет требованию. Если это не так, то было бы интересно услышать, где кратность корня требуется вычислить.
Будут, конечно. Попробуйте свой алгоритм для многочлена (x - 1/3)^2 начиная с отрезка 0,1 и используя просто середину отрезка как функцию.
У Вас получится отрезок (1 - eps)/3, (1 + 2*eps)/3 , значения в обоих концах положительные. А значит по алгоритму
Как и ранее, итерационный процесс прекращается, если
становится меньше заранее заданного
. Далее ступень считается относящейся к вещественному корню полинома, если на границах интервала значение полинома имеет разные знаки, или если разность значений индекса точки на границах интервала нечетна.
он будет выкинут.
Правда, такой корень он пропустит в любом случае. И скорее всего пропустит любой корень четной кратности. Справедливое замечание. Я не придавал и по-прежнему не придаю этому большого значения, по тем причинам, которые рассмотрены в разделе об ограниченной точности вычислений. Стоит немного изменить коэффициенты, и такой корень уже будет найден. А невозможность находить корни кратной четности это не проблема в практических применениях, так как в таких точках полином не пересекает нулевую прямую, а лишь касается ее. Нужно будет подумать над коррекцией указанного критерия.
Метод интересный, но слишком переусложненный. Этот ваш индекс весьма интересный объект, и его свойства стоит исследовать. Но есть более простой метод нахождения корней не уступающий вашему.
Раз уж вы все-равно вычисляете все производные и бинарный поиск делаете, то можно делать проще: сначала найдем все корни производной тем же методом. Потом между любыми двумя точками из этого списка ищем корень полинома бинарным поиском (если значения в концах имеют разный знак).
Еще могут быть корни левее первой и правее последней экстремальной точки. Вы этот вопрос обсуждаете в секции "на всей числовой прямой". Но тут тоже есть более простое решение. Надо посмотреть на знак производной и полинома снаружи от крайней точки. Если они разные - то там есть корень. Потом надо взять интервал длиной, например, 1 снаружи от крайней точки. Потом надо удваивать длину этого интервала, пока полином не поменяет знак на второй границе и потом опять же бинпоиском найти корень. Вот так за O(log (M/eps)) вы найдете сначала интервал, где корень есть, а потом и сам корень. M тут - максимальное значение корня по модулю. Такая же оценка как в вашем алгоритме.
Технически, метод бинарного поиска последовательно по производным (как и метод, предложенный автором) плохо работает с многочленами вида (x-a)^(2n) + b, где b очень мало по модулю (по сравнению с a). И это, в общем-то понятно - очень малое изменение свободного члена (с a^2n - b до a^2n + b) сильно влияет на корни (либо они есть, либо их нет). Собственно, никакой численный метод, кажется, за заранее фиксированное число шагов не может ответить на вопрос, есть ли корень, или нет (и это то, зачем нужен метод Штурма).
Заинтересовали. Собственно как метод Штурма поможет здесь ? Даже чисто по определению оба метода заключаются в вычислении числа перемен знака в ряду значений полиномов, только в методе Штурма каждый следующий полином зависит от предыдущего, что вносит (по опыту, часто очень большие) отклонения в сами коэффициенты полиномов ряда, не говоря уже о корнях. Не совсем понятно, что значит "плохо работает с многочленами вида" ?
У многочленов (x-a)^(2n) + b и (x-a)^(2n) - b совпадают и все производные, и знаки значений почти всюду (кроме очень маленького интервала a +- b^{1/(2n)} ).
В итоге, любому численному методу нужно попасть в этот интервал (который может быть гораздо меньше eps), чтобы отличить эти два случая.
А метод Штурма даёт точку из этого интервала (собственно a) когда делит многочлен на производную, там как раз (x-a) выносится.
Вы этот вопрос обсуждаете в секции "на всей числовой прямой".
Не совсем, в этой секции я рассматриваю, как правильнее действовать в ситуации, когда полином имеет, например, вид: 
. А именно, когда допустимое значение переменной
это в том числе некоторая окрестность нуля. Чтобы был яснее, откуда такие ситуации - попробуйте составить уравнение пересечения луча, по направлению
из точки
с неявно заданной поверхностью
, то есть уравнение 
. Если поверхность имеет направленные в бесконечность части, как раз возникает указанная проблема. Есть и другие геом. задачи, в которых присутствует такая проблема, пути решения которых как раз обсуждаются в указанной секции.
Потом надо взять интервал длиной, например, 1 снаружи от крайней точки. Потом надо удваивать длину этого интервала, пока полином не поменяет знак на второй границе и потом опять же бинпоиском найти корень.
Такой подход вполне может быть таким же эффективным, как и указанные мной, однако, на мой взгляд такая реализация сложнее в программировании. Можно сказать вкусовщина.
Раз уж вы все-равно вычисляете все производные и бинарный поиск делаете, то можно делать проще: сначала найдем все корни производной тем же методом. Потом между любыми двумя точками из этого списка ищем корень полинома бинарным поиском (если значения в концах имеют разный знак).
На самом деле да, отличия не столь значительны, однако они имеются (До какой то степени это тоже можно назвать вкусовщиной). В частности, в разделе о полиномиальной оценке указаны очень эффективные способы вычисления полинома и всех его производных в точке. Для указанного вами метода приходится пользоваться либо простой схемой Горнера, либо схемой Эстрина. В разделе о частичном индексе также указан один из способов оптимизации. Насколько такие методы отличаются в вычислительной эффективности - вероятно, не очень сильно.
В частности, в разделе о полиномиальной оценке указаны очень эффективные способы вычисления полинома и всех его производных в точке
Ваш алгоритм вычисляет все производные и полином в O(N^2 Log M/eps) точках. Суммарно O(N^3 log M/eps) вычислений полином. Мой же в N раз меньше: каждую производную отдельно в O(N Log M eps) точках. Чуть сложнее оценки, ведь степени производной уменьшаются на один каждый раз, но в целом все тоже самое будет. Поэтому даже без хитрых схем там одинаковое количество вычислений.
Я имел ввиду не асимптотическую сложность, а абсолютную. Значения полинома и всех его производных вычисляются без лишних домножений на коэффициенты, при этом с высокой локальной параллельностью, которой не достичь с обычными производными.
Я тоже имел ввиду количество вычислений, но мне лень было точные значения приводить, поэтому я привел оценки.
Ваш индекс меняет знак максимум в n(n+1)/2 точках. И вы для каждой точке в итоге прогоняете бинпоиск. Т.е. вы считаете все полинномы в n(n+1)/2 log(M/eps) точках.
В моем методе надо найти лишь n корней самого полинома, поэтому там n бинпоисков, и весь полином считается лишь n log(M/eps) раз. Но до этого сначала надо найти корни производной. Там n-1 отрезок и этот полином считается в (n-1) log(M/eps) раз. И т.д. в итоге мой метод считает один полином в n(n+1)/2 log(M/eps) точках, вы же считаете в таком же количестве точек ВСЕ полиномы. Даже если ваша схема сокращает вычисления в n раз, по сравнению с тупо вычислением всех полиномов отдельно, то она только-только сравняется с предложенным методом по количеству вычислений.
Хабр - торт!
В первом приближении очень интересно, и, хотя мне оно на практике не требуется, но вкурить в эти методы, однозначно, стоит. Завтра на свежую голову в рабочее время вместо работы :))
Надежный алгоритм поиска корней полиномиальных уравнений