Список всех простых чисел
Sergei M. Abramov abram_AT_botik.ru
Простите, что не про Рефал. Просто не могу вспомнить то, что (как мне
помнится) случилось в обстабовке рефал-компании на заре моей
профессиональной деятельности (практика в НИЦЭВТе или работа в нем).
По каким-то причинам я задумался тогда над построением списка всех
простых чисел. И в некий момент сам придумал алгоритм, которрый
опирался вот на какое определение последовательности простых чисел:
1. Очередное натуральное число n будет простым, если оно
взаимнопросто со всеми ранее найдеными простыми числами.
или (что то же самое, но уже совсем близко к коду)
2. Очередное натуральное число n -- простое, если gcd(pp,n)==1, где
pp -- произведение всех ранее найденых простых чисел.
Понятно, что нужна сверхдлинная целочисленная арифметика. Но,
gcd(pp,n) обещал посчитаться быстрее, чем деление n по очереди на все
сомножители из pp.
Я этот алгоритм (точнее подход, поскольку сам алгоритм похитрее, о нем
в конце) придумал сам и очень этим гордился. Но потом в
рефал-компании кто-то мне дал книжку. Серегей Романенко? Андрей или
Алик Климов? Колля Кондратьев? -- не помню.
Книжка была кого-то из великих... Дейскстра? Вирт? И, почему-то,
мне кажется, что называлась она "Записки из башни слоновой кости". А
в книжке, как мне сегодня кажется, были программистсткие побасенки.
Ну, и был коротенький шутливый этюд (близко к тексту):
================ Цитата по памяти ===============================
Для эффективной генерации всех простых чисел мнe достаточно всего двух
ячеек памяти:
integer pp, n;
n := 2;
pp := 1;
while ( True ) do begin
if ( gcd(pp, n)=1 ) then begin
println(n);
pp := pp * n;
end;
n := n + 1;
end;
==================================================================
С тех пор я узнал, что хотя я и самостоятельно нашел такой подход к
генерации простых, но я был не единственным (и, скорее всего, не
первым)
ПРОСЬБА О ПОМОЩИ:
1. Коллеги, если кто-то может мне дать имя автора, название книги,
которая мне помнится, то я буду благодарен.
2. И если есть еще какие-либо ссылки на такой подход к проверке
простоты -- gcd(pp, n)=1,-- то тоже полезно.
Почему вспомнил? Вспомнил я про эту историю во время очередной своей
лекции по Хаскеллю в ЯрГУ. Ну, и сегодня написал немного строк на
Хаскелле -- приложены. Взгляните, для интереса.
primesAD -- это генератор, прямо в лоб переписан приведенный выше
паскале-подобный текст.
А теперь поговорим про решето. Решето вообще (не обязательно
Эратосфена) -- это когда берем начальное пространство поиска --
например, [2..]. И дальше всегда делаем такую пару операций:
1. Некие первые элементы (один или больше) в текущем прастранстве
поиска -- заведомо простые числа -- выдаем их в результат.
2. Затем при помощи только что найденых простых фильтруем (сужаем)
пространство поиска.
Если брать только одно простое n из пространства поиска и фильтровать
предикатом ((/=0).(`mod`n)) -- то это решето Эратосфена.
А если брать много заведомо простых из пространства поиска, вычислять
их произведение pp и фильтровать предикатом ((==1).(gcd pp)) -- то это
решето ... назову решетом Абрамова ;-) Пока Вы мне не подскажете
другого автора такого же подхода (или где искать).
Три остальные программки в Хаскелле -- это разные варианты решета:
primesE, primesE1 -- решето Эратосфена (отличии в изобразительных
средствах).
primesА -- решето Абрамова.
А заканчивается файл статистикой (время счета и память) прогона
некоторых тестов в HUGSе и в GHCi.
Посмотрите статистику... Шуточный этюд "достаточно две ячейки памяти"
имеет нешуточное приложение, как мне кажется ;-)
Я уж молчу про то, что gcd по Евклиду (как оно описано в прелюдии
хаскелля) не самая эффективная реализация gcd (по Кнуту).
Всего доброго,
Сергей Абрамов
Александр Коновалов a.v.konovalov87_AT_mail.ru
Добрый день всем!
Нужно спасать рассылку от офтопа.
Переписал исходник Сергея Михайловича на Рефал-5. Функции primesE, primesAD и primesA переписаны максимально близко к тексту, ради этого пришлось ввести функции filter, o (композиция функций), bind-right (чтобы связать Mod со вторым аргументом), eq и neq.
Рефал не поддерживает ленивые вычисления и бесконечные списки, поэтому все функции на входе принимают последовательность чисел от 2 до 40 000.
Понятно, что вызов <filter (o ((neq 0) (bind-right (Mod s.n)))) e.ns> на Рефале будет выполняться медленно. Поэтому я также написал оптимизированные версии primesE-opt и primesA-opt, в которых фильтрация выполняется не функцией высшего порядка, а специально написанной функцией. Функция primesE в этом случае ускорилась на порядок, primesA — незначительно.
Замеры времени есть в исходнике. Масштаб величин примерно такой же: primesE на порядок медленнее primesAD, primesAD на порядок медленнее primesA. Функции primesE-opt и primesAD требуют одинакового времени работы, но первая делает гораздо больше шагов. Причина этого расхождения, по-видимому, в длинной арифметике, в вычислении gcd(pp, n), поскольку pp — длинное число.
С уважением,
Александр Коновалов
-----Original Message-----
From: Sergei M. Abramov abram_AT_botik.ru <re...@botik.ru>
Sent: Tuesday, October 8, 2019 10:29 PM
To: re...@botik.ru
Subject: Список всех простых чисел
День добрый, всем!
Простите, что не про Рефал. Просто не могу вспомнить то, что (как мне помнится) случилось в обстановке рефал-компании на заре моей профессиональной деятельности (практика в НИЦЭВТе или работа в нем).
Arkady Klimov arkady.klimov_AT_gmail.com
Sergei M. Abramov abram_AT_botik.ru
> То есть для последнего "поколения" эта работа лишняя.
чисел?
> Наверно, лучше накапливать в виде списка, а при подаче на Filter
> для следующего шага все перемножить.
> Интересно, для Хаскела это тоже актуально, или там умножение pp*n ленивое?
> А еще хотелось бы посмотреть статистику по числу (доле) отсеянных
> (прошедших) через каждый из фильтров (это мне уже как математику
> интересно).
есть ли хоть какие-то работы (а если есть, то какие) в этой отрасли,
где используется (\ k -> (gcd pp k) == 1)?
> Нужно спасать рассылку от офтопа.
> Понятно, что вызов <filter (o ((neq 0) (bind-right (Mod s.n))))
> e.ns> на Рефале будет выполняться медленно.
ns' = filter (\ k -> (gcd pp k) == 1) ns
работает немного медленнее, чем
ns' = filter ((== 1).(gcd pp)) ns
Вот нифига ж себе?
Всего доброго,
Сергей Абрамов
Andrei Klimov andrei_AT_klimov.net
Вот, хочу сказать, что меня удивило:
ns' = filter (\ k -> (gcd pp k) == 1) ns
работает немного медленнее, чем
ns' = filter ((== 1).(gcd pp)) ns
Вот нифига ж себе?
Sergei M. Abramov abram_AT_botik.ru
Замечу, что если взять несколько первых простых чисел p1...pK, их
произведение pp и профильтровать (\ k -> (gcd pp k) == 1) весь
бесконечный список [(pk+1).. ], то результат фильтрации имеет очень
простое, эффективное, конструктивное представление на Haskell-е. Q
это в своих этюдах использовал.
То есть, несколько поколеий можно прогнать в статике, на этапе
написания программы.
>> А еще хотелось бы посмотреть статистику по числу (доле) отсеянных
>> (прошедших) через каждый из фильтров (это мне уже как математику
>> интересно).
знакомое с детсва 6n-1, 6n+1:
[ n+r | n <- [6, 12..], r <- [-1, 1]]
Если взять 2, 3, 5 (pp=30), то:
[ n+r | n<-[30, 60 ..], r<-[-23,-19,-17,-13,-11,-7,-1,1]]
Если взять 2, 3, 5, 7 (pp=210), то:
[ n+r | n <- [210, 420 ..], r <- bs]
bs = [-199,-197,-193,-191,-187,-181,-179,-173,-169,-167,
-163,-157,-151,-149,-143,-139,-137,-131,-127,-121,
-113,-109,-107,-103,-101,-97,-89,-83,-79,-73,-71,-67,
-61,-59,-53,-47,-43,-41,-37,-31,-29,-23,-19,-17,-13,
-11,-1,1] :: [Integer]
Сжатие серч-спэйса (доля неотсортированного), о чем спросил Алик,
состaвит для этих трех случаев: 2/6, 8/30, 48/210
Программка, которая по заданному начальному списку [2, 3,..pK] строит
bs, для выражения:
[ n+r | n <- [pp, 2*pp ..], r <- bs]
очень простая.
Всего доброго,
Сергей Абрамов
Arkady Klimov arkady.klimov_AT_gmail.com
День добрый, всем!
> То есть для последнего "поколения" эта работа лишняя.
А что такое последнее поколение, если я строю список всех простых
чисел?
> Наверно, лучше накапливать в виде списка, а при подаче на Filter
> для следующего шага все перемножить.
> Интересно, для Хаскела это тоже актуально, или там умножение pp*n ленивое?
Ленивое, несомненно.
> А еще хотелось бы посмотреть статистику по числу (доле) отсеянных
> (прошедших) через каждый из фильтров (это мне уже как математику
> интересно).
Алик, а мне интересно от математиков (может знакомые есть?) узнать,
есть ли хоть какие-то работы (а если есть, то какие) в этой отрасли,
где используется (\ k -> (gcd pp k) == 1)?
> Нужно спасать рассылку от офтопа.
Спамера забанить!
> Понятно, что вызов <filter (o ((neq 0) (bind-right (Mod s.n))))
> e.ns> на Рефале будет выполняться медленно.
Вот, хочу сказать, что меня удивило:
ns' = filter (\ k -> (gcd pp k) == 1) ns
работает немного медленнее, чем
ns' = filter ((== 1).(gcd pp)) ns
Вот нифига ж себе?
Всего доброго,
Сергей Абрамов
Sergei M. Abramov abram_AT_botik.ru
> следующая граница (lim) получается из предыдущей примерно
> возведением в квадрат, так ведь?
> Поколение - то, что выходит после очередной фильтрации как кусок
> списка простых и служит сомножителями в следующем pp. Меня
> интересует сколько (какую долю) отсеивает каждая следующая
> фильтрация. Понятно, что сначала это 1/2, далее 1-(2/3)*(4/5),
> далее 1-(6/7)*(10/11)*...*(22/23) и т.д. Произведение - это доля
> остающихся. Просто хотел бы увидеть десятичные числа, как они
> меняются, наверно, убывают к 0. Насколько быстро?
Хотя и не напрямую.
>> Интересно, для Хаскела это тоже актуально, или там умножение pp*n ленивое?
>
> Ленивое, несомненно.
> Неужели?
> Значит на Хаскеле будет автоматически накапливаться в форме терма 7*11*13*..., да?
прохождение через lim и начнется первый gcd с накопленным, но не
посчитанным, произведением. Вот там случится need результата.
Все жадные операции перечислены в Haskell-е. И я не вижу среди них
умножения, тем более Integer.
Всего доброго,
Сергей Абрамов
Anton Orlov orlovan_AT_gmail.com
where sieve (n:ns) lim ~ps@(p1:p2:ps3)
| n < lim = n : sieve ns lim ps
| otherwise = sieve (filter ((/=0).(`mod`p1)) ns) (p2*p2) (p2:ps3)
176087
(0.16 secs, 168,810,848 bytes)
*Main> primesE2!!16000
176087
(0.33 secs, 169,225,768 bytes)
*Main> primes!!16000
176087
(0.03 secs, 54,413,416 bytes)
| l1 < r1 = l1:(minus ls r)
| l1 > r1 = minus l rs
| otherwise = minus ls rs
unionAll ((l1:ls):xs) = l1 : (unionAll $ f ls xs)
where f l@(l1:ls) ((r@(r1:rs)):xs)
| l1 < r1 = l:r:xs
| otherwise = r:(f l xs)
176087
(0.75 secs, 642,165,016 bytes)
Sergei M. Abramov abram_AT_botik.ru
Огромное спасибо за письмо и интересные материалы.
> 1) primesA сильно быстрее, чем primesE, за счёт прыжков до
> следующего квадрата.
> 2) Решето Эратосфена -- оно же не про делимость, а про вычёркивание
> с определённым шагом. Эффективная запись этого на Хаскеле выглядит так:
Особенно то, что в коде нет ни одного деления (ни mod, но gcd)
> Я подсмотрел определение primes выше вот на этой страничке:
> https://wiki.haskell.org/Prime_numbers Там эффективная реализация
> решета выводится последовательно, и можно понять, что происходит.
упоминается, а вот gcd -- нет.
> Эта unionAll -- довольно хитрая штука. С наивной реализацией
> (0.75 secs, 642,165,016 bytes)
Тоже не плохо.
Всего доброго,
Сергей Абрамов
Александр Коновалов a.v.konovalov87_AT_mail.ru
Добрый день, Аркадий!
Прошу прощения за некропостинг, но на письмо отвечу.
Я остановился на 40 000, поскольку решето Эратосфена на Рефале оказалось сильно медленным, и его время растёт нелинейно от длины исходного списка. Причём, мне показалось, что даже квадратично. Поэтому результата в 176 100 я просто не дождусь.
Время я замерял для переславской реализации Рефала-5. Рефал-5λ сильно медленнее (в разы), поэтому мне стыдно было показывать эти числа 😳.
В последующей переписке Антон Орлов продолжил оффтопик, мне уже лень пытаться переписывать его примеры на Рефал.
С уважением,
Александр Коновалов
Sergei M. Abramov abram_AT_botik.ru
> Прошу прощения за некропостинг, но на письмо отвечу.
> Я остановился на 40 000, поскольку решето Эратосфена на Рефале
> оказалось сильно медленным, и его время растёт нелинейно от длины
> исходного списка. Причём, мне показалось, что даже квадратично.
> Поэтому результата в 176 100 я просто не дождусь.
> Время я замерял для переславской реализации Рефала-5. Рефал-5λ
> сильно медленнее (в разы), поэтому мне стыдно было показывать эти числа .
> В последующей переписке Антон Орлов продолжил оффтопик, мне уже лень
> пытаться переписывать его примеры на Рефал.
То, что нарыл Антон я ввел в свой курс по FP (Хаскелль) под заголовком
"Эратосфен это про вычеркивание, а не про делимость".
Действительно, Эратосфен говорил: оставь P как очередное простое и иди
с шагом P дальше и вычеркивай. То есть, он не про делимость, а про
удаление из изходного множества арифметической прогрессии:
[2P, 3P ..]
С этой точки зрения, чистый Эротосфен выглядит так (это все можно
описать на Рефале, но у меня нет Рефала):
primesE1 :: [Integer]
primesE1 = sieve [2..]
where sieve (p:xs) = let c1 = 2*p
c2 = c1 + p
in p : sieve (minusOS xs [c1, c2 ..])
> primesE1!!2000
> 17389
> (14831010 reductions, 27298448 cells, 2 garbage collections)
здесь minusOS -- вычитание "упорядоченных множеств", представленных
отсортированными ленивыми бесконечными списками:
minusOS :: [Integer] -> [Integer] -> [Integer]
minusOS (x:xs) (y:ys)
| x < y = x : minusOS xs (y:ys)
| x == y = minusOS xs ys
| x > y = x : minusOS xs ys
То, что нашел Антон, это следующих шаг: давайте вычеркивать будем не
как сказал Эратосфен, а вычеркнем ОБъЕДИНЕНИЕ всего, что советовал
вычеркивать Эратосфен. Появляется функция объединения бесконечного
множетсва "упорядоченных множеств", причем они перечислены в порядке
роста первого элемента:
primesE :: [Integer]
primesE = 2 : minusOS [3 ..] cs
where
cs = unionOS
[ [c1, c2 ..] | p <- primesE,
let c1 = 2*p,
let c2 = c1+p
]
unionOS :: [[Integer]] -> [Integer]
unionOS ((x:xs) : xss) = x : unionOS (ins xs xss)
where ins (x:xs) ((y:ys) : yss)
| x < y = (x:xs) : (y:ys) : yss
| x > y = (y:ys) : ins (x:xs) yss
| x == y = (y:ys) : ins xs yss
и этот вариант не лучше, а хуже (в 3+ раз) обычного Эратосфена:
> primesE!!2000
> 17393
> (47292483 reductions, 83601402 cells, 8 garbage collections)
А вот дальше можно сообразить, что первое число, которое удастся
простому P действительно вычеркнуть своей арифметической прогрессией
[2P, 3P..] будет P*P -- все меньшие члены прогрессии уже вычеркнуты до
него. Поэтому прогрессию можно заменить на [P*P, P*P+P ..]
В коде выше это замена всех c1 = 2*p на c1 = p*P
И это ускоряет обычный Эратосфен совсем немного, а Эратосфен с
объединением unionOS -- радикально и он на порядок (раз в 6) обгоняет
обычного!
> primesE1!!2000
> 17387
> (14692824 reductions, 27114235 cells, 2 garbage collections)
> primesE!!2000
> 17393
> (2731729 reductions, 4790869 cells)
Тупое прозрачное опредленеие "N -- простое, если не делется ни на одно
простое P, которое P*P<=K" дает лишь втрое хуже результат:
primes1 = 2 : filter (isPrime1) [3..]
where isPrime1 n = null(filter (\ p -> n`mod`p == 0)
(takeWhile (\ k -> (k*k <= n)) primes1))
> primes1!!2000
> 17393
> (7394774 reductions, 10355414 cells, 1 garbage collection)
И, наконец, мное описанное решето Абрамова-Дейкстры(?), которое не про
вычеркивание, не про делимость, а про взаимную простоту с
произведением простых (gcd pp n)==1 дает такой результат:
primesA = sieve 1 1 2 [2..]
where sieve pp n' lim ns@(n:ns1)
| n <= lim = n : sieve (pp*n) n lim ns1
| otherwise = sieve 1 1 (n'*n') (filter ((== 1).(gcd pp)) ns)
> primesA!!2000
> 17393
> (4254611 reductions, 7465819 cells)
что хуже ("всего" в 1.55 раз) Эратосфена вычеркивания с объединением
unionOS.
Это все.
Есть еще одна идея ускорения, но она применимы ко всем. И если
исследовать ее эффект, то применяя ее ко всем. В коде, который Антон
Орлов нашек на Haskell.org эта идея применена для No=1.
Идея состоит в том, что можно "в уме" вычеркнуть из [2..] все числа,
которые делятся на 2, 3, 5, ... p (взято No первых простых чисел).
Пусть ns результат вычеркивания. Тогда он определяется вот так:
Пусть pp произведение 2, 3, 5, ... p. Тогда:
ns = [ n+r | n <- [pp, 2*pp ...], r <- rs ]
where rs = filter ((==1).(gcd pp)) [p+1-pp .. p]
Если взять No простых, то после фильтрации в rs останется не pp
элементов, a меньше -- length rs. И от исходного множества у нас
останется вот такая доля: (length rs) / pp. Как-то так:
No length rs pp Доля
1 1 2 50.0%
2 2 6 33.3%
3 8 30 26.7%
4 48 210 22.9%
5 480 2310 20.8%
6 5760 30030 19.2%
7 92160 510510 18.1%
На мой вкус, разумно внедрять случай No=4 или No=5. Ускорение должно
быть вполне себе.
И, последнее, жалко, что не имею рефала. Можно все упомянутое
написать на рефале. Даже в условиях отсутствия ленивости, ее
эммуляция в этих алгоритмах на сложна.
С уважением,
Абрамов С.М.
ab...@botik.ru
мобильный: +7(903)2928308
Boyko Bantchev boykobb_AT_gmail.com
наследника APL, считает, что разницу в скорости меньше чем в два раза
не следует считать разностью. И действительно, она вполне может быть
вызвана не самими алгоритмами, а другими обстоятельствами.
Александр Коновалов a.v.konovalov87_AT_mail.ru
Добрый день всем!
Нет, я когда-нибудь напишу конвертор Рефала/2 в Рефал-5. Это будет быстрее, чем научиться бегло читать эту шифровку.
С уважением,
Александр Коновалов
From: Василий Стеллецкий swi_AT_cnshb.ru [mailto:re...@botik.ru]
Sent: Sunday, November 10, 2019 1:16 PM
To: re...@botik.ru
Subject: Re: Список всех простых чисел
Добрый день!
Странно! Попробовал реализовать "решето Эратосфена", правда, свой алгорим. Для поля 176100 около 7 сек.
Все потроха приложил в архиве...
Основная часть алгоритма (сложение и умножение чисел в символьном виде не привел, можно посмотреть в прилож.архиве.):
SWAP ^n
a =/0/k/^n/k/bi/k/GI/...k/erit/()k/i/('0')(k/RDR//^n/.)..
bi e1' 'ee=e1
i w1w1=
w1w2='1'k/i/(k/add/w1'1'.)w2. * забиваем массив единичками
erit ()s1s2ee=k/f0/(s1s2)ee.
(e0)e1'1'ee=k/f0/(e0e1'1')ee.
f0 (e1)ee=k/fc/(k/LENGW/e1.)(k/LENGW/k/RDR//^n/..)ee.
fc (sne1)ee=k/f1/(e1)(k/LENGW/k/mul/(k/CVDN/sn.)k/CVDN/sn...)ee. *часть проверки, что квадрат простого числа не больше размера массива
f1 w0(s1e2)(s1e4)ee=k/f2/k/cmp/(e2)e4.w0ee.
w0(s1e2)(s3e4)ee=k/f2/k/cmp/(s1)s3.w0ee.
f2 '>'(e0)ee=k/e/('1')e0ee.
ssw0ee=k/fe/k/f/()w0ee..
fe ee(e0)(e1)=k/erit/(e0e1)ee.
f (e0)(s1)s2ee='0'k/f/()(e0s1)ee. * заменяем дырку на нуль
(e0)(s1e2)ssee=ssk/f/(e0s1)(e2)ee.
(e0)(e1)=(e0)(e1)
e w0'0'ee=k/e/(k/add/w0'1'.)ee.
(e0)'1'ee=k/P/e0.k/e/(k/add/(e0)'1'.)ee. *вывод полученного списка
w0=
--
С уважением,
--
Василий Стеллецкий
10.11.2019, 00:50, "Александр Коновалов a.v.konovalov87_AT_mail.ru" <re...@botik.ru>:
Александр Коновалов a.v.konovalov87_AT_mail.ru
Речь идёт о сравнении различных алгоритмов, скажем пузырьковая сортировка vs быстрая сортировка? Или о чём речь?
Когда я занимаюсь оптимизациями в своём компиляторе, я замеряю и пятипроцентное ускорение. По нескольким запускам и с учётом доверительного интервала.
-----Original Message-----
From: Boyko Bantchev boykobb_AT_gmail.com [mailto:re...@botik.ru]
Sent: Sunday, November 10, 2019 2:16 PM
To: re...@botik.ru
Subject: Re: Список всех простых чисел
Boyko Bantchev boykobb_AT_gmail.com
> Речь идёт о сравнении различных алгоритмов, скажем пузырьковая сортировка vs быстрая сортировка? Или о чём речь?
> Когда я занимаюсь оптимизациями в своём компиляторе, я замеряю и пятипроцентное ускорение. По нескольким запускам и с учётом доверительного интервала.
незначительного ускорения. Что считать незначительным, полагаю, нужно
решать по контексту.
Нужно иметь ввиду и то, что ускорение может быть только кажущимся –
вдруг оно исчезнет или даже окажется замедлением при переходе к
другому компилятору того же языка, к другой версии тоже же самого
компилятора или к машине с другими параметрами.
Василий Стеллецкий swi_AT_cnshb.ru
> Нет, я когда-нибудь напишу конвертор Рефала/2 в Рефал-5. Это будет быстрее, чем научиться бегло читать эту шифровку.
Шифровка на мой взгляд еще круче:
С уважением,
Александр Коновалов a.v.konovalov87_AT_mail.ru
Добрый вечер, Василий!
Если расставить отступы в более привычном для Рефала-5 стиле, то в программе разобраться проще.
А если переименовать некоторые встроенные функции, заменить ящик на копилку и добавить недостающие add, mul и cmp, работающие со строками цифр (литер цифр), то программа откомпилируется и даже заработает.
$ENTRY Go {
= <Br n '=' <bi <Card>>> <erit () <i ('0') (<Cp n>)>>;
}
bi {
e.1 ' ' e.e = e.1;
e.1 = e.1;
}
* забиваем массив единичками
i {
t.1 t.1 = ;
t.1 t.2 = '1' <i ( <add t.1 '1'>) t.2>;
}
erit {
() s.1 s.2 e.e = <f0 (s.1 s.2) e.e>;
(e.0) e.1 '1' e.e = <f0 (e.0 e.1 '1') e.e>;
}
f0 { (e.1) e.e = <fc (<Lenw e.1>)(<Lenw <Cp n>>) e.e> }
fc {
* часть проверки, что квадрат простого числа не больше размера массива
(s.n e.1) e.e = <f1 (e.1) (<Lenw <mul (<Symb s.n>) <Symb s.n>>>) e.e>;
}
f1 {
t.0 (s.1 e.2) (s.1 e.4) e.e = <f2 <cmp (e.2) e.4> t.0 e.e>;
t.0 (s.1 e.2) (s.3 e.4) e.e = <f2 <cmp (s.1) s.3> t.0 e.e>;
}
f2 {
'>' (e.0) e.e = <e ('1') e.0 e.e >;
s.s t.0 e.e = <fe <f () t.0 e.e >>;
}
fe { e.e (e.0) (e.1) = <erit (e.0 e.1) e.e> }
f {
(e.0) (s.1) s.2 e.e ='0' <f () (e.0 s.1) e.e>;
* заменяем дырку на нуль
(e.0) (s.1 e.2) s.s e.e = s.s <f (e.0 s.1) (e.2) e.e>;
(e.0) (e.1) = (e.0) (e.1);
}
e {
t.0 '0' e.e = <e ( <add t.0 '1'>) e.e >;
(e.0) '1' e.e = <Prout e.0> <e (<add (e.0) '1'>) e.e>;
* вывод полученного списка
t.0 = ;
}
cmp {
(e.L) e.R
, <Compare (<Numb e.L>) <Numb e.R>>
: {
'+' = '>';
'0' = '=';
'-' = '<';
}
}
mul { (e.L) e.R = <Symb <* (<Numb e.L>) <Numb e.R>>> }
add { (e.L) e.R = <Symb <+ (<Numb e.L>) <Numb e.R>>> }
Замер для классического Рефала-5:
D:\Mazdaywik\Documents\MEGA\primes>refc primes-swi.ref
Refal-5 Compiler. Version PZ Oct 29 2004
Copyright: Refal Systems Inc.
D:\Mazdaywik\Documents\MEGA\primes>echo 176400 | refgo -a primes-swi.rsl > primes.txt
Refal-5 System. Version PZ Oct 29 2004
*** The View Field:
<STOP$ >
*** Number of Steps = 16408006
Elapsed system time: 14.406 seconds
*** Buried:
((n '=176400' ))
Memory allocated = 2826240 Bytes
Замер для Рефала-5λ:
D:\Mazdaywik\Documents\MEGA\primes>srefc primes-swi.ref
*Compiling primes-swi.ref:
** Compilation succeeded **
D:\Mazdaywik\Documents\MEGA\primes>echo 176400 | primes-swi.exe > primes.txt
Total program time: 36.00000 seconds (100.0 %).
(Total refal time): 29.62200 seconds (82.3 %).
Linear result time: 18.31600 seconds (50.9 %).
Linear pattern time: 10.88800 seconds (30.2 %).
Runtime time overhead: 3.92800 seconds (10.9 %).
Native time: 2.45000 seconds (6.8 %).
t- and e-var copy time: 0.40200 seconds (1.1 %).
Repeated t-var match time (outside e-loops): 0.01600 seconds (0.0 %).
Step count 53533728
Identifiers allocated: 89
Memory used 177000 nodes, 177000 * 16 = 2832000 bytes
Рефал-5 оказался медленнее Рефала/2, как мне кажется, из-за конвертирования чисел в строки и наоборот в функциях add и mul. Возможно, если переписать на обычную длинную арифметику, будет быстрее.
А у меня стыдоба: медленнее более чем в два раза (36/14,4 = 2,5). Причины: кривой промежуточный код, длинная арифметика написана пополам на Рефале и Си++, в то время как в Рефале-5 вся арифметика написана на Си. Максимальная оптимизация (-OdDPRS) даёт лишь 21,625 секунд, что тоже медленнее Рефала-5. В общем, у меня ещё много работы.
Хотел сделать замер для Рефала/2 на своей машине, но не получилось. Версия на сайте Василия вылетает для чисел, больших 32665. На сайте лежит исполнимый файл refalb.exe, но в архиве, который выкладывал Василий, упоминается refald. Видимо, эта какая-то более свежая версия, ещё не выложенная на сайте.
Мы с Василием спасаем ветку от офтопика!
Александр Коновалов a.v.konovalov87_AT_mail.ru
Василий!
У меня в папке лежат test.bat (см. ниже), REFAL.CM (скачан с сайта), REFALB.EXE (скачан с сайта), Refald.exe (из предыдущего письма). Батник имеет следующее содержание:
@setlocal
@cls
set PROMPT=$P $T$G
REFALB.EXE REFAL.CM erit.ref erit.cm
echo %1 > input.txt
REFALB.EXE erit.cm < input.txt > primes.txt
rem OK
@endlocal
Когда я в нём использую REFALB.EXE, всё работает как надо, но до 32 665. Если я заменяю REFALB.EXE на Refald.exe, то ничего не работает, даже не создаётся erit.cm. Похоже, что лежащий на сайте REFAL.CM не совместим с новой версией интерпретатора.
Написал такой сценарий, который использует erit.cm из Вашего письма:
@setlocal
@cls
set PROMPT=$P $T$G
echo %1 > input.txt
Refald.exe erit.cm < input.txt > primes.txt
rem OK
@endlocal
Вывод:
D:\Mazdaywik\Documents\MEGA\primes\swi>set PROMPT=$P $T$G
D:\Mazdaywik\Documents\MEGA\primes\swi 21:01:54,81>echo 176400 1>input.txt
D:\Mazdaywik\Documents\MEGA\primes\swi 21:01:54,81>Refald.exe erit.cm 0<input.txt 1>primes.txt
D:\Mazdaywik\Documents\MEGA\primes\swi 21:02:03,98>rem OK
Отработало за 9,17 секунд. Другие замеры дали 9,36, 9,33, 9,24 секунды. Так что у меня не около 7, а около 9 секунд получается. Но всё равно быстрее, чем в классическом Рефале-5, и существенно быстрее, чем в Рефале-5λ.
С уважением,
Александр Коновалов
From: Василий Стеллецкий swi_AT_cnshb.ru [mailto:re...@botik.ru]
Sent: Sunday, November 10, 2019 8:30 PM
To: re...@botik.ru
Subject: Re: Список всех простых чисел
Добрый вечер, Александр!
в приложенном архиве два файла с двойным расширением (иначе не проходит) второе надо убрать.
Думаю, разберетесь!
Успехов!
--
С уважением,
--
Василий Стеллецкий
10.11.2019, 20:05, "Александр Коновалов a.v.konovalov87_AT_mail.ru" <re...@botik.ru>:
Василий Стеллецкий swi_AT_cnshb.ru
Рад, что файл получили, и у Вас получилось (с доставкой на gmail опять проблемы как не крутись :( ). Мне казалось, что просили 176100, а не 176400, но разница, наверное, небольшая. Видимо разница 7 и 9 сек из-за скорости компа ;)
А вообще результаты интересные.
На каком же алгоритме Вы не могли дождаться результата? Наверное, проверяли все полученные простые? а не до квадратного корня?
С уважением,
Александр Коновалов a.v.konovalov87_AT_mail.ru
Василий!
Я не мог дождаться результата на программе, переписанной втупую из Хаскеля:
primesE {
s.n e.ns = s.n <primesE <filter (o ((neq 0) (bind-right (Mod s.n)))) e.ns>>;
/* empty */ = /* empty */;
}
Оригинал на Хаскеле:
primesE = sieve [2..]
where sieve :: [Integer] -> [Integer]
sieve (n:ns) = n : sieve (filter ((/=0).(`mod`n)) ns)
Функции высшего порядка в Рефале-5 — не самая быстрая вещь. Их вообще нет в Рефале-5 как таковых. Вместо них есть вызов по имени при помощи «метафункции» Mu, через неё эмулируется косвенный вызов «по указателю». Понятно, что эмулирование этим путём каррирования, композиции и прочей фильтрации получается довольно медленной.
В первом письме ветки Абрамов написал несколько функций вычисления последовательности простых чисел, при этом primesE была самой медленной. Я их всех переписал на Рефал, на сколько это возможно, и тоже замерил производительность. В детали алгоритмов я вообще не вдавался — где вычисляется квадрат, а где нет.
«Приятно, что моя реализация даёт неплохие результаты...»
Да, неплохие. Но есть реализация и побыстрее. У меня есть компилятор Рефал-05, компилирующий подмножество Рефала-5 в Си. На нём тот же пример вычислился быстрее:
D:\Mazdaywik\Documents\MEGA\primes>set R05CCOMP=bcc32 -w -DR05_SHOW_STAT=1
D:\Mazdaywik\Documents\MEGA\primes>refal05c primes-swi.ref Library refal05rts
*Compiling primes-swi.ref:
+Linking D:\Mazdaywik\Documents\Refal-05\lib/Library.c
+Linking D:\Mazdaywik\Documents\Refal-05\lib/refal05rts.c
Borland C++ 5.5.1 for Win32 Copyright (c) 1993, 2000 Borland
primes-swi.c:
d:\mazdaywik\documents\refal-05\lib/library.c:
d:\mazdaywik\documents\refal-05\lib/refal05rts.c:
Turbo Incremental Link 5.00 Copyright (c) 1997, 2000 Borland
*** Compilation successed ***
Builtin time: 2.594 seconds (100.0 %).
Total program time: 2.594 seconds (100.0 %).
Step count 38622
Memory used 34889 nodes, 34889 * 16 = 558224 bytes
D:\Mazdaywik\Documents\MEGA\primes>echo 176400 | primes-swi.exe > primes.txt
Total program time: 6.422 seconds (100.0 %).
(Total refal time): 5.937 seconds (92.4 %).
Linear result time: 3.654 seconds (56.9 %).
Linear pattern time: 2.283 seconds (35.5 %).
Builtin time: 0.485 seconds (7.6 %).
Step count 16408006
Memory used 176704 nodes, 176704 * 16 = 2827264 bytes
Он не поддерживает длинную арифметику, но для этого теста вполне достаточно 32-разрядных целых чисел. Наверное, этот прирост как раз из-за отсутствия длинной арифметики.
«Топунов Владимир Леонидович нам-студентам говорил, что если ожидаемый выигрыш от изменения (оптимизации) алгоритма менее 10%, то за эту работу вообще не стоит браться…»
Проблема в том, что ожидаемый выигрыш порой трудно определить умозрительно, поэтому можно ошибиться в обе стороны. Или сделать оптимизацию, которая не даст ожидаемых 10 %, или наоборот отказаться от оптимизации, которая может дать значительный прирост.
С уважением,
Александр Коновалов
Александр Коновалов a.v.konovalov87_AT_mail.ru
Доброе утро, Василий!
Во-первых, я улучшил код конвертора Рефал/2→Рефал-5. Полученный результат по-прежнему сразу не совместим с Рефалом-5, но ручных правок требуется меньше. В частности, правок по стилю. Например, если раньше k/f/e1. конвертировалось в < f e.1 > (с пробелами вокруг скобок), то сейчас избыточные пробелы удаляются: <f e.1>. А в классическом Рефале-5 пробелы после < запрещены.
Во-вторых, перенёс обновлённый вариант на Рефал-5.
Замеры производтельности (для 176400):
Рефал/2: 1,12; 1,16; 1,17; 1,18; 1,23. Медиана — 1,17 секунд.
Рефал-5 PZ: 1,343; 1,344; 1,359; 1,360; 1,375. Медиана — 1,359 секунд.
Рефал-5λ (без оптимизаций): 21,234; 21,234; 21,282; 21,265; 21,407. Медиана — 21,282. (Два раза 21,234 — не ошибка, поскольку время квантуется.)
Рефал-5λ (максимальная оптимизация): 20,313; 20,328; 20,343; 20,438; 20,359. Медиана — 20,343. Оптимизация почти ничего не дала.
Рефал-05: 1,047; 1,047; 1,063; 1,063; 1,078. Медиана — 1,063. Числа повторяются из-за шага квантования 1/64 секунды ≈ 0,015–0,016 секунд.
Сравнивать Рефал/2 с другими реализациями напрямую некорректно, поскольку функции add и mul в них реализованы сильно по-разному. Но время выполнения сопоставимо с временем Рефала-5 PZ.
Почему Рефал-5λ показывает такие безобразные результаты — буду разбираться. Профилировка показывает, что более 90 % времени уходит на функции стандартной библиотеки (First, арифметика), а они написаны на смеси Рефала и Си++.
Рефал-05 опережает другие реализации только из-за того, что он компилируется сразу в Си, а не в интерпретируемый код.
С уважением,
Александр Коновалов
From: Василий Стеллецкий swi_AT_cnshb.ru [mailto:re...@botik.ru]
Sent: Monday, November 11, 2019 2:37 AM
To: refal@botik.ru
Subject: Re: Список всех простых чисел
Александр!
Вот как раз подумал про оптимизацию...
я функцию f писал в расчете на строки любой длины. В данном же случае можно воспользоваться встроенной функцией FIRST в которой длина задается одним символом.
Получилась замена:
*f (e0)(s1)s2ee='0'k/f/()(e0s1)ee.
* (e0)(s1e2)ssee=ssk/f/(e0s1)(e2)ee.
* (e0)(e1)=(e0)(e1)
f snw0(e1s2)ee=e1'0'k/f/snw0k/FIRST/snee..
sn w0'*'ee=ee()w0 * w0 протаскиваю для совместимости с уже написанной системой управления
в синтаксисе рефал-5 примерно так:
f { s.n t.0( e.1 s.2) e.e= e.1'0' < f s.n t.0 < FIRST s.n e.e > >;
s.n t.0'*' e.e= e.e() t.0;
}
Результат:
C:\temp>set PROMPT=$P $T$G
C:\temp 1:52:59,98>echo 176400 1>input.txt
C:\temp 1:53:00,00>Refald.exe eritf.cm 0<input.txt 1>primes.txt
C:\temp 1:53:01,18>rem OK
т.е. 1,18 сек :) ;)
Значит, основное время работы всё же давала не длинная целочисленная арифметика...
Приложил батник, текст и скомпилированный вариант, а также refald.exe, для не получивших его ранее (кое-где с двойным расширением, которое следует убрать ;) )
--
С уважением,
--
Василий Стеллецкий
10.11.2019, 22:28, "Александр Коновалов a.v.konovalov87_AT_mail.ru" <re...@botik.ru>: