Уравнение Пелля. Поиск корней

Автор темы a62 
ОбъявленияПоследний пост
ОбъявлениеРаботодателям и кадровым агентствам: Размещение вакансий26.03.2008 03:07
ОбъявлениеЗапущен новый раздел «Задачки и головоломки»29.08.2019 00:42
ОбъявлениеТеХнический редактор - LATEX18.01.2020 21:57
30.12.2019 22:21
Уравнение Пелля. Поиск корней
Чтобы не заставлять читать лишнее и ненужное, объясню, о чем тема.
В рамках школьного «индивидуального проекта» сын писал программу, а я немного помогал ему с алгоритмом. По ходу работы возникли некоторые мысли, которые мы реализовали в программе, а часть не реализовали.
Часть из этих "мыслей" удалось найти в теории и дать ссылки. Часть осталась неопределенной. Конкретно: есть несколько утверждений и формул, которые в поиске в явном виде найти не удалось. В связи с этим ниже их излагаю (максимально кратко). При этом надеюсь в ходе обсуждения темы, возможно, выяснится, что:
- мы «изобрели велосипед»
- или ошиблись
- или не ошиблись, но есть более эффективные пути, т.е. был сделан «напрасный труд».
Ну или, уж на крайний случай, возможно, окажется, что что-то из всего этого «не напрасно».
Заранее прошу не судить строго оформление, с LaTex познакомился совсем недавно.
Текст большой, я его много раз вычитывал и проверял формулы, но, все равно, возможны ошибки. Если укажете, постараюсь исправить

