Как писать быстрые программы

[Dmitriy40]

Dmitriy40 в сообщении #1718536 писал(а):

Выигрыш GPU 130 раз уменьшится до 130/15=8.7 раза если сравнивать с avx(2) версией кода выше. А на поток - CPU быстрее в 30 раз. т.е. чтобы обогнать 64 потока CPU надо GPU с 2000 потоками минимум. Хотя это сравнение несколько разных кодов не слишком адекватно.

Прислали результаты запуска OpenCL кода на RTX4060, получилось она быстрее моего CPU в один поток на 3.5ГГц в 70 раз. А в расчёте на поток GPU медленнее в 45 раз.
Если учесть и частоты, то в тактах получается GPU в 30 раз медленнее CPU.
Значит мой сервер где-то в 2.5 раза проигрывает RTX4060.

[Yadryara]

Dmitriy40 в сообщении #1719329 писал(а):

заводим вектор флагов s[] что число цепочки точно не подходит или разложено полностью;
в простых фильтрах указываем какие числа полностью разложились чтобы не раскладывать повторно;
после простых фильтров запускаем разложение каждого недоразложенного числа на секунду (меньше нельзя) командой alarm() и проверяем разложилось ли оно;
если разложилось и не дало нужное количество делителей или и не может уже дать - тоже выставляем флаг;
если набралось много недопустимых чисел (т.е. valids точно мал и такие не интересны) - обрываем проверку;
повторяем разложение оставшихся чисел уже с ограничением на 3 секунды;
повторяем разложение оставшихся чисел уже с ограничением на 15 секунд;
повторяем разложение оставшихся чисел уже с ограничением на 90 секунд;
если valids оказался достаточно большим - запускаем полное разложение оставшихся чисел.

Спасибо. Ну, вы видели: актуальные числа нынче 62-значные. И очень вряд ли надо лезть выше. И это итоговые числа цепочки, а недоразложенные-то обычно на 3-5 порядков меньше.

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

Dmitriy40 в сообщении #1719329 писал(а):

f=factor(n+d-1,2^15); fn=matsize(f)[1];\\Не успело, проверим все ли делители найдены, до 2^15 проверяется очень быстро

Ну вот этим кодом очень здорово проверяется делимость на все предвычисленные простые до 2^(19-20). Так что в этой части вроде ничего не ускорить.

Счётчики успеха пока оставил.

