Dmitriy40Спасибо, что откликнулись. Вот это сильно обнадеживает:
Для не слишком больших чисел, скажем до
, извлечение целочисленного квадратного корня не проблема, не быстро, но миллионы в секунду вполне реально (например до
корень извлекается 15 млн.раз/с), это у меня уже давно написано.
Мне нужно аккуратно прикинуть, насколько велики будут коэффициенты решаемых квадратных уравнений, если исходить из границы
. Надеюсь уложиться в Ваши размеры. Позже я напишу подробное изложение алгоритма; он в самом деле крайне примитивный, ничего общего с (действительно) очень сложным алгоритмом отыскания целых точек на эллиптической кривой общего вида. А сейчас выложу код на Maple, чтобы Вы могли чисто визуально (а может, и содержательно) понять, что с точки зрения математики здесь все примитивно.
Используются две вспомогательные процедуры. Первая находит целую часть квадратного корня из натурального
:
Код:
FloorSqrt:=proc(N)
local x,x_next;
if N=0 then return 0 end if;
x:=iquo(N+1,2); x_next:=iquo(x+iquo(N,x),2);
while x>x_next do
x:=x_next; x_next:=iquo(x+iquo(N,x),2)
end do;
return x
end proc:
Вторая решает квадратные уравнения
с целыми коэффициентами в целых числах:
Код:
solve_int:=proc(a,b,c)
local S,d,d2;
if a=0 and b=0 then return {} end if;
if a=0 then if irem(c,b)=0 then return {-c/b} else return {} end if end if;
S:={};
d:=b^2-4*a*c;
if d>=0 then d2:=FloorSqrt(d);
if d2^2=d then
if irem(-b+d2,2*a)=0 then S:=S union {(-b+d2)/2/a} end if;
if irem(-b-d2,2*a)=0 then S:=S union {(-b-d2)/2/a} end if;
end if;
end if;
return S;
end proc:
Ну и, собственно, сама программа:
Код:
runge:=proc(H)
local S,m,Q,k,res,x,y;
S:={[0,-1]};
m:=1+floor(evalf(abs(H)^(1/4)));
Q:=floor(evalf(m/(m^2-2)+sqrt((m+abs(H))/(m^2-2)+2/(m^2-2)^2)));
for k from -m to m do
res:=solve_int(k^2-2,2*k,H-k+1);
for x in res do S:=S union {[x,-k*x-1]} end do
end do;
for x from 1 to Q do
res:=solve_int(x,1,-2*x^3+H*x+1);
for y in res do S:=S union {[x,y]} end do
end do;
for x from -Q to -1 do
res:=solve_int(x,1,-2*x^3+H*x+1);
for y in res do S:=S union {[x,y]} end do
end do;
return S
end proc:
По существу, здесь решается две порции квадратных уравнения (выделены вертикальными пробелами), причем вторая порция --- это исходное уравнение, решаемое относительно
(а
в роли параметра, который перебирается в заранее вычисленных границах от
до
, исключая
). Первая порция квадратных уравнений зависит от вспомогательного параметра
(который перебирается в заранее вычисленных границах от
до
), а неизвестным в уравнении является
(по известному
и найденному корню
затем вычисляется
).
Если бы Вы заменили мои вульгарно написанные (с точки зрения программирования; математически они корректны) вспомогательные процедуры Вашими оптимизированными, можно было бы покуситься на рекорд.