Re: Список всех простых чисел

[Anton Orlov orlovan_AT_gmail . com]

Добрый день!

Прошу прощения, продолжу оффтоп про Хаскель.

Два небольших соображения:

1) primesA сильно быстрее, чем primesE, за счёт прыжков до следующего
квадрата. Если это реализовать для primesE, то разница в производительности
между ними уже не такая драматическая:

primesE2 = sieve [2..] 4 primesE2
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)

(Пришлось поизворачиваться, чтобы раскрутиться)

*Main> primesA!!16000
176087
(0.16 secs, 168,810,848 bytes)
*Main> primesE2!!16000
176087
(0.33 secs, 169,225,768 bytes)

2) Решето Эратосфена -- оно же не про делимость, а про вычёркивание с
определённым шагом. Эффективная запись этого на Хаскеле выглядит так:

primes = 2 : 3 : Data.List.Ordered.minus [5,7..]
(Data.List.Ordered.unionAll [[p*p, p*p+2*p..] | p <- tail primes])

*Main> primes!!16000
176087
(0.03 secs, 54,413,416 bytes)

minus -- ленивая разница возрастающих списков
unionAll -- ленивое объединение возрастающих списков, которые между собой
отсортированны по возрастанию первых элементов.

Эта unionAll -- довольно хитрая штука. Я подсмотрел определение primes выше
вот на этой страничке: https://wiki.haskell.org/Prime_numbers
Там эффективная реализация решета выводится последовательно, и можно
понять, что происходит.

С наивной реализацией unionAll всё, конечно, гораздо хуже:

minus l@(l1:ls) r@(r1:rs)
  | 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)

primes = 2 : 3 : minus [5,7..] (unionAll [[p*p, p*p+2*p..] | p <- tail
primes])

*Main> primes!!16000
176087
(0.75 secs, 642,165,016 bytes)

Антон

On Thu, Oct 10, 2019 at 1:52 PM Sergei M. Abramov abram_AT_botik.ru <
refal@botik.ru> wrote:

> > Мог ошибиться в границах, но смысл должен быть понятен: каждая
> > следующая граница (lim) получается из предыдущей примерно
> > возведением в квадрат, так ведь?
>
> Так, так, именно!
>
> > Поколение - то, что выходит после очередной фильтрации как кусок
> > списка простых и служит сомножителями в следующем pp.  Меня
> > интересует сколько (какую долю) отсеивает каждая следующая
> > фильтрация.  Понятно, что сначала это 1/2, далее 1-(2/3)*(4/5),
> > далее 1-(6/7)*(10/11)*...*(22/23) и т.д.  Произведение - это доля
> > остающихся.  Просто хотел бы увидеть десятичные числа, как они
> > меняются, наверно, убывают к 0.  Насколько быстро?
>
> Алик, я думаю, мое второе письмо (про bs) отвечает на этот вопрос.
> Хотя и не напрямую.
>
> >> Интересно, для Хаскела это тоже актуально, или там умножение pp*n
> ленивое?
> >
> >  Ленивое, несомненно.
>
> > Неужели?
>
> Ну да.  Пока кому-то не понадобиться это значение...
>
> > Значит на Хаскеле будет автоматически накапливаться в форме терма
> 7*11*13*..., да?
>
> Да.  А понадобиться это только, когда действительно случится
> прохождение через lim и начнется первый gcd с накопленным, но не
> посчитанным, произведением.  Вот там случится need результата.
>
> Все жадные операции перечислены в Haskell-е.  И я не вижу среди них
> умножения, тем более Integer.
>
> Всего доброго,
>
> Сергей Абрамов
>
>

Re: Список всех простых чисел

[Sergei M. Abramov]

> Мог ошибиться в границах, но смысл должен быть понятен: каждая
> следующая граница (lim) получается из предыдущей примерно
> возведением в квадрат, так ведь?

Так, так, именно!

> Поколение - то, что выходит после очередной фильтрации как кусок
> списка простых и служит сомножителями в следующем pp.  Меня
> интересует сколько (какую долю) отсеивает каждая следующая
> фильтрация.  Понятно, что сначала это 1/2, далее 1-(2/3)*(4/5),
> далее 1-(6/7)*(10/11)*...*(22/23) и т.д.  Произведение - это доля
> остающихся.  Просто хотел бы увидеть десятичные числа, как они
> меняются, наверно, убывают к 0.  Насколько быстро?

Алик, я думаю, мое второе письмо (про bs) отвечает на этот вопрос.
Хотя и не напрямую.

