Решение упражнения 1.28 из SICP
Это упражнение немного заковыристо сформулировано, но сложным не является.
Вот код модифицированной процедуры expmod, которая сигнализирует о нахождении нетривиального квадратного корня из 1, возвращая 0 (как и рекомендовано в условии).
(define (square x) (* x x))
(define (apply-trivial-check k m r) (if (and (not (= k 1)) (not (= k (- m 1))) (= r 1)) 0 r))
(define (remainder-or-trivial k m) (apply-trivial-check k m (remainder (square k) m)))
(define (modified-expmod base exp m) (cond ((= exp 0) 1) ((even? exp) (remainder-or-trivial (modified-expmod base (/ exp 2) m) m)) (else (remainder (* base (modified-expmod base (- exp 1) m)) m))))
Теперь запишем всю процедуру, реализующую тест Миллера-Рабина. Она прогоняет для заданного числа n модифицированный тест Ферма на всех числах a<n и подсчитывает количество успехов. Если их больше половины, то число в соответствии с тестом Миллера-Рабина идентифицируется как простое. В противном случае оно составное.
(define (miller-rabin-test n) (define (miller-rabin-iteration a t n) (define (try-it a) (= (modified-expmod a (- n 1) n) 1)) (if (= a n) (> t (/ n 2)) (if (try-it a) (miller-rabin-iteration (+ a 1) (+ t 1) n) (miller-rabin-iteration (+ a 1) t n)))) (miller-rabin-iteration 1 0 n))
Проверим написанную процедуру на нескольких простых и составных числах:
> (miller-rabin-test 5) #t > (miller-rabin-test 15) #f > (miller-rabin-test 97) #t > (miller-rabin-test 121) #f > (miller-rabin-test 1003) #f > (miller-rabin-test 1009) #t > (miller-rabin-test 100003) #t > (miller-rabin-test 100005) #f
Теперь проверим первые 6 чисел Кармайкла:
> (miller-rabin-test 561) #f > (miller-rabin-test 1105) #f > (miller-rabin-test 1729) #f > (miller-rabin-test 2465) #f > (miller-rabin-test 2821) #f > (miller-rabin-test 6601) #f
Как видим, определение простоты проходит корректно во всех случаях.
Comments
Comment from thror
Date: February 14, 2008, 11:45 am
И все же, несмотря на корректность работы процедуры на всех входных данных…
Общая схема работы теста Миллера-Рабина для числа n такова:
[1] Выбирается s случайных значений величины a из {1,…,n-1}
[2] Если одно из них свидетельствует о том, что n—составное, то n признается составным. Такой вывод всегда верный.
[3] Если в ходе s попыток не было таких свидетельств, то предполагается, что нет причин считать число n составным, так что оно считается простым.
Справедливы следующие теоремы (Кормен, Лейзерсон, Ривест, Штайн. Алгоритмы—построение и анализ, 2-е издание. Теоремы 31.38, 31.39)
Если n—нечетное составное число, то количество свидетельств того, что n—составное, не меньше (n-1)/2.
При любом нечетном n > 2 и положительном целом s вероятность того, что тест Миллера-Рабина выдаст неправильный результат не превышает 1/2**s.
Comment from thror
Date: February 14, 2008, 11:48 am
В дополнение два варианта, реализующих этот тест.
;; SICP
(define (square x)
(* x x))
(define (sqrmod x m)
(let ((y (remainder (square x) m)))
(if (and (= y 1) (not (= x 1)) (not (= x (- m 1))))
0
y)))
(define (expmod base exp m)
(cond ((= exp 0) 1)
((even? exp) (sqrmod (expmod base (/ exp 2) m) m))
(else (remainder (* base (expmod base (- exp 1) m)) m))))
(define (witness? a n)
(not (= (expmod a (- n 1) n) 1)))
(define (miller-rabin-test n)
(not (witness? (+ 1 (random (- n 1))) n)))
(define (fast-prime? n times)
(cond ((= times 0) true)
((miller-rabin-test n) (fast-prime? n (- times 1)))
(else false)))
Comment from thror
Date: February 14, 2008, 11:49 am
;; CLRS
(define (square x)
(* x x))
(define (sqrmod x m)
(remainder (square x) m))
(define (ntz n)
(cond ((or (= n 0) (odd? n)) 0)
(else (+ 1 (ntz (/ n 2))))))
(define (ntz n)
(define (iter m r)
(cond ((or (= m 0) (odd? m)) r)
(else (iter (/ m 2) (+ r 1)))))
(iter n 0))
(define (rtz n)
(cond ((or (= n 0) (odd? n)) n)
(else (rtz (/ n 2)))))
(define (expmod base exp m)
(cond ((= exp 0) 1)
((even? exp) (sqrmod (expmod base (/ exp 2) m) m))
(else (remainder (* base (expmod base (- exp 1) m)) m))))
(define (witness? a n)
(define (iter x t)
(let ((y (remainder (square x) n)))
(cond ((= t 0) (not (= x 1)))
((and (= y 1) (not (= x 1)) (not (= x (- n 1)))) true)
(else (iter y (- t 1))))))
(iter (expmod a (rtz (- n 1)) n) (ntz (- n 1))))
(define (miller-rabin-test n)
(not (witness? (+ 1 (random (- n 1))) n)))
(define (fast-prime? n times)
(cond ((= times 0) true)
((miller-rabin-test n) (fast-prime? n (- times 1)))
(else false)))
Comment from thror
Date: February 14, 2008, 11:54 am
Один из них, как ясно из комментариев, “навеян” SICP, другой—CLRS =)
Оба приведены к общей схеме, и соответствующие функции названы одинаково, для большей ясности.
–thror
Comment from Sergey Khenkin
Date: February 21, 2008, 11:52 pm
Спасибо за комментарии. К сожалению, долго не возвращался к этому упражнению. Вернувшись, перечитал соответствующий раздел SICP, а также главу из замечательной книги Кормена, Лейзерсона и Ривеста “Алгоритмы: построение и анализ”.
В итоге хочу сказать, что наши подходы во многом сходны, хотя есть и одно очень существенное различие. Мне больше нравится разбивка на процедуры, сделанная не мной, а throrом (в первом варианте по мотивам SICP), поэтому я буду объяснять в ее терминах.
Действительно, в книге Кормена и др. тест Миллера-Рабина дается как вероятностный тест. Запущенный s раз на случайных значениях a<n он ошибается с вероятностью не более 1/2s в силу того, что если n не простое, то модифицированная процедура expmod, описанная выше, обнаружит нетривиальный квадратный корень из 1 по модулю n по крайней мере для половины a<n. Именно этот подход положен в основу проверки на простоту, предложенной throrом. И, как мне сейчас кажется, именно так на практике и используется тест Миллера-Рабина.
Но есть и другой ход мыслей. В условии упражнения сказано, что в нем идет речь о тесте “который невозможно обмануть”. Из этого можно заключить, что речь идет не о вероятностном тесте, а о тесте, который достоверно дает ответ на вопрос, является ли число простым. Мы знаем, что что если n не простое, то модифицированная процедура expmod, описанная выше, обнаружит нетривиальный квадратный корень из 1 по модулю n по крайней мере для половины a В терминах процедур, введенных thror-ом, этот подход можно записать так: Я полагаю, что этот вариант решения также имеет право на жизнь и является, по сути, правильным.
(define (fast-prime? n)
(define (count-witnesses from)
(if (= from n)
0
(if (witness? a n)
(+ 1 (count-witnesses (+ from 1)))
(count-witnesses (+ from 1)))))
(< (count-witnesses 1) (/ (- n 1) 2)))
Хотя с практической точки зрения, скорее всего, он менее полезен и применим.
Comment from Sergey Khenkin
Date: February 22, 2008, 9:03 am
В сущности, можно и не считать количество свидетелей того, что n является составным числом. Ведь достаточно найти хотя бы одно свидетельство для того, чтобы опровергнуть гипотезу о простоте числа. Таким образом, процедуру можно упростить и записать так:
(define (fast-prime? n)
(define (any-witness? from)
(cond ((= from n) false)
((witness? from n) true)
(else (any-witness? (+ from 1)))))
(not (any-witness? 1)))
Однако это все же далеко от практического применения.
Ведь определение простоты числа этим методом требует порядка Θ(n log n) шагов, что для больших чисел n неприемлемо. В то же время вероятностный тест, предложенный thror-ом, для s попыток требует только Θ(s log n) шагов и дает корректный результат с вероятностью не меньше, чем 1-2-s.
Comment from gorilych
Date: September 27, 2008, 7:45 am
кхм. Я нигде не вижу сравнения результата modified-expmod с нулём. Таким образом, вся ваша модификация коту под хвост - решение неправильное!
Общий алгоритм решения таков:
1) Берём n/2+1 чисел от 1 до n.
2) Последовательно вычисляем для каждого числа a из этого списка (modified-expmod a (- n 1) n)=res
Если res=0 - это значит, что мы нашли “нетривиальный квадратный корень из 1 по модулю n” - n - составное (однозначно)
Если res!=1 - это значит, что необходимое условие простоты n не удовлетворено, то есть опять n - составное (однозначно)
Если res=1 - продолжаем тестировать оставшиеся числа.
3) Если ни одно из чисел не выявило “нетривиальный квадратный корень из 1 по модулю n” - n - простое.
Comment from gorilych
Date: September 27, 2008, 8:16 am
а, ну да, оно и не должно сравниваться с нулём. однако у вас в miller-rabin-test подсчитывается, сколько раз результат modified-expmod был не равен единице. Достаточно одного раза - это уже гарантирует, что число n - составное.
Comment from трор
Date: November 2, 2008, 4:05 pm
и все же, уважаемый горилич, уверяю вас, пока что природа происходящего для вас не ясна. Самой простой рекомендацией было бы-читать цлрс до просветления, или самому придумать алгоритм логарифмической сложности для подобного анализа. Однако это показалось бы слишком грубым. Поэтому я попробую еще раз растолковать.
Итак, задача в том, чтоб потратив к*лог(н) действий выяснить-составное число, или простое?
Для этого мы выбираем к чисел случайно, и каждого используем как свидетеля составной природы числа. Если хоть один свидетель найдем, то число составное.
Витнесс, или на латинском-свидетель, в моих процедурах возвращает правду, когда число не простое. В этом вы легко убедитесь и сами. Для составного числа модифицированный экспмод вернет либо ноль, либо какое нибудь не равное единице число. Но не единицу.
А далее, следуем простой схеме-если к свидетелей не утверждают, что число составное, то считаем его простым. Как в жизни.
Откуда же взялось это загадочное н/2, которое способно даже грамотных (судя по 1.45) людей ввести в заблуждение? Зачем нужна эта громадная цифра? А для чисел в четыре тясячи бит длиной она действительно потрясает воображение.
Ответ кроется в простой идее, которую не трудно доказать-для нечетного составного числа н не менее (н-1)/2 свидетелей того, что оно составное.
Это обеспечивает экспоненциальное убывание вероятности ошибки при росте к.
Для сергея-на практике именно так он и реализован. Вы можете, к примеру, обратиться к старым реализациям библиотек шифрования из 4.4бсд-не помню сейчас точно, или План9-а там очень хорошая идея как модифицируют алгоритм на практике.
Вообщем-на самом деле то не слишком важно, правильное решение или нет-важно то, на какие идеи оно нас наталкивает, и о чем заставляет задуматься.
Comment from трор
Date: November 4, 2008, 10:05 am
вы можете взглянуть на код полной криптосистемы РСА (включающей свозможность генерировать длинный ключ случайно-но прежде чем утверждать что ключ получается не слишком случайный прочтите пожалуйста кнута 3.2.2 и упражнение 18 к ней, а также 3.3.2-критерии серий, и не говорите, что я не в курсе=), позволяет шифровать и расшифровывать сообщения, создавать цифровую подпись, создавать подписанные сообщения, и проверять подписи, и предоставлющую функции взлома приватного ключа по публичному, основанные как на классическом методе разложения, так и на эвристическом ро-методе полларда. На моей странице http://www.uic.nnov.ru/~kuaa29
–трор
Write a comment