18.05.2019 15:26
wataru 43 6800 Источник
В комментариях к одному из прошлых постов о решете Эратосфена был упомянут этот короткий алгоритм из Википедии:

Алгоритм 1:

1: для i := 2, 3, 4, ..., до n: 
2:  если lp[i] = 0:
3:       lp[i] := i
4:       pr[] += {i}
5:   для p из pr пока p ? lp[i] и p*i ? n:
6:       lp[p*i] := p
Результат:
 lp - минимальный простой делитель для кажого числа до n
 pr - список всех простых до n.

Алгоритм простой, но не всем он показался очевидным. Главная же проблема в том, что на Википедии нет доказательства, а ссылка на первоисточник (pdf) содержит довольно сильно отличающийся от приведенного выше алгоритм.

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

Примечание о сложности алгоритма
Хоть этот алгоритм и асимптотически быстрее стандартного решета Эратосфена за O(n log log n), ему требуется гораздо больше памяти. Поэтому для по-настоящему больших n, где бы этот алгоритм засиял во всей красе, он не применим. Однако, он очень полезен даже для небольших n, если вам надо не только найти все простые числа, но и быстро раскладывать на множители все числа до n.

Определения


$\mathcal{P}$ — множество простых чисел.
$lp(x)$ — минимальный простой делитель числа: $lp(x) = min \{p \in \mathcal{P}\ \vert\ p\vert x \}$
$rd(x)$ — остальные множители: $rd(x) = x / lp(x)$

Обратите внимание, все определения выше даны для $x \ge 2$.

Некоторые очевидные свойства введенных выше функций, которые будут использоваться дальше:
$lp(a\times b) = min(lp(a), lp(b))$
$rd(x) < x$
$p \in \mathcal{P} => lp(p) = p$

Доказательство


Лемма 1: $lp(x) \le lp(rd(x)), \forall x \notin \mathcal{P}, x>1$
Доказательство: Т.к. любой делитель $rd(x)$ также является делителем $x$, а $lp(x)$ не превосходит любой делитель $x$, то $lp(x)$ не превосходит и любой делитель $rd(x)$, включая наименьший из них. Единственная проблема в этом утверждении, если $rd(x)$ не имеет простых делителей, но это возможно только в случае $x \in \mathcal{P}$, который исключен в условии леммы.
$\square$

Пусть $E = \{(lp(x), rd(x)) \vert \forall x \notin \mathcal{P}, 2\le x \le n\}$

Т.к. $lp(x)\times rd(x) = x$ (по определению $rd()$), если бы нам было дано множество $E$, то мы смогли бы вычислить $lp()$ для всех составных чисел до n. Это делает, например, следующий алгоритм:

Алгоритм 2:

1: Для всех (l,r) из E:
2:    lp[l*r] := l;

Заметим, что $\vert E \vert \le n$, поэтому алгоритм 2 выше работает за линейную сложность.

Далее я докажу, что алгоритм 1 из Википедии на самом деле просто перебирает все элементы этого множества, ведь его можно параметризовать и по-другому.

Пусть $E' = \{(p, r) \vert p \in \mathcal{P}, 2\le r < n, p \le lp(r), p \times r \le n\}$

Лемма 2: $E = E'$
Доказательство:

Пусть $(a,b)\in E$ => $\exists x \notin \mathcal{P},\ 2\le x\le n \mid a=lp(x), b=rd(x)$.

По определению $lp(), rd()$: $a \in \mathcal{P}$, $a\times b = x$. По лемме 1, $a \le lp(b)$.
т.к. $rd(x) < x$, то $b < n$.
Поскольку $x \notin \mathcal{P}$, $b \ge 2$.
Так же, по определению $E$, $x \le n$, следовательно, $a \times b \le n$.
Все 4 условия $E'$ выполнены, значит, $(a,b) \in E'$ => $E \subset E'$.