Уравнение Пелля $ y^2=Ax^2+1$. Рассматривать будем только натуральные числа $A,x, y$. Ищем нетривиальные корни. В первую очередь наименьший из них.
Здесь и дале $x_0$ и $ y_0$ – наименьшее нетривиальное решение уравнения Пелля.
И еще далее везде $ A=n_0^2+c_0$ где $n_0$ наибольшее число, квадрат которого меньше $A$.
Чтобы не путать участников форума, - обозначения: подходящие дроби ${q_i, P_i}$ разложения в цепную дробь числа корня из числа ${\sqrt{A}= \overline{a_0, a_1 , a_2 … }}$.
По классическому определению, нужные значения ${x_0=q_{i } ; y_0=P_{i }}$ надо искать при таком $i$, которое соответствует предпоследнему значению в периоде $a_i$.
Но мы использовали другой критерий нахождения нужной дроби.
Схематично, «цепные» вычисления образуют последовательность действий для получения массивов ($a_i,q_i,P_i$), которая сводится к такой схеме:
$a_i=[\frac{n_0+n_{i-1}}{c_{i-1}}] → n_i=a_i c_{i-1}-n_{i-1} → c_i=\frac{A-n_i^2}{c_{i-1}}$ при этом ($A=n_0^{2}+c_0; a_0=n_0$) ------------------------------ (1)
Остается только «досчитать» по Эйлеру:
$q_i =a_i q_{i-1}+q_{i-2} → P_i=a_i P_{i-1}+P_{i-2}$ при этом $(q_{-1}=0; P_{-1}=1; q_0=1; P_0=n_0)$ ------------------------------ (2)
Здесь и далее $i$ индекс в расчетах в алгоритме цепных дробей. Расчет начинается с индекса $0$. В частности, $a_0=n_0$, но период $a_i$ начинается с $a_1=[\frac{2n_0}{c_0} ] $ и заканчивается $a_N=2n_0$ Само значение $a_0$ в период не входит.
Утверждение $I$ : $с_i$ лучше показывает ход поиска корня, чем $a_i$ .
Важные свойства $c_i $:
Свойство 1:
$с_i, q_i$ и $P_i$ связаны соотношением: $ P_i^2-Aq_i^2 =(-1)^{i+1} с_i$ ------------------------------ (3)
Это означает, что если $(-1)^{i+1}с_i=+1$ то $q_i$ и $P_i$ – корни уравнения Пелля.
Свойство 2:
$ c_i$ периодичен и длина $N$ его периода такая же, как и для $a_i$. Числа периода $c_i$ отличаются от чисел периода $a_i$ и начинается период $c_i$ всегда со значения $c_0$, а заканчивается всегда значением $с_{N-1}=1$
Свойство 3:
Если из периода $c_i$ исключить последнюю цифру (единицу), то в остальной части он (период) обладает центральной симметрией: $с_i= с_{N-i-2}$.
Свойство 4
Если в ходе вычислений оказалось, что $с_i=с_{i-1}$, то это значит, что $N=2i+1$ и:
1) $с_{(N-3)/2}= с_{(N-1)/2}$
2) $(-1)^{N} с_{N-1}=-1$
3) период $N$ нечетный, а ($q_{2N-1}; P_{2N-1}$) – наименьшее нетривиальное решение уравнения Пелля для данного $A$ находится в конце второго периода и $(-1)^{2N} с_{2N-1}=+1$;
4) $ q_{(N-3)/2}^2+q_{(N-1)/2}^2=q_{N-1}; P_{N-1}=\sqrt{Aq_{N-1}^{2}-1} $; ------------------------------ (4)
5) $ x_0=q_{2N-1}=2q_{N-1} P_{N-1}$; $ y_0=P_{2N-1}= P_{N-1}^2+1$ ------------------------------ (5)
В этом случае $A$ принадлежит последовательности A003814 (по международной классификации) и соответственно всегда представимо в виде суммы квадратов. Это условие необходимое, но не достаточное. Условно, назовем эти числа первым типом. Их в первой тысяче насчитывается $152$.
Если $A= m^2+n^2$, то число $с_{(N-3)/2}= с_{(N-1)/2}$ - всегда одно из этих двух ($m$ и $n$), причем всегда нечетное.
Все остальные числа $A$ условно назовем вторым типом. Среди них числа, представимые в виде суммы квадратов встречаются, но крайне редко. Например, в первой сотне ко второму типу относятся $18$ и $34$.
Для второго типа чисел $A $ справедливы другие свойства:
Свойство 5:
Если в ходе вычислений оказалось, что $с_{i-2}=с_i$, то это значит, что $N=2i$ и:
1) $с_{N/2-2}= с_{N/2}$
2) период четный, $(q_{N-1}; P_{N-1})$ – наименьшее нетривиальное решение уравнения Пелля для данного $A$ находится в конце первого периода, $(-1)^{N} с_{N-1}=+1$;
4) $ P_{N/2-1}=\frac{с_{N/2-1}}{2}(q_{N/2-2}+q_{N/2})$; ------------------------------ (6)
5) $x_0=q_{N-1}=\frac{2q_{N/2-1} P_{N/2-1}}{c_{N/2-1}} ; y_0=P_{N-1}= P_{N/2-1}^2+(-1)^{N/2} $ ------------------------------ (7)
Определение: Условно назовем уравнение $y^2=Ax^2±s$ «предшествующим» по отношению к уравнению Пелля, если его наименьшее нетривиальное решение $(x_s,y_s)$ связано с решением самого уравнения Пелля при одном и том же $A$ соотношением:
$x_0=\frac{2x_s y_s}{s}$, $s \in \mathbb{N}$ ------------------------------ (8)
Утверждение $II$.
Для любого значения $A$ у уравнения Пелля есть «предшествующее» уравнение, причем $s$ – делитель $A$ или удвоенный делитель $A$, а корни $x_s$ и $y_s$также являются парой ($q_i$ и $P_i$) подходящей дроби.
Значение $s$ равно: если $ N$ четное, то $s=c_{N/2-1}$; если $N$ нечетное, то $s=-1$.
При этом для простых $A$ значение $s$ равно строго одной из трех величин: $ -1$, $-2$ или $+2$.


