Вероятность появления взаимно простых чисел

[wrest]

B3LYP в сообщении #1726193 писал(а):

мне хочется понять что-то в математике, чтобы не было “за деревьями не видно леса”.

B3LYP в сообщении #1726174 писал(а):

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

Тогда почитайте например про функцию Эйлера и её свойства.
https://ru.wikipedia.org/wiki/%D0%A4%D1 ... 1%80%D0%B0

[wrest]

B3LYP в сообщении #1726193 писал(а):

если никому тут не интересно обсуждать общие приёмы алгоритмической оптимизации, жаль...

Неинтересно тут другое: играть с вами в угадайку что вы там написали и как ускорили, не видя вашего кода

B3LYP в сообщении #1726193 писал(а):

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

Вот код который считает очень быстро, код очень простой, одна строчка на pari/gp:
pt(n)=my(s=sum(i=1,n,eulerphi(i)));return((2*s-1)/(n^2))
Но вот разобраться в нём для вас было бы неплохо, с математической стороны вопроса.
Это код считает искомую вероятность по формуле

$$P(n)=\dfrac{1}{n^2} \left( 2 \sum \limits _{k=1}^n \varphi(k) -1 \right)$$

где $\varphi(k)$ функция Эйлера
Почему это работает?

[B3LYP]

Мне хочется продолжить свои изыски, чтобы в итоге “отвелосипедить” функцию Эйлера; я статью в Вики про неё ещё не читал, и пока даже не просмотрел внимательно ответы в теме.
Переписал алгоритм и думал что он должен ускориться, но нет – чуть даже замедлился. Попробовал немного провести “среднеуровневую оптимизацию” (поварьировал размеры динамических массивов и пр.) – не помогло. Тогда я упростил код, убрав дополнительные массивы – стало чуть быстрее и заодно уменьшилась потребность в памяти. Итого новый алгоритм всего на 20% быстрее считает график выше, и ясно что без функции Эйлера тут никак.
В новом алгоритме изначально быстрее стал заполнятся массив со списком простых чисел для каждого перебираемого числа, но это ничего не изменило, поскольку лимитирующей стадией является перебор N*N пар чисел. Для каждой итерации этого перебора, разные алгоритмы ещё запускали один цикл или два вложенных цикла (по спискам простых делителей), но особой разницы не было. Если мы считаем P(100000), то итераций будет 10 миллиардов и на моём компьютере всё это считается примерно 5 минут.

wrest в сообщении #1726201 писал(а):

Вот код который считает очень быстро, код очень простой, одна строчка на pari/gp:

Как-то уж очень просто. Как я понимаю, вы с вашими высокоуровневыми языками программирования можете найти готовую библиотеку для чего угодно, а я привык изобретать велосипеды со своим Delphi. И вдруг я когда-то довелосипедируюсь до того, что придумаю что-то новое, чего в ваших библиотеках нет?

[Dmitriy40]

B3LYP в сообщении #1726231 писал(а):

лимитирующей стадией является перебор N*N пар чисел.

Так вот этого вот перебора и не нужно если знаем результат для N-1, достаточно перебора только N чисел (или одного вычисления функции Эйлера, что ещё быстрее).
Если не знаем, то плюс перебор всех чисел до N.
Но никаких внутренних циклов (ни одного ни двух) не надо.

[wrest]

B3LYP в сообщении #1726231 писал(а):

И вдруг я когда-то довелосипедируюсь до того, что придумаю что-то новое, чего в ваших библиотеках нет?

Непременно!
Но функцию Эйлера можно иногда и в уме вычислить, не то что в билиотеке или на дельфи
Например $\varphi(1000)=1000*(1-1/2)(1-1/5)=1000\cdot 1/2 \cdot 4/5 =400$ - чисел меньших 1000 взаимно просты с ним

-- добавлено через 9 минут --

B3LYP в сообщении #1726231 писал(а):

но это ничего не изменило, поскольку лимитирующей стадией является перебор N*N пар чисел.

Ну так если пара скажем (1000,1001) взаимно проста, то и пара (1001,1000) -- тоже. Зачем же перебирать обе?

[wrest]

B3LYP в сообщении #1726231 писал(а):