for(i:small = starti, starti + koli - 1, 

for(j = 1, #prter, if( bo[j] == i % prter[j], next(2)));
\\usp[1]++;

for(j = 1, #T,     if(  T[j] == i % dopp[j] , next(2)));
\\usp[2]++;

kanp = n0 + m * i;
if(!ispseudoprime(kanp,1), next);
n = v[ip+1] * kanp - ip;
\\usp[3]++;

ostmes = ostmesbig;


kolfac=vectorsmall(lenv,mes,1); pro=vector(lenv,mes,1);

for(nompr=1, kolpredp,
mes = lenv - (ostn0a[nompr] + ostma[nompr]*i - ip + lenv - 1) % predp[nompr];
if(mes<1, next); kolfac[mes]++; pro[mes] *= predp[nompr]; 
\\usp[4]++;

if( kolfac[mes] == nu[mes], 
\\usp[5]++;

if(!ispseudoprime((n+mes-1)/(v[mes]*pro[mes]),1), next(2),
ostmes = setminus( ostmes, [mes] ))));
\\usp[6]++;

for(j=1, #ostmes, if(ispseudoprime((n+ostmes[j]-1)/(v[ostmes[j]]*pro[ostmes[j]])),next(2)));

hormes = [];
for(j = 1, #ostmes,
if(nu[ostmes[j]] - kolfac[ostmes[j]] == 1,
if( #factor( (n+ostmes[j]-1)/(v[ostmes[j]]*pro[ostmes[j]]) )[2] > 2, next(2),
hormes = concat( hormes, ostmes[j]);
\\print(hormes, "   ", #hormes);
)));

ostmes = setminus( ostmes, hormes );

print1(#ostmes, "   " );

for(j = 1, #ostmes,
if(nu[ostmes[j]] - kolfac[ostmes[j]] == 3,
matfac = factor( (n+ostmes[j]-1)/(v[ostmes[j]]*pro[ostmes[j]]) );
fac2 = #matfac[2];
if( fac2 < 4 && vecprod(matfac[2]) == 1,  nedobor++;print("nedobor = ",nedobor);  next(2));
if( fac2 > 4, prev++;print("prev = ",prev); next(2))));

И дальше уже идёт проверка на valids по всему полю и показ найденных цепочек.

[Dmitriy40]

Yadryara в сообщении #1719337 писал(а):

Ну вот этим кодом очень здорово проверяется делимость на все предвычисленные простые до 2^(19-20). Так что в этой части вроде ничего не ускорить.

Здесь не ускорение, здесь накопление в f[] всех найденных делителей, не только до 2^15, но и больших если они были найдены предыдущим factor() под alarm(). Тут вообще лучше поставить 2^24 чтобы получились все найденные делители (больше 2^24 сохраняются из-за factor_add_primes=1, а меньшие найдём тут). Но это не принципиально: если за секунду не разложилось, то вряд ли там есть делители в диапазоне 2^15...2^24, и даже если и есть, то эта позиция будет ещё раз проверена следующим alarm(factor()) уже с учётом всех найденных больших делителей, а factor() без ограничения (под alarm) малые делители находит моментально (не перебором как при указании ограничения, а методом ECM) и они не станут тормозом.

Yadryara в сообщении #1719337 писал(а):

И дальше уже идёт проверка на valids по всему полю и показ найденных цепочек.

Ну вот alarm+factor (или кстати сразу numdiv, это без разницы) вставлять перед проверкой valids по всему полю.

[Yadryara]

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

Ну вот вставил засечку времени перед первым factor:

if(nu[ostmes[j]] - kolfac[ostmes[j]] == 1,
t3=getwalltime();
if( #factor( (n+ostmes[j]-1)/(v[ostmes[j]]*pro[ostmes[j]]) )[2] > 2, print(">2  ",strtime(getwalltime()-t3));next(2),
print("   ",strtime(getwalltime()-t3));hormes = concat( hormes, ostmes[j]);

   1 ms
>2  242 ms
   80 ms
   1,683 ms
   1,504 ms
>2  148 ms
   2 ms
   3 ms
   44 ms
>2  119 ms
   3 ms
   9 ms
>2  136 ms
>2  14 ms
>2  184 ms
>2  32 ms
   1 ms
>2  117 ms
   2,004 ms
   1,138 ms
>2  85 ms
   125 ms
>2  718 ms
>2  141 ms
>2  277 ms
>2  143 ms
>2  56 ms

Когда факторов больше 2-х, то работает быстро.

[Dmitriy40]

Yadryara
А каким образом у Вас вообще работает код #factor()[2]?! У меня он выдаёт ошибку:

? #factor(113*5)[2]
  ***   at top-level: #factor(113*5)[2]
  ***                               ^---
  ***   incorrect type in _[_] OCcompo1 [not a vector] (t_MAT).
  ***   Break loop: type 'break' to go back to GP prompt
break> break

? factor(113*5)[2]
  ***   at top-level: factor(113*5)[2]
  ***                              ^---
  ***   incorrect type in _[_] OCcompo1 [not a vector] (t_MAT).
  ***   Break loop: type 'break' to go back to GP prompt
break> break

И это правильно:
если порядок выполнения #(factor()[2]), то [2] неприменимо к результату factor, который матрица, а не вектор (как и сказано в ошибке выше);
если порядок выполнения (#factor())[2], то [2] неприменимо к результату #, который число, а не вектор.

[Yadryara]

Да, у меня в интерпретаторе PARI тоже не работает:

(18:49) gp > factor( 123456 )
%2 =
[  2 6]

[  3 1]

[643 1]

(18:49) gp > #factor( 123456 )
%3 = 2
(18:49) gp > #factor( 123456 )[2]
  ***   at top-level: #factor(123456)[2]
  ***                                ^---
  ***   incorrect type in _[_] OCcompo1 [not a vector] (t_MAT).
(18:51) gp >

Ну вот вроде здорово получается: Транслятор и Убунта понимают меня так как надо. Я пошагово перепроверял сравнительно недавно. То есть для числа 123456 имеется три фактора и вся цепочка будет отброшена, потому что 3 > 2.

-- 03.03.2026, 19:06 --

Кстати, попытался alarm применить, но он не сработал как раз на том числе, для которого факторизация длилась дольше секунды.

alarm(1,#factor( (n+ostmes[j]-1)/(v[ostmes[j]]*pro[ostmes[j]]) )[2])

   2 ms
>2  142 ms
   60 ms
  ***   at top-level: init_Rab_243_20_0_test()
  ***                 ^------------------------
  *** init_Rab_243_20_0_test: forbidden comparison t_INT , t_ERROR.
  ***   Break loop: type 'break' to go back to GP prompt
break>

-- 03.03.2026, 19:10 --

А ежели переделать на 2 секунды, то этот аларм пока не вмешивается, видать потому что все времена меньше 2-х секунд:

(Оффтоп)

   1 ms
>2  176 ms
   95 ms
   1,499 ms
   1,703 ms
>2  146 ms
   1 ms
   2 ms
   56 ms
>2  123 ms
   3 ms
   6 ms
>2  109 ms
>2  20 ms
>2  190 ms
>2  32 ms
   1 ms
>2  91 ms
   1,816 ms
   983 ms
>2  75 ms
   82 ms
>2  774 ms
>2  138 ms
>2  308 ms
>2  119 ms
>2  43 ms
   92 ms
   1,140 ms
   452 ms
   64 ms
   17 ms
   1,806 ms
   1,266 ms

[Dmitriy40]

DemISdx в сообщении #1719354 писал(а):

У Вас в строке:
"if(type(alarm(1,factor(n+d-1)))!="t_ERROR",\\Дадим 1с на разложение"
Вроде как не совсем по канонам получается...

Почему не по канонам? Ведь проверяется тип результата, он или матрица от factor, или ошибка от alarm, я проверяю на <не_ошибку>, имею право. Ну и второй пример по ссылке ровно как у меня, только в две строки вместо одной и с лишней переменной. Вот дальше да, в примере лучше чем у меня, хотя работает и так и так.

-- 03.03.2026, 19:34 --

Yadryara в сообщении #1719353 писал(а):

Кстати, попытался alarm применить, но он не сработал как раз на том числе, для которого факторизация длилась дольше секунды.

Он как раз сработал - и выдал не число, а e_ALARM с типом t_ERROR, которые Вы попытались сравнить с двойкой и разумеется это ошибка (объект типа t_ERROR нельзя сравнивать с объектом типа t_INT). О чём и сказано в выводе.

-- 03.03.2026, 19:39 --

Yadryara в сообщении #1719353 писал(а):

Ну вот вроде здорово получается: Транслятор и Убунта понимают меня так как надо.

Это ПЛОХО! Они должны понимать согласно нормам языка, по которым выражение #factor()[2] недопустимо:

Dmitriy40 в сообщении #1719351 писал(а):

если порядок выполнения #(factor()[2]), то [2] неприменимо к результату factor, который матрица, а не вектор (как и сказано в ошибке выше);
если порядок выполнения (#factor())[2], то [2] неприменимо к результату #, который число, а не вектор.

И так писать программы не стоит.

[Yadryara]

Да, невнимательно посмотрел на alarm.

Задача.
Когда комп пытается факторизовать число и потратил уже секунду, прервать факторизацию, то есть считать что факторов не более 2-х и цепочку пока не отбрасывать, продолжить проверку других чисел цепочки.

[Dmitriy40]

Yadryara в сообщении #1719365 писал(а):

Задача.
Когда комп пытается факторизовать число и потратил уже секунду, прервать факторизацию, то есть считать что факторов не более 2-х и цепочку пока не отбрасывать, продолжить проверку других чисел цепочки.

Так приведённый мной цикл именно так и делает. Т.е. задача решена ещё 4 года назад и код выкладывал в облако, например вот для M48n21, там применён alarm. И показанный код практически этот и есть, с незначительными доработками.

-- 03.03.2026, 21:03 --

В нём только removeprimes вызывается неправильно, правильно так: removeprimes(addprimes()).

-- 03.03.2026, 21:06 --

Yadryara в сообщении #1719365 писал(а):

то есть считать что факторов не более 2-х

Это некорректно, делителей может быть и три и больше, но все относительно большие и за секунду ECM их не нашёл, ну бывает, ECM не гарантирует нахождение делителей за конечное количество итераций (и не гарантирует нахождение делителей строго в порядке возрастания, тому море примеров, один из последних: ECM нашёл делители из 36 и 39 цифр, но не нашёл из 34 цифры).

-- 03.03.2026, 21:25 --

Yadryara в сообщении #1719365 писал(а):

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

Впрочем, если так упростить задачу, то всё ещё проще:

for(d=0,#v-1, if(type(nd=alarm(1,numdiv(n+d))) <> "t_ERROR" && nd<>96, next(2)); );\\Если за секунду успели разложить, то проверяем подходит ли и если нет, то переход к следующей цепочке

[DemISdx]

Dmitriy40 в сообщении #1719360 писал(а):

DemISdx в сообщении #1719354 писал(а):

У Вас в строке:
"if(type(alarm(1,factor(n+d-1)))!="t_ERROR",\\Дадим 1с на разложение"
Вроде как не совсем по канонам получается...

Почему не по канонам? Ведь проверяется тип результата, он или матрица от factor, или ошибка от alarm, я проверяю на <не_ошибку>, имею право. Ну и второй пример по ссылке ровно как у меня, только в две строки вместо одной и с лишней переменной. Вот дальше да, в примере лучше чем у меня, хотя работает и так и так.

Я имел ввиду именно как "не совсем" по канонам.
Другими словами когда есть "однострочник", то бывает крайне сложно понять в какой именно точке "прилетает" нечто, что все ломает.
В таких случаях начинаю раскладывать все "на условные примитивы".
Тогда это позволяет определить уже конкретное место проблемы.

(Оффтоп)

Не понимаю пока, как именно, донести свою мысль.
Как обратный реверс:
1.

if(type(alarm(1,factor(n+d-1)))!="t_ERROR",\\Дадим 1с на разложение
      s[d]=numdiv(n+d-1); if(s[d]!=nd, break);\\Если успело разложиться проверим подходит ли (numdiv выполнится моментально)
   ,
      f=factor(n+d-1,2^15); fn=matsize(f)[1];\\Не успело, проверим все ли делители найдены, до 2^15 проверяется очень быстро

2.

if(type(alarm(1,factor(n+d-1)))!="t_ERROR", 
s[d]=numdiv(n+d-1); if(s[d]!=nd, break); 
, f=factor(n+d-1,2^15); fn=matsize(f)[1];

3.

if(type(alarm(1,d(n+d-1)))!="t_ERROR", 
s[d]=c(n+d-1); if(s[d]!=nd, break); 
, f=b(n+d-1,2^15); fn=a(f)[1];

И т.д. в обратную сторону.
При этом в п.3 сознательно заменены matsize на a, равно как factor на b, и т.д. в обратную сторону...
Это еще не доходя до уровня раскладывания самой строчки
if(type(alarm(1,factor(n+d-1)))!="t_ERROR"

Конечно было-бы легче пощупать некий предельно простой примитив кода, на котором валится сам factor() (это предположительно, для данного случая).

[Dmitriy40]

DemISdx в сообщении #1719378 писал(а):

Конечно было-бы легче пощупать некий предельно простой примитив кода, на котором валится сам factor() (это предположительно, для данного случая).

Factor() ни при чём, валится alarm(), при дальнейшем разбирательстве в теме про PARI даже нашли в каком месте исходника вероятна ошибка. Что можно было легко увидеть если пройти по сразу данной мною ссылке и почитать пяток сообщений с обсуждением далее ...

[Yadryara]

Dmitriy40 в сообщении #1719371 писал(а):

Впрочем, если так упростить задачу, то всё ещё проще:

for(d=0,#v-1, if(type(nd=alarm(1,numdiv(n+d))) <> "t_ERROR" && nd<>96, next(2)); );

Да, вот эту версию удалось адаптировать под мой код:

for(j = 1, #ostmes,

if(nu[ostmes[j]] - kolfac[ostmes[j]] == 1,

t3=getwalltime();

if( type(kfnew = alarm(1,factor( (n+ostmes[j]-1)/(v[ostmes[j]]*pro[ostmes[j]]) ) ))

 <> "t_ERROR" && #kfnew[2] > 2, print(">2   ",strtime(getwalltime()-t3));next(2),

print("   ",strtime(getwalltime()-t3));hormes = concat( hormes, ostmes[j]);
);

));

И в результате максимальное время проверки теперь 1,001 ms:

   1 ms
>2   123 ms
   52 ms
   1,001 ms
   1,001 ms
>2   102 ms
   1 ms
   1 ms
   30 ms
>2   62 ms
   1 ms
   4 ms
>2   75 ms
>2   12 ms
>2   106 ms
>2   29 ms
   1 ms
>2   69 ms
   1,001 ms
   656 ms

Убрал засечку и вывод времени, протестировал по известным результатам и вставил в рабочую программу.

Пусть поработает, сколько-то процентов времени, видимо, удастся выиграть. Дальше надо думать, как сделать порог ниже секунды.

[Yadryara]

Никакого заметного ускорения выявить не удалось, даже на 5%.

Кроме того, видимо и у меня возник глюк с алармом к исходу первого часа работы.

  *** init_Rab_243_20_2: Warning: increasing stack size to 16000000.
  *** init_Rab_243_20_2: Warning: increasing stack size to 32000000.
  *** init_Rab_243_20_2: Warning: increasing stack size to 64000000.
  *** init_Rab_243_20_2: Warning: increasing stack size to 128000000.
14   nedobor = 1
12   nedobor = 2
  *** init_Rab_243_20_2: Warning: increasing stack size to 256000000.
13   61127909902100883545614516363964389591223990348728555069474808              111 1 1111    1  1        10
16   nedobor = 3
10   nedobor = 4
12   nedobor = 5
  ***   at top-level: init_Rab_243_20_2()
  ***                 ^-------------------
  *** init_Rab_243_20_2: alarm interrupt after 1,000 ms.
  ***   Break loop: type 'break' to go back to GP prompt
break>

И в результате один поток остановился. После двукратного Cntrl+D прога в этом окне перезапустилась со следующего юнита.

[Yadryara]

Ну вот без сбоев проработали потоки достаточное время, теперь можно сравнить.

8 потоков, старая версия:

14.   0.          26min, 57,251 ms
89.   0.          26min, 30,424 ms
44.   0.          27min, 5,002 ms
74.   0.          26min, 47,907 ms
59.   0.          28min, 2,333 ms
29.   0.          28min, 7,194 ms
119.  0.          26min, 55,583 ms
104.  0.          26min, 38,369 ms

8 потоков, новая версия:

31.   0.          27min, 44,691 ms
16.   1.          27min, 50,117 ms
61.   1.          28min, 4,407 ms
76.   1.          27min, 30,654 ms
91.   1.          26min, 57,898 ms
46.   1.          27min, 51,822 ms
106.  1.          27min, 24,864 ms
2.    0.          28min, 34,948 ms

То есть не то что не быстрее, а где-то в среднем на 2% медленнее считались. Да, юниты не ровно те же самые, но они однтотипные.

[DemISdx]

Dmitriy40 в сообщении #1719380 писал(а):

Factor() ни при чём, валится alarm(), при дальнейшем разбирательстве в теме про PARI даже нашли в каком месте исходника вероятна ошибка. Что можно было легко увидеть если пройти по сразу данной мною ссылке и почитать пяток сообщений с обсуждением далее ...

О!
Извините.
Реально не дочитал дальше...
Интересно конечно.
aitap
Там высказал вполне интересное соображение.
Единственное, что мне пока не совсем понятно, воспроизводится ли подобное на чистом линуксе (не WSL).