Частные решения.
Введем функцию$ f(A)$, у которой аргумент - число $A$ из уравнения $y^2=Ax^2+1$, а значение функции - величина $x_0$ наименьшего нетривиального решения $ f(A) = x_0$ (для этого $A$)
Параметр $A$ можно представить в виде $A =n_0^2±с_0$, $(n_0, с_0) \in \mathbb{N}$.
Тогда для большинства частных решений уравнения Пелля есть простая формула:
$f(n_0^2±с_0 )= \frac{2n_0}{с_0}$ , если $\frac{2n_0}{с_0}$ – целое число ------------------------------ (9)
Единственное исключение $f(n_0^2-1)=1$
Есть еще 2 частных случая решений для отдельных значений $A=n_0^2+с_0$:
если $с_0=4$, то ответ: $x_0=\frac{n_0}{2}$ (для четного $n_0$) или $x_0=\frac{n_0 (An_0^2+3)}{2}$ (для нечетного $n_0$) ------------------------------ (10)
если$ с_0=2n_0-3$, то ответ: $x_0= \frac{n_0+1}{2}$ (для нечетного $n_0$) или $x_0=\frac{n_0^2}{2}+n_0$ (для четного $ n_0$) ------------------------------ (11)
Частные решения устраняют проблему «коротких периодов», которые усложняют алгоритм поиска корня по критериям $с_{i-1}= с_i$ и $с_{i-1}= с_{i+1}$1
Частный случай $c_i=4$
Если $с_i=4$, можно ускорить поиск корня $x_0$ по такой схеме:
Производим вычисления по основному алгоритму, доходим до $с_i=4$, а значит получаем соответствующие ему $q_i$ и $P_i$ и окончательный ответ уже вычисляем сразу по формулам:
- для $(-1)^{i+1} c_i=+4$: Если $q_i P_i$ четно, то $x_0=\frac{q_i P_i}{2}$ , иначе $x_0=\frac{q_i (P_i^2-1)}{2}$ ------------------------------ (12)
- для $(-1)^{i+1} c_i=-4$: Если $q_i P_i$ четно, то $x_0=\frac{q_i P_i}{2}$, иначе $x_0=\frac{q_i P_i (Aq_i^2 P_i^2+3)}{2}$ ------------------------------ (13)
Подходящие дроби, для которых $с_i=4$, всегда находятся в цепочке вычислений раньше всех других описанных критериев.
Для примера: $A=157$ длина периода $N=17$, нужная подходящая дробь на позиции $i=2\cdot17-1=33$, т.е. в базовых вычислениях $33$ цикла, а $(-1)^{i+1} c_i=-4 $ находится на позиции $i=4$, т.е. в $5$м цикле. Для $A=421$ требуется $74$ цикла, а $(-1)^{i+1} c_i=-4$ находится в $13$м цикле.

https://c.radikal.ru/c18/1912/8e/de73db14ca70t.jpg
Возможно есть еще частные решения в явном виде, но мне они не известны. Говорят, есть методичка по уравнению Пелля и там есть все частные случаи, в т.ч. все вышеописанные. Но я эту методичку в открытом доступе не нашел.
Для максимально эффективного применения предложенных критериев нужно было бы заранее определить тип числа $A$ и факт того, можно ли встретить среди $c_i=4$. Ни для того, ни для другого применимых способов в литературе найти не удалось.

Множество всех решений уравнения Пелля.

В теории рекомендовано находить все корни путем возведения в степень. Мы решили, что удобнее будет использовать рекуррентные соотношения (если известно наименьшее нетривиальное решение $x_0$ и $y_0$) :
$y_i = 2y_0 y_{i-1} – y_{i-2}; x_i = 2y_0 x_{i-1} – x_{i-2}$ при этом $(x_{-1}=0; y_{-1}=1)$ ------------------------------ (14)

Пробегая по всем натуральным $i$ начиная с $1$, мы получаем все решения строго по возрастанию. Предшествующие уравнения (и только они) имеют такие же рекуррентные соотношения. Но у них другие значения ($x_{-1}$ и $y_{-1}$).

Из того, что оказалось интересным, но найти примерение не получилось:

"Максимальные" значения.