Если мы считаем P(100000), то итераций будет 10 миллиардов и на моём компьютере всё это считается примерно 5 минут.

Насчёт Delphi не знаю, а на Си без всяких библиотечных функций, то есть с честным решетом, P(10^5) на планшете считается за 0,03с. Потому что сложность должна быть не O(n^2), а окололинейная, т.е. в примерно n раз (в 10^5 в этом случае) быстрее.

[B3LYP]

Честно говоря, я пока даже стараюсь по возможности не читать ответы мне в теме (потом перечитаю конечно), потому что хочется самому отвелосипедить функцию Эйлера или что-нибудь ещё. Итак, как я бы ускорил мой алгоритм. Предположим, надо найти P(10). Соответственно старым алгоритмом это 100 итераций. Но мы можем перебирать все простые числа, и то что им кратно. Вначале взяли число 2, и перебрали 4, 6, 8, 10. Поскольку все эти числа кратны два, мы знаем что они заведомо не являются взаимно-простыми. Эти 5 чисел имеют 1+2+3+4=10=4*5/2 комбинаций, а если менять местами первое и второе число в комбинациях, то комбинаций 20. Т.е. мы знаем что из 100 всех комбинаций уже 20% это не взаимно простые числа. Далее надо сделать такой же перебор для 3, 6, 9, но тут уже надо обращать внимание на такой момент - число 6 уже входило в комбинации для чисел, кратных двум, важно не учитывать это дважды. Возможно, если всё это правильно учесть и сделать чёткий алгоритм, то можно будет очень даже всё ускорить. Я в верном направлении двигаюсь? Это подступ к функции Эйлера?

[Dmitriy40]

Да, если перебирать не все простые, а только делители N, то это фактически и получится вычисление через факторизацию N.

А решетом можно за практически линейное время вычислить сразу все P(1...N).

[wrest]

B3LYP в сообщении #1726517 писал(а):

но тут уже надо обращать внимание на такой момент - число 6 уже входило в комбинации для чисел, кратных двум, важно не учитывать это дважды. Возможно, если всё это правильно учесть и сделать чёткий алгоритм, то можно будет очень даже всё ускорить. Я в верном направлении двигаюсь? Это подступ к функции Эйлера?

Тут скорее подступ к принципу включений-исключений, ну и функции Эйлера заодно
https://ru.wikipedia.org/wiki/%D0%A4%D0 ... 0%B8%D0%B9
Пока стоит на запасном пути ещё функция Мёбиуса, ждёт когда освоите функцию Эйлера

[B3LYP]

Извиняюсь что очень медленно продвигаюсь с этой темой, типа совсем тугодум, но упомяну что я до этой темы не слышал ничего про функцию Эйлера.
Итак, что мне пришло в голову. Предположим, мы ищем P(12). Перебираем все простые числа, начиная с двух, и для каждого вторым циклом перебираем все делимые на них числа. Т.е. сначала перебираем 2, 4, 6, 8,10, 12 – всего 6 чисел. Начиная с четвёрки, помечаем эти числа как не-простые, чтобы потом повторно не перебирать. И ещё создаём массив с булевыми значениями, для которых проставляем true – это маркеры что мы уже посчитали суммы для этих чисел. Итого для первого простого числа – двойки – имеем 6 чисел, это 6*6=36 комбинаций, добавляем это в сумму. Потом перебираем все делимые числа для следующего простого числа 3: 3, 6, 9, 12, т.е. 4 числа. И тут надо иметь в виду, что два из этих чисел – 6 и 12 – уже помечены на предыдущей итерации как обработанные. Итого к сумме надо добавить 4*4-2*2=12, общая сумма будет 48. Далее берём простое число 5, т.е. 5 и 10 – это два числа, но одно из них (10) уже помечено раньше, итого к сумме прибавляем 2*2-1*1=3, сумма будет 51. Наконец, для 7 к сумме можно добавить только 1, и для 11 то же самое, итого сумма будет 53. Всё сошлось.
Алгоритм работает до N=14, а начиная с 15 делает ошибку. Раньше я думал что ошибка начнётся на N=30 - это первое число, являющееся произведением сразу трёх простых чисел (2*3*5). Чуть позже разберусь в чём дело. В целом мой вопрос – я пока двигаюсь в верном направлении?