>> Интересно, для Хаскела это тоже актуально, или там умножение pp*n ленивое?
>  
>  Ленивое, несомненно.

> Неужели?

Ну да.  Пока кому-то не понадобиться это значение...

> Значит на Хаскеле будет автоматически накапливаться в форме терма 
> 7*11*13*..., да?

Да.  А понадобиться это только, когда действительно случится
прохождение через lim и начнется первый gcd с накопленным, но не
посчитанным, произведением.  Вот там случится need результата.

Все жадные операции перечислены в Haskell-е.  И я не вижу среди них
умножения, тем более Integer.

Всего доброго,

Сергей Абрамов

Re: Список всех простых чисел

[Arkady Klimov arkady . klimov_AT_gmail . com]

чт, 10 окт. 2019 г. в 18:32, Sergei M. Abramov abram_AT_botik.ru <
refal@botik.ru>:

> День добрый, всем!
>
> > То есть для последнего "поколения" эта работа лишняя.
>
> А что такое последнее поколение, если я строю список всех простых
> чисел?
>
> Да, всех, но не сразу ведь.
Сначала весь начальный список [2..] фильтруется на   ((== 1).(gcd pp)) для
pp = 2,
потом полученный список снова фильтруется для pp=3*5,
потом для 7*11*13*17*19*23,
потом для 29*31*...
Мог ошибиться в границах, но смысл должен быть понятен: каждая следующая
граница (lim) получается из предыдущей примерно возведением в квадрат,
так ведь?
Поколение - то, что выходит после очередной фильтрации как кусок списка
простых и служит сомножителями в следующем pp.
Меня интересует сколько (какую долю) отсеивает каждая следующая фильтрация.
Понятно, что сначала это 1/2, далее 1-(2/3)*(4/5), далее
1-(6/7)*(10/11)*...*(22/23) и т.д. Произведение - это доля остающихся.
Просто хотел бы увидеть десятичные числа, как они меняются, наверно,
убывают к 0. Насколько быстро?
А самому кодить недосуг :)

> > Наверно, лучше накапливать в виде списка, а при подаче на Filter
> > для следующего шага все перемножить.
>
> > Интересно, для Хаскела это тоже актуально, или там умножение pp*n
> ленивое?
>
> Ленивое, несомненно.
>
Неужели? Значит на Хаскеле будет автоматически накапливаться в форме терма
7*11*13*..., да?

>
> > А еще хотелось бы посмотреть статистику по числу (доле) отсеянных
> > (прошедших) через каждый из фильтров (это мне уже как математику
> > интересно).
>
> Алик, а мне интересно от математиков (может знакомые есть?) узнать,
> есть ли хоть какие-то работы (а если есть, то какие) в этой отрасли,
> где используется (\ k -> (gcd pp k) == 1)?
>
> Не знаю, и спросить кого тоже пока не приходит в голову.


> > Нужно спасать рассылку от офтопа.
>
> Спамера забанить!
>
> > Понятно, что вызов  > e.ns> на Рефале будет выполняться медленно.
>
> Вот, хочу сказать, что меня удивило:
>
>  ns' = filter (\ k ->  (gcd pp k) == 1) ns
>
> работает немного медленнее, чем
>
>  ns' = filter ((== 1).(gcd pp)) ns
>
> Вот нифига ж себе?
>
Да, забавно.

>
> Всего доброго,
>
> Сергей Абрамов
>
>
-- 
___
*С уважением, *
*Аркадий Климов,*
*с.н.с. ИППМ РАН,*
*+7(499)135-32-95*
*+7(916)072-81-48*

Re: Список всех простых чисел

[Sergei M. Abramov]

>> То есть для последнего "поколения" эта работа лишняя.

Про поколения.

Замечу, что если взять несколько первых простых чисел p1...pK, их
произведение pp и профильтровать (\ k -> (gcd pp k) == 1) весь
бесконечный список [(pk+1)..  ], то результат фильтрации имеет очень
простое, эффективное, конструктивное представление на Haskell-е.  Q
это в своих этюдах использовал.

То есть, несколько поколеий можно прогнать в статике, на этапе
написания программы.

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

Например, если взять 2 и 3.  pp=6, то результатом фильтрации будет
знакомое с детсва 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]

очень простая.

Всего доброго,

Сергей Абрамов

Re: Список всех простых чисел

[Andrei Klimov andrei_AT_klimov . net]

On Thu, Oct 10, 2019 at 6:32 PM Sergei M. Abramov abram_AT_botik.ru <
refal@botik.ru> wrote:

>
> Вот, хочу сказать, что меня удивило:
>
>  ns' = filter (\ k ->  (gcd pp k) == 1) ns
>
> работает немного медленнее, чем
>
>  ns' = filter ((== 1).(gcd pp)) ns
>
> Вот нифига ж себе?
>

В голову приходит такая гипотеза: ко второму варианту успешно применяются
какая-нибудь дефорестация, а на первом она ломается, так как второй
удовлетворяет некоторому шаблону (в аргументе fold'а композиция закарренных
функций), а первый воспринимается как общий вид и оставляется как есть.

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

Кстати, меня все-таки удивляет, что Simon P-J со
своим аспирантом Maximilian Bolingbroke не смогли подобрать вариант
суперкомпиляции, который работал бы в GHC за разумное время, как хотел
Simon (он об этом говорил). Подозрение такое: если остановится на
достаточно простых свистке и обобщении, то не набирался диссер
Bolingbroke'у. А когда увлеклись "научностью" для диссера (он весьма
солидный) , то не
получилось устойчиво (по времени) работающего свистка.

А может они потом сделали в GHS более богатую версию дефростации, которая
впитала какие-то идеи суперкомпиляции. Не знаю. Не натыкался на
соответствующие работы.

Всего наилучшего,
Андрей

Re: Список всех простых чисел

[Sergei M. Abramov]

День добрый, всем!

> То есть для последнего "поколения" эта работа лишняя.

А что такое последнее поколение, если я строю список всех простых
чисел?

> Наверно, лучше накапливать в виде списка, а при подаче на Filter
> для следующего шага все перемножить.

> Интересно, для Хаскела это тоже актуально, или там умножение pp*n ленивое?

Ленивое, несомненно.

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

Алик, а мне интересно от математиков (может знакомые есть?) узнать,
есть ли хоть какие-то работы (а если есть, то какие) в этой отрасли,
где используется (\ k -> (gcd pp k) == 1)?

> Нужно спасать рассылку от офтопа.

Спамера забанить!

> Понятно, что вызов  e.ns> на Рефале будет выполняться медленно.

Вот, хочу сказать, что меня удивило:

 ns' = filter (\ k ->  (gcd pp k) == 1) ns

работает немного медленнее, чем

 ns' = filter ((== 1).(gcd pp)) ns

Вот нифига ж себе?

Всего доброго,

Сергей Абрамов

Re: Список всех простых чисел

[Arkady Klimov arkady . klimov_AT_gmail . com]

Александр, а можно попросить вас пропустить вариант для 176100 вместо 4
- тогда результаты (времена) будут сопоставимы.
Я подумал еще об одной оптимизации. Сейчас произведение простых формируется
сразу, а ведь оно может и не понадобиться.
То есть для последнего "поколения" эта работа лишняя.
Наверно, лучше накапливать в виде списка, а при подаче на Filter для
следующего шага все перемножить.
Интересно, для Хаскела это тоже актуально, или там умножение pp*n ленивое?
А еще хотелось бы посмотреть статистику по числу (доле) отсеянных
(прошедших) через каждый из фильтров (это мне уже как математику интересно).
И соответствующие значения lim.
Я так понимаю, что количество поколений здесь 4-5, да?
Аркадий

ср, 9 окт. 2019 г. в 13:13, Александр Коновалов a.v.konovalov87_AT_mail.ru <
refal@botik.ru>:

> Добрый день всем!
>
> Нужно спасать рассылку от офтопа.
>
> Переписал исходник Сергея Михайловича на Рефал-5. Функции primesE,
> primesAD и primesA переписаны максимально близко к тексту, ради этого
> пришлось ввести функции filter, o (композиция функций), bind-right (чтобы
> связать Mod со вторым аргументом), eq и neq.
>
> Рефал не поддерживает ленивые вычисления и бесконечные списки, поэтому все
> функции на входе принимают последовательность чисел от 2 до 40 000.
>
> Понятно, что вызов 
> на Рефале будет выполняться медленно. Поэтому я также написал
> оптимизированные версии 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 
> Sent: Tuesday, October 8, 2019 10:29 PM
> To: refal@botik.ru
> Subject: Список всех простых чисел
>
>
>
> День добрый, всем!
>
>
>
> Простите, что не про Рефал.  Просто не могу вспомнить то, что (как мне
> помнится) случилось в обстановке рефал-компании на заре моей
> профессиональной деятельности (практика в НИЦЭВТе или работа в нем).
>
>
>
> По каким-то причинам я задумался тогда над построением списка всех простых
> чисел.  И в некий момент сам придумал алгоритм, которрый опирался вот на
> какое определение последовательности простых чисел:
>
>
>
> 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..].  И