Пусть $A$ растет последовательно и для каждого $A$ мы вычисляем $f(A)$. Начнем процедуру с $A=61$. Вычисляем $f(61)=226 153 980$. Это первый член последовательности. Следующим будет минимальное из чисел, $A$, таких, что $f(A)>f(61)$ и т.д. Т.е. все значения функции между двумя максимальными меньше их обоих.
Вычисления показывают, что логарифм максимальных значений очень близок к величине $3\sqrt{A}$, вплоть до значений порядка $10^8$.
http://c.radikal.ru/c01/1912/c2/afbfd9393b3b.jpg
В дальнейшем коэффициент при корне квадратном растет, но очень медленно. Т.е. для аппроксимации вполне можно использовать оценку $\lgx_{max}\approx3\sqrt{A_{max}}$
Этой же формулой, наверное, можно пользоваться для оценки сверху функции $x_0=f(A)$
Все, без исключения, максимальные числа, начиная с $61$, представляют из себя простые числа типа $12k+1$, где $k\in \mathbb{N}$. И все они первого типа ( $ s=-1$), соответственно, период у всех нечетный.
Все "максимальные" числа от $1$ до $1 000 000$:
http://b.radikal.ru/b11/1912/5a/477ef40a093c.jpg
Интересно, что все "максимальные" числа заканчиваются только на $1$ или на $9$. Забавно,также, что, начиная с числа $421$, остаток от деления на $120$ дает только $61$ и $109$, которые сами являются "максимальными" числами.

Буквально несколько слов о времени вычислений:
При расчете всех корней для значений $A$ от $1$ до $1 000 000$ (1 млн.) время вычислений в одинаковых условиях (один и тот же компьютер Intel Pentium CPU G4620 @ 3.70 GHz):
Если считать просто по методу цепных дробей по рекуррентным формулам подходящие дроби $q_i$ и $P_i$ : $215$ секунд
Если считать только $q_i$ : $180$ секунд.
Если применить критерий $с_i=с_{i-2}$ : $90$ секунд.
На всякий случай: программа написана на Питоне, использована его библиотечная "длинная арифметика", использовалась только целочисленная арифметика.

Итого мы пришли к возможности сокращения вычислений по алгоритму «цепных дробей» при помощи использования свойств $c_i$:
1. Частные решения при отдельных значениях $c_0$ Они выполняются до начала циклов и носят «разовый характер». Благодаря им сокращается время вычислений и решается неприятная проблема «коротких периодов». Мы успешно вписали формулу частных решений в реализацию программы ;
2. Критерий $c_{i-2}=c_i$ позволяет сократить вдвое количество циклов и для первого типа (нечетный период) и для второго (четный период). Позволяет вообще не считать $P_i$. Естественно, был успешно реализован;
3. Критерий $c_{i-1}=c_i$ позволяет сократить еще вдвое количество циклов, но только для первого типа (нечетный период), а это примерно $15%$ (по крайненй мере в первой тысяче это так). Для второго типа (четный период) является лишним действием во всех циклах и, соответственно замедляет вычисления. Если я правильно провел поиск, пока никто не умеет эффективно определять заранее, до вычислений, тип числа $A$. Формулы реализации критерия $c_{i-1}=c_i$ требуют возведения в квадрат и извлечения корня, хотя и всего один раз, но в длинной арифметике этот фактор мы оценивать не стали и от использования критерия отказались;
4. Критерий $c_{i-1}=4$ Позволяет очень существенно сократить количество вычислений. Но тоже годится лишь для чисел $A$, в вычислениях которых есть $c_{i-1}=4$. Для остальных это лишние действия в каждом цикле. Сколько таких чисел – не очень понятно как оценивать. Поэтому не стали анализировать оптимизационные тонкости с точки зрения информатики, просто не стали и этот критерий использовать.
5. Рекуррентные формулы для расчета всех корней, когда наименьшие найдены, нам показались гораздо эффективнее, чем общепринятая схема возведения в степень.
Возможно мы что-то упустили и с критериями можно было поступить более рационально? Также остались все вопросы, которые изложены в самом начале темы.
Были, на самом деле, и другие вопросы, но о них, возможно, в другой раз.



Редактировалось 1 раз(а). Последний 31.12.2019 03:16.
Извините, только зарегистрированные пользователи могут публиковать сообщения в этом форуме.

Кликните здесь, чтобы войти