[wrest]

B3LYP в сообщении #1726670 писал(а):

итого сумма будет 53. Всё сошлось.

Что с чем сошлось?

[B3LYP]

wrest в сообщении #1726691 писал(а):

Что с чем сошлось?

Количество комбинаций чисел, имеющих общий делитель, для N=12. В данном случае $P(12)=(12 \cdot 12-53)/(12 \cdot 12)=91/144=0.63194$.
Ошибка на 15 у меня возникает потому, что когда подсчитываются числа делимые на 5, программа вспоминает, что числа 10 и 15 были ранее уже помечены; но не учитывается, что они относились к разным группам - 10 к группе для двойки, а 15 к группе для тройки. Соответственно, по описанной выше формуле, программа добавляет к сумме $3 \cdot 3-2 \cdot 2$, а надо $3 \cdot 3-1 \cdot 1-1 \cdot 1$ (хотя и 10 и 15 были раньше помечены, к сумме комбинаций чисел, имеющих общий делитель, не была добавлена пара 10-15 и 15-10). Буду думать как это исправить, но что-то уж очень тяжко скрипеть мозгами(

[wrest]

B3LYP
Я бы предлжил идти не с начала, а с конца. Для 12 как-то так:

=== Подсчёт не взаимно простых пар для n = 12 ===  
Шаг 1: Инициализация pairs[d] для d <= n/2
pairs[2] = (12\2)^2 = 6^2 = 36
pairs[3] = (12\3)^2 = 4^2 = 16
pairs[4] = (12\4)^2 = 3^2 = 9
pairs[5] = (12\5)^2 = 2^2 = 4
pairs[6] = (12\6)^2 = 2^2 = 4
Шаг 2: Для всех d > n/2 pairs[d] = 1 (т.к. n\d = 1)
pairs[7] = 1
pairs[8] = 1
pairs[9] = 1
pairs[10] = 1
pairs[11] = 1
pairs[12] = 1

Шаг 3: Вычитаем кратные (проход с pairs[n\2] до pairs[2])
Внутренний цикл:  Для каждого делителя d 
  вычитаем из pairs[d] все бОльшие кратные 
  то есть все от pairs[2*d] до pairs[n] с шагом d

Обработка d = 6, pairs[6] = 4
  Вычитаем pairs[12] = 1
  После вычитания: pairs[6] = 3  
  
Обработка d = 5, pairs[5] = 4
  Вычитаем pairs[10] = 1
  После вычитания: pairs[5] = 3 
  
Обработка d = 4, pairs[4] = 9
 Вычитаем pairs[8] = 1
 Вычитаем pairs[12] = 1
 После вычитания: pairs[4] = 7

Обработка d = 3, pairs[3] = 16
 Вычитаем pairs[6] = 3
 Вычитаем pairs[9] = 1
 Вычитаем pairs[12] = 1
 После вычитания: pairs[3] = 11
  
Обработка d = 2, pairs[2] = 36
 Вычитаем pairs[4] = 7
 Вычитаем pairs[6] = 3
 Вычитаем pairs[8] = 1
 Вычитаем pairs[10] = 1
 Вычитаем pairs[12] = 1
 После вычитания: pairs[2] = 23 
 
Шаг 4: Суммируем все pairs[d] для d >= 2  
 result += pairs[2] = 23 
 result += pairs[3] = 11 
 result += pairs[4] = 7  
 result += pairs[5] = 3
 result += pairs[6] = 3  
 result += pairs[7] = 1
 result += pairs[8] = 1  
 result += pairs[9] = 1
 result += pairs[10] = 1 
 result += pairs[11] = 1
 result += pairs[12] = 1 
=== ИТОГО: 53 ===

То есть: в этом случае вы ходите всё время по одному массиву размером N (в примере выше это массив pairs , допонительных структур данных нет. Никакие билиотечные функции не используются, только целочисленное деление (в ексте выше обозрачено обратным слэшом \, которое должно быть и в Дельфи.

Алгоритмическая сложность по времени получается $O(n \log n)$, по памяти $O(n)$ - не лучшие показатели, ну так вы же хотите сами оптимизировать вот вам простор.