Пусть $(a,b)\in E'$. Пусть $x = a \times b$.
По определению $E'$, $a \in \mathcal{P}$. Значит, $a$ — простой делитеть $x$.
$lp(x) = lp(a\times b) = min(lp(a), lp(b)) = min(a, lp(b)).$
Т.к $a \le lp(b)$, то $lp(x) = a.$ Значит, $b = rd(x).$
По определению, $x = a \times b \le n.$ Также, очевидно, $x = a \times b \ge 2.$
Все условия для $E$ выполнены => $E' \subset E.$

Следовательно, $E = E'.$
$\square$

Теперь осталось перебрать правильные $r$ и $p$ из определения множества $E'$ и применить алгоритм 2. Именно это и делает Алгоритм 1 (только вместо r используется переменная i). Но есть проблема! Что бы перебрать все элементы множества $E'$ для вычисления $lp,$ нам надо знать $lp,$ ведь в определении есть условие $p \le lp(r)$.

К счастью, эта проблема обходится правильным порядком перебора. Надо перебирать $r$ во внешнем цикле, а $p$ — во внутреннем. Тогда $lp(r)$ уже будет вычислено. Этот факт доказывает следующая теорема.

Теорема 1:
Алгоритм 1 поддерживает следующий инвариант: После выполнения итерации внешнего цикла при i=k, все простые числа до k включительно будут выделены в список pr. Также будет подсчитано $lp()$ для всех $x \notin \mathcal{P} \mid x\le n,\ rd(x) \le k$ в массиве lp.

Доказательство:
Докажем по индукции. Для k=2 инвариант проверяется вручную. Единственное простое число 2 будет добавлено в список pr, потому что массив lp[] изначально заполнен нулями. Также, единственное составное число у которого $rd(x) \le 2$ — это 4. Несложно убедиться, что внутренний цикл выполнится ровно один раз (при n>3) и правильно выполнит lp[4] := 2.

Теперь допустим, что инвариант выполняется для итерации i=k-1. Докажем, что он будет выполнятся и для итерации i=k.

Для этого достаточно проверить что число i, если оно простое, будет добавлено в список pr и что $lp()$ будет подсчитано для всех составных чисел $x \le n,$ т.ч. $rd(x) = i.$ Именно эти числа из инварианта для k не покрыты уже инвариантом для k-1.

Если i простое, то lp[i] будет равно 0, ведь единственная операция записи в массив lp, которая теоретически могла бы подсчитать lp[i] (в строке 6), всегда пишет в составные индексы, ведь p*i (для i > 1) — всегда составное. Поэтому число i будет добавлено в список простых. Также, в строке 3 будет подсчитано lp[i].

Если же i составное, то на начало итерации lp[i] уже будет подсчитано по инварианту для i=k-1, ведь $rd(i) < i$ или $rd(i) \le i-1,$ следовательно i попадает под условие инварианта в предыдущей итерации. Поэтому составное число i не будет добавлено в список простых чисел.

Далее, имея корректное lp[i] и все простые числа до i включительно цикл в строках 5-6 переберет все элементы $(p, i) \in E'$, у которых вторая часть равна i:

  • $p \in P,$ потому что оно из списка pr
  • $p \le lp(i),$ по условию останова цикла
  • $p\times i \le n,$ по условию останова цикла
  • $i < n,$ следует из $p\times i \le n,\ p > 1$

Все нужные простые числа в списке pr есть, т.к. нужны только числа до $lp(i) \le i$. Поэтому будут подсчитаны значения lp[] для всех составных чисел $x$, у которых $rd(x)=i$. Это ровно все те числа, которых не хватало при переходе от инварианта для k-1 к инварианту для k.

Следовательно, инавриант выполняется для любых i = 2..n.

$\square$

По инварианту из Теоремы 1 для i=n получается, что все простые числа до n и все lp[] будут подсчитаны алгоритмом 1.

Более того, поскольку в строках 5-6 перебираются различные элементы множества $E$, то суммарно внутренний цикл выполнит не более $\vert E \vert < n$ операций. Операция проверки в цикле выполнится ровно $\vert E \vert + n - 1 < 2n$ раз ($\vert E \vert$ раз вернет true и n-1 раз, для каждого i, вернет false). Все оставшиеся операции вложены в один цикл по i от 2 до n.
Отсюда следует, что сложность алгоритма 1 — O(n).

Комментарии (43)


  1. pohjalainen
    18.05.2019 18:44
    #20170408
    +1

    Поясните последний абзац, почему из двух вложенных циклов не более чем по n получается O(n), а не O(n^2)?


    1. wataru Автор
      18.05.2019 18:58
      #20170468
      +1

      Вложенный цикл суммарно выполняет O(n) операций. Это как, допустим 


      step = n;
for (i = 0; i < n; ++i) {
    step /= 2;
    for (j = 0; j < step; ++j) {
         do_somthing();
    }
}

      Вложенный цикл выполнит n/2 + n/4 + n/8 +… операций. Cуммарно — меньше n.


      Еще одна аналогия — обход в ширину графа c n ребрами (если ребра заданы списком смежности). Вроже бы опять 2 вложенных цикла, но внутренняя итерация суммарно выполнится на более n раз.


      1. pohjalainen
        18.05.2019 20:33
        #20170698
        +1

        n/2 + n/4 + n/8 +… = n(1-(1/2)^([log_2 n]))(1-1/2)
        Асимптотически получаем все-таки O(n)
        Внутренний цикл — O(n), внешний цикл — n…
        Интересно, что в англоязычной вики этого алгоритма нет.


        1. wataru Автор
          18.05.2019 20:40
          #20170722
          +2

          Внутренний цикл — O(n), внешний цикл — n…

          Ну да, большинство итераций внешнего цикла не сделают почти ничего. Цикл for посмотрит на условие и не выполнит ни одну итерацию. Также и алгоритм в этой статье.


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


          Я подозреваю, что этот алгоритм — плод народного творчества. Поэтому нигде не удалось найти объяснение алгоритма с доказательством.


  1. Rsa97
    18.05.2019 18:53
    #20170452
    +1

    На мой взгляд, здесь модифицированное решето Эратосфена, и сложность будет примерно та же. Внешний цикл перебирает все числа, внутренний — только простые, меньшие найденного.
    Для n верхняя граница количества простых чисел в интервале [1,n] в первом приближении оценивается как ln(n). То есть, общая оценка алгоритма — O(n*log(n)).


    1. wataru Автор
      18.05.2019 19:02
      #20170480
      +1

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


      Внутренний цикл перебирает только простые меньше lp(i), что, например, для всех четных чисел — всего 2. Т.е. перебираются сильно меньше чисел.


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


      Для n верхняя граница количества простых чисел в интервале [1,n] в первом приближении оценивается как ln(n). То есть, общая оценка алгоритма — O(n*log(n)).

      Кстати, верхняя оценка количества простых — n/log n, а не log n, как вы написали.


      1. Rsa97
        18.05.2019 21:08
        #20170796
        +1

        Действительно, я невнимательно посмотрел. Каждое непростое число вычёркивается только один раз, значит сложность O(n).
        А КМП и обход дерева в ширину реализуются без вложенных циклов, так что тут ваши примеры несколько неудачны.


        1. wataru Автор
          18.05.2019 22:26
          #20170964

          А КМП и обход дерева в ширину реализуются без вложенных циклов, так что тут ваши примеры несколько неудачны.

          Ну, КМП, допустим, можно реализовать через конечный автомат (кстати, моя вторая статья как раз про это в том числе). Но это немного другой алгоритм, да и я сомневаюсь, что вы это имели ввиду. Как обход в ширину делается без вложенных циклов — я вообще не в курсе.


          Очень интересно. Могли бы вы привести наброски кода? И для кмп или хотя бы только для обхода в ширину.


          1. Rsa97
            19.05.2019 01:04
            #20171348

            С КМП реализуется легко и без конечного автомата:

            Заголовок спойлера
            var needle = "abracadabra";
var needle_pref = [0];
var i = 0;
var j = 1;
while (i < needle.length) {
  if (needle[i] == needle[j]) {
    needle_pref[i++] = ++j;
  } else if (j == 0) {
    needle_pref[i++] = 0
  } else {
    j = needle_pref[j-1];
  }
}

var haystack = "abyrvalgabramabracadabracadabrass";
i = 0;
j = 0;
while (i < haystack.length) {
  if (haystack[i] == needle[j]) {
    if (j == needle.length - 1) {
      console.log(i-j, haystack.substr(i-j, j+1));
      i++;
      j = needle_pref[j];
    } else {
      i++;
      j++;
    }
  } else if (j == 0) {
    i++;
  } else {
    j = needle_pref[j];
  }
}


            1. wataru Автор
              19.05.2019 09:42
              #20171784

              Действительно, можно было и так реализовать КМП, спасибо!


      1. printercu
        19.05.2019 08:58
        #20171732
        +1

        Ну покажите же, пожалуйста, цифры :)


        Искать неточность в доказательстве сложнее, чем замерить производительность. У меня получилась нелинейная зависимость времени от n для вашего алгоритма:


        Calculating -------------------------------------
                1000      2.755k (± 3.9%) i/s -     13.860k in   5.038808s
               10000    273.275  (± 2.9%) i/s -      1.378k in   5.046871s
              100000     26.691  (± 3.7%) i/s -    134.000  in   5.027327s
             1000000      2.590  (± 0.0%) i/s -     13.000  in   5.021461s
             2000000      1.309  (± 0.0%) i/s -      7.000  in   5.349981s

Comparison:
                1000:     2755.2 i/s
               10000:      273.3 i/s - 10.08x  slower
              100000:       26.7 i/s - 103.23x  slower
             1000000:        2.6 i/s - 1063.97x  slower
             2000000:        1.3 i/s - 2105.49x  slower

        код
        # ruby
def fn(n)
  pr = []
  lp = Array.new(n)
  2.upto(n) do |i|
    unless lp[i]
      lp[i] = i
      pr.push i
    end
    pr.each do |p|
      break unless p <= lp[i] && p * i <= n
      lp[p * i] = p
    end
  end
  pr
end

puts fn(100).join(' ')

require 'benchmark/ips'
puts Benchmark.ips { |x|
  [1_000, 10_000, 100_000, 1_000_000, 2_000_000].each do |n|
    x.report(n.to_s) { fn(n) }
  end
  x.compare!
}


        1. adukhounik
          19.05.2019 09:42
          #20171788

          стороной проходил, но вроде бы линейная по вашему отчету)


        1. wataru Автор
          19.05.2019 11:00
          #20171970

          Eсли точно подсчитать, то сложность будет ?((c1 — c2/logn) * n).


          Это следует из того, что количество составных чисел ~ n — n/log n, а алгоритм выполняет действия для каждого из них.


          Обратите внимание, эта сложность остается O(n), потому что с ростом n слагаемое c1*n доминирует над c2*n/logn.


          Ваши же цифры вполне подходят под эту теорию — даже n log log n растет сильно быстрее, чем время работы вашего решения.


          Плюс к этому, на больших массивах оно перестает влезать в кэш и работает медленнее.
          Мерить точное время работы при скачках по 4-64Мб массиву — дело неблагодарное.


          Давайте лучше подсчитаем количество операций:


          Код
          long long Erato(unsigned int n) {
    int cnt1 = 0, cnt2 = 0, cnt3 = 0;
    vector<int> pr;
    pr.reserve(n-1);
    vector<int> lp(n+1);
    unsigned int i, j;
    long long sum = 0;
    for (i = 2; i <= n; ++i) {
      if (lp[i] == 0) {
        lp[i] = i;
        sum += i;
        pr.push_back(i);
        ++cnt1;
      }
      int maxp = n / i;
      for (j = 0; ++cnt3 && j < pr.size() && pr[j] <= lp[i] && pr[j] <= maxp; ++j) {
        ++cnt2;
        lp[i*pr[j]] = pr[j];
      }
    }
    cout << n << ' ' << cnt1 << ' ' << cnt2 << ' ' << cnt3 << '\n';
    return sum;
}


  1. Rsa97
    18.05.2019 21:07
    #20170790
    +1

    del


  1. mikhaelkh
    18.05.2019 23:17
    #20171080

    А в какой вычислительной модели оценка O(N)?

    Тут есть ?(N) умножений двух чисел длины log(N), то есть ?(N * log(N)) битовых операций.


    1. wataru Автор
      18.05.2019 23:56
      #20171174

      Модель такая же, как и для обычного решета — умножения происходят за константу.


      Так-то, если битовую сложность считать, то за n надо брать не границу до которой ищутся простые числа, а ее длину. В такой модели сложность вообще будет O(2^n*n*log n). Еще и памяти n*2^n надо.


      1. mikhaelkh
        19.05.2019 01:03
        #20171346

        за n надо брать не границу до которой ищутся простые числа, а ее длину

        А почему не длину этой длины, степенные башенки должны быть выше!
        Если серьёзно, то длина это как раз граница до которой ищутся простые числа — n. Или КМП у вас работает за 2^n, а числа перемножаются за n*2^n?


        1. wataru Автор
          19.05.2019 09:36
          #20171776

          А почему не длину этой длины, степенные башенки должны быть выше!

          Потому что берется длина входных данных в битах, если вы битовую сложность хотите.
          В кмп как раз длина 8*n. И при перемножении чисел в k знаков сложность считается от k, а не значения числа.


      1. mikhaelkh
        19.05.2019 22:36
        #20173716

        У вас n — это граница до которой ищутся простые числа, зачем называть этой буквой ещё и её длину?


        1. wataru Автор
          19.05.2019 22:55
          #20173764

          Что бы использовать привычную нотацию вида O(f(n)). Так-то можно обозначить такой буквой, какой мы захотим. Но тут вы правы, я зря Вас запутал этим.


          В любом случае, если вы хотите считать битовые операции, то и стандартное решето будет не O(n*log log n), как везде написано, а O(n*log n*log log n), ведь там числа до n складываются.


          1. mikhaelkh
            21.05.2019 08:14
            #20178524

            Есть решето Аткина с временем работы ?( n / log(log(n)) ) — статья
            Там нет умножений чисел размера O( log(n) ), только сложения-вычитания и прочие линейные операции, т.е. битовая сложность ?( n * log(n) / log(log(n)) )
            При сложности умножения ?(n * log(n)) у вас битовая сложность ?(n * log n * log log n), такая же как и у обычного решета Эратосфена.
            Т.е. «сэкономиленные» log(log(n)) вы на самом деле спрятали в сложность умножения.


            1. wataru Автор
              21.05.2019 10:43
              #20179120

              такая же как и у обычного решета Эратосфена.

              Да, в этой модели вычислений — вы правы.


              1. mikhaelkh
                21.05.2019 11:40
                #20179486

                От log(log(n)) легко избавиться, если заметить, что мы вычисляем p*i, для фиксированного i, а p бежит по простым. Пользуясь тем, что (pr[k]-pr[k-1]) = O(pr[k] ^ (2/3)), можно предпосчитать q*i для нужных q, а потом умножение свести к сложению предыдущего с предпосчитанным.


                1. wataru Автор
                  21.05.2019 12:45
                  #20179816

                  Была такая мысль, но я так и не придумал, как это сделать быстро.


                  Ведь для разных i использутеся разный префикс pr. Так что, мы можем много-много итераций не трогать конец списка. А потом должны как-то быстро подсчитать (pr[k]-pr[k-1])*i.


                  Или я вас не правльно понял?


                  1. mikhaelkh
                    21.05.2019 16:37
                    #20181018
                    +1

                    Я же зафиксировал i, значит для разных i считается независимо.
                    Так как (pr[k]-pr[k-1]) = O(pr[k] ^ (2/3)) асимптотически мало, мы можем предпосчитать все нужные произведения (pr[k]-pr[k-1])*i и это будет асимптотически меньше чем последующие сложения pr[k]*i = pr[k-1]*i + (pr[k]-pr[k-1])*i.


                    1. wataru Автор
                      21.05.2019 16:46
                      #20181052

                      Так как (pr[k]-pr[k-1]) = O(pr[k] ^ (2/3)) асимптотически мало,

                      Это все-равно 2/3 log n бит. Вы таким образом ускорите умножения, но не ассимптотически.


                      1. mikhaelkh
                        21.05.2019 16:51
                        #20181076
                        +1

                        Ну будет O(k^(2/3)*log(k)*log(log(k))), это асимптотически меньше чем O(k*log(k))


                      1. mikhaelkh
                        21.05.2019 17:35
                        #20181248
                        +1

                        Точнее O(k^(2/3)*log(k)*log(n)) у предпосчёта против O(k*log(n))
                        k — количество простых в цикле для данного i.
                        Предподсчёт тоже будем не в лоб умножать, а складывать.


                        1. wataru Автор
                          21.05.2019 18:10
                          #20181476

                          Почему вы количество простых-то берете в степени? все теже итерации сделать предется, все столько-же i*p значений подсчитать. Да, теперь вместо умножения двух чисел в log n бит вы умножаете число в 2/3 log n бит (pr[k]-pr[k-1]) на число в log n бит (i). Все те же n умножений во всем алгоритме, но числа на 33% короче. На асимтотику не влияет.


                          1. mikhaelkh
                            21.05.2019 19:00
                            #20181714
                            +1

                            Потому что смысл предподсчёта — посчитать i*j только для j до k^(2/3)*log(k), так как pr[k] = O(k*log(k)), а разности между соседними простыми — O(pr[k]^(2/3)). А потом основной цикл — предыдущее складываем с предподсчитанным.


                            1. wataru Автор
                              21.05.2019 19:16
                              #20181786

                              Я думал k у вас — сколько раз внутренний цикл для данного i выполнится. Туда по определению записано ограничение i*p <= n. А теперь мне кажется, что k у вас количество простых в списке pr для данного i.


                              Но это все не важно. Суть в том, что то же самое количество значений i*p вам все-равно придется вычислить, ведь это все вычеркнутые по одному разу числа и их надо обязательно вычеркнуть.


                              Теперь вместо каждого умножения ip вы предлагаете делать prev + i(pr[j]-pr[j-1]), так? Количество умножений не изменилось, потому что мы считаем все те же значения. Сами умножения работают на на ~15% быстрее, потому что один из множителей на треть короче, но на асимптотику это не влияет.


                              1. mikhaelkh
                                21.05.2019 19:33
                                #20181842
                                +1

                                Все i(pr[j]-pr[j-1]) мы предпосдчитали заранее в отдельный массив и заново умножать не надо, только складывать. Всего O(k) сложений, но так как (pr[j]-pr[j-1]) = O( k^(2/3)*log(k) ), то умножений в предподсчёте гораздо меньше.


                                1. wataru Автор
                                  21.05.2019 19:41
                                  #20181866

                                  Все i(pr[j]-pr[j-1]) мы предпосдчитали заранее

                                  Вот этих комбинаций, реально используемых O(n) и будет. Ровно столько, сколько умножений делает наивный алгоритм.


                                  Я не понимаю, как из размера разности между простыми числами следует, что понадобится меньше умножений. Ведь каждое используемое значение, которое мы предподсчитываем, делается умножением. А нужно нам в каждом внутреннем цикле все те же k, суммарно O(n).


                                  Если вас не затруднит, хотя бы псевдокод вашего решения приведите, пожалуйста.


      1. mikhaelkh
        21.05.2019 20:21
        #20182062
        +1

        for (i=2; i<=n; ++i) {
	if (lp[i]==0)
		pr[++np] = lp[i] = i;
	m = min(lp[i], n/i), d=0;
	for (k=0; k<np && pr[k+1]<=m; ++k)
		d = max(d, pr[k+1]-pr[k]);
	for (j=1; j<=d; ++j)
		a[j] = a[j-1]+i;
	pi=0;
	for (j=1; j<=k; ++j)
		pi += a[pr[j]-pr[j-1]], lp[pi] = pr[j];
}

        От n/i можно избавиться, но тогда код станет сложнее.


        1. wataru Автор
          21.05.2019 21:19
          #20182280

          Шикарно вообще! Ни единого умножения. Вы — гений. Я вообще вас не правильно сначала понял.


          Единственное, надо подумать, а не станет ли d сильно больше k? И не сломает ли это асимптотику? Ведь если для каких-то определенных значений k Ваш алгоритм может сделать сильно больше операций и если это k довольно часто встречается, то все может поломаться.


          Кажется, что не должно, ведь разница между простыми до n в среднем log n. Но это в среднем. И судя по Вики эта разница всегда меньше количества простых. Надо только убедится, что так и будет для очень больших простых.


          P.S. жаль, не могу вам плюс в карму поставить!


          1. mikhaelkh
            21.05.2019 21:50
            #20182366

            > а не станет ли d сильно больше k?

            d = O(k ^ (2/3)), и это не самая лучшая известная оценка.

            Есть Гипотеза Крамера, утверждающая, что d = O(log(m)^2) = O(log(k)^2)

            Так что всё, кроме вычисления n/i, выполняется в сумме за O(n*log(n)) битовых операций. Ну а как справиться с n/i подумайте сами.


  1. Shador
    19.05.2019 09:22
    #20171760

    rd(x) < x ведь всегда, а не только если x ? P. Если x — простое, то rd(x) = 1 < x.


    1. wataru Автор
      19.05.2019 09:38
      #20171782

      Исправил, спасибо.


  1. NotThatEasy
    19.05.2019 12:31
    #20172214

    Линейная скорость, это, конечно, замечательно, да вот только, вроде бы, в проджект Эйлере, например, в первой же десятке заданий есть решето и там, насколько я помню, линейная скорость не то чтоб слишком годится, я не готов минуту ждать расчёта против варианта за 100-300 ms.


    1. wataru Автор
      19.05.2019 12:42
      #20172246

      В примечании сразу после хабраката я указал, что на практике этот алгоритм работает плохо из-за огромных потребностей в памяти. Этот алгоритм используется прежде всего, когда надо иметь возможность раскладывать все числа до n на простые множители.


      Но все же, если вы про 10-ое задание на Эйлере, то там надо подсчитать простые всего до 2 миллионов. Тут что линейное решето требует 0,02с, что стандартное 0.015с на моем компьютере. И вся разница объясняется необходимостью работать с 8mb памяти вместо 250kb.


      1. NotThatEasy
        19.05.2019 13:11
        #20172320

        Мог тупануть и сравнить с временем исполнения poorly-written хаскеля и с++, ручками оптимизированную, в таком случае прошу прощения.
        Всё же 20мс на линейном решете мне кажутся не то чтоб самым реальным результатом, наверняка улучшали его, как минимум, если все числа > 2 крутить в цикле, то будут отжираться огромные куски времени проца на молочение впустую составных чисел.


        1. wataru Автор
          19.05.2019 14:08
          #20172470

          Ну помилуйте, цикл до 2 миллионов срабатывает мгновенно. На олимпиадах всегда был такой прием — просто считаем операции, вроде присвоения, if-ов и т.д. В секунду таких 400-600 миллионов точно влезет. Тут же выполнится, дай бог 6-8 миллионов. Как раз 1/50-ая


          Код и счетчики операций есть в этом комментарии. Никаких оптимизаций почти нет.


          Потом, по-моему, вы что-то путаете. То что O(n log log n) работает за 15мс у вас вопросов не вызывает, а в то что O(n) работает за 20мс вы поверить не можете.


          1. NotThatEasy
            21.05.2019 07:01
            #20178432

            С готовностью признаю, что никогда не оценивал сложность алгоритмов, помимо пары пар в универе.
            Ежели Ваша оценка сложности вычислена верно, то подобные значения вполне правдоподобны, я же просто приводил цифры, что меня удивили; в любом случае, благодарю за занятный диалог.