lion137
2018-11-02 20:41

Pakiety w Pythonie

Dawno nie było nic o Pythonie, to jest. Znowu dość praktycznie: pakiety. Jak tworzymy, jak są izolowane, jak wygląda w nich praca.
Jeśli chcemy napisać coś większego w Pythonie, dobrze jest zacząć od stworzenia środowiska wirtualnego:

➜  test git:(master) ✗ python3.7 -m virtualenv --no-site-packages venv
Using base prefix '/usr/local'
New python executable in /home/lion/code/workspace/python/tests/test/venv/bin/python3.7
Not overwriting existing python script /home/lion/Code/workspace/python/tests/test/venv/bin/python (you must use /home/lion/Code/workspace/python/tests/test/venv/bin/python3.7)
Installing setuptools, pip, wheel...
done.

Aktywujemy je:

➜  test git:(master) ✗ source venv/bin/activate 

I jak teraz sprawdzimy, jakie mamy zainstalowane pakiety, lista jest pusta (pip, jako nakładka na setuptools służy do instalacji, listowania, itp...):

(venv) ➜  venv git:(master) ✗ pip freeze
(venv) ➜  venv git:(master) ✗ 

Lista jest pusta. Kilka uwag, Komenda --no-site-packages słuźy do kompletnego odizolowania środowiska (tak dla pewności, Python domyślnie może tego nie robić). W naszym mini unixowym drzewie katalogów mamy plik wykonawczy Pythona (i kilka innych), jest to dokładnie ta sama wersja, którą użyliśmy do stworzenia środowiska(i będzie do niego przypisana na zawsze). Jesteśmy odseparowani od reszty, cokolwiek nie robimy jesteśmy w naszym venv.
Instalujemy coś:

(venv) ➜  venv git:(master) ✗ pip install pytest
Collecting pytest
 Successfully installed atomicwrites-1.2.1 attrs-18.2.0 more-itertools-4.3.0 pluggy-0.8.0 py-1.7.0 pytest-3.9.3 six-1.11.0
(venv) ➜  venv git:(master) ✗ python
Python 3.7.0 (default, Sep 15 2018, 00:51:17) 
[GCC 5.4.0 20160609] on linux
Type "help", "copyright", "credits" or "license" for more information.
>>> import pytest
>>> pytest.__file__
'/home/lion/code/workspace/python/tests/test/venv/lib/python3.7/site-packages/pytest.py'
>>> 

Jak widać pytest jest w naszym lokalnym (nie globalnie) venv.
Teraz pip freeze:

(venv) ➜  venv git:(master) ✗ pip freeze
atomicwrites==1.2.1
attrs==18.2.0
more-itertools==4.3.0
pluggy==0.8.0
py==1.7.0
pytest==3.9.3
six==1.11.0
(venv) ➜  venv git:(master) ✗ 

Zależności, jeżeli przekierujemy to wyjście do pliku:

(venv) ➜  venv git:(master) ✗ pip freeze > requirements.txt && cat requirements.txt
atomicwrites==1.2.1
attrs==18.2.0
more-itertools==4.3.0
pluggy==0.8.0
py==1.7.0
pytest==3.9.3
six==1.11.0
(venv) ➜  venv git:(master) ✗ 

To mamy gotową listę zależności, teraz możemy ją spakować (plus cokolwiek co mamy w projekcie) i w nowym środowisku wirtualnym możemy my (lub ktokolwiek inny), zainstalować nasz projekt:

(venv) ➜  venv git:(master) ✗ deactivate 
➜  venv git:(master) ✗ python3.7 -m virtualenv --no-site-packages venv2     
Using base prefix '/usr/local'
New python executable in /home/lion/code/workspace/python/tests/test/venv/venv2/bin/python3.7
Also creating executable in /home/lion/code/workspace/python/tests/test/venv/venv2/bin/python
Installing setuptools, pip, wheel...
done.
➜  venv git:(master) ✗ ls
bin  include  lib  requirements.txt  venv2
➜  venv git:(master) ✗ source venv2/bin/activate
(venv2) ➜  venv git:(master) ✗ pip install -r requirements.txt
Successfully installed atomicwrites-1.2.1 attrs-18.2.0 more-itertools-4.3.0 pluggy-0.8.0 py-1.7.0 pytest-3.9.3 six-1.11.0
(venv2) ➜  venv git:(master) ✗ 

Ciekawostka: jeśli zedytujemy requirements.txt, zmienimy wersję jakiegoś pakietu i ponownie odpalimy pip install -r req..., to pip zainstaluje nową i usunie poprzednią, czyli nie mamy opcji zrobienia chaosu z różnymi wersjami pakietu w projekcie:)
To tyle, jak wszystko poszło dobrze:), to powinno śmigać:). Pochwalcie się czymś na githubie:)
#Python #python #theory

KageYoshi

W windowsie ręczna instalacja jest tragedią :/ setki dziwnych błędów, u każdego spowodowane czym innym. Dla własnego zdrowia psychicznego używam już tylko condy.

lion137
2018-09-27 11:29

Typowany Python?

Witam, po przerwie, jak w tytule, co tam z tzw. "type annotations" w Pythonie, otóz mają się nieźle; mypy coraz lepiej zintegrowane z trójką (najlepiej działa 3.5>=). Praca polega (przynajmnije na Linux) na zainstalowaniu mypy:

$ python -m pip install mypy

Co ono robi? Sprawdza typy zgodnie z naszymi annotacjami, które [te annotacje] nie mają wpływu na wykonanie programu, tzn.:

$ python my_file.py

będzie działać, chociaż:

$ mypy my_file.py

zwróci błędy. Czyli do istniejącego kodu, można je dodawać stopniowo i sprawdzać, testować, etc... Kilka przykładów, jak ktoś chce to znajdzie resztę:)

Deklaracja wygląda tak:

from typing import List
x: List[int] = []
x.append(1)

mypy zaprotestuje, jeśli dodamy do listy cokolwiek innego niż int.

Mamy generyki i type variables :

from typing import Generic, TypeVar

A = TypeVar('A')

class Tree(Generic[A]):
    def __init__(self, value: A, left: 'Tree[A]', right: 'Tree[A]') -> None:
        self.payload = value
        self.left_tree = left
        self.right_tree = right

Żeby interpreter się nie wysypał, parametry typów w __init__ umieszczamy w apostrofach. Funkcja zwracająca None, to odpowiednik void z Javy, czy C/C++.

Kolejna funkcja to mapa (w sensie jak w programowaniu funkcyjnym), zwracająca stream - w Pytonie są to generatory. Callable[[A], A] - przyjmuje argument typu A i zwracająca takiż element : (A) -> A. Drugi argument Iterable[A] jest supertypem dla listy (może to być, również Dict, Tuple, Range):

from typing import Callable, Iterable, Generator

def my_map_generator(f: Callable[[A], A], xs: Iterable[A]) -> Generator:
    for x in xs:
        yield f(x)

Jak akurat nie chcemy strumienia, to proszę bardzo:

def my_map(f: Callable[[A], A], xs: List[A]) -> List[A]:
    for ind, x in enumerate(xs):
        xs[ind] = f(x)
    return xs

Ta wersja bierze taką samą funkcję, jak poprzednia, i listę typu A, a zwraca listę typu A.
To tyle, opcja warta rozważenia w dużych projektach, przy skryptach, to typowy "overkill". Pozdrawiam.
#Python #python #theory #programming

Leroy

Pracowałem w korpo pythonowym projekcie na Python'ie 2 gdzie w sumie wdrazylem MyPy.
Samo narzedzie jest fajne, ale dalej, Python byl, jest i bedzie dynamicznie typowany.
Spielismy to z CI, napisalismy wlasne toole, ktore pozwalaly na ustawianie roznych poziomow skanowania i walidacji w roznych podprojektach i pakietach.
Ogolnie wg mnie MyPy dziala, ale to tylko proteza, ktora mimo ze pomaga, ma swoje wady i czesto natrafialismy na sciane gdzie musielismy pisac wlasne stub'y.
Najwiekszym plusem bylo to, ze PyCharm po prostu zaczal dzialac.
O wiele lepszym rozwiazaniem jest po prostu statycznie typowany jezyk.

lion137
2018-07-24 12:32

Python, Dwie Sztuczki i Algorytm

Dawno nic nie pisałem, tak to już jest, czasem człowiek jest busy, albo lazy:). Do rzeczy, dwa małe, praktyczne, triki (jakby ktoś nie znał). Pierwszy, otwieranie pliku przy pomocy iteratora(załóżmy, że mamy otwarty plik, na przykład z logami, jako f):

def good_lines(f):
    for line in f:
        line = line.strip("\n ")
        if line.startswith('#'):
            continue
        if not line:
            continue
        yield line

Utworzony iterator, omija linie komentarzy, puste linie, obcina spacje i znaki końca linii, teraz można już robić coś z interesującymi nas partiami pliku:

for line in good_lines(f):
    do_stuff(line)

Argumentem do funkcji good_lines(f) niekoniecznie musi być plik, może to być lista, zbiór, generalnie każda struktura wspierająca iterację.

Dalej,
W Pythonie nie ma łatwej możliwości przerwania podwójnej pętli, trzeba tworzyć dodatkowe zmienne kontrolne, warunki, co raczej zaciemnia kod, utrudnia testowanie, itd... Np.:

for y in range(height):
    for x in range(width):
        if sth_found(x, y):
            break # nie przerwie obu pętli

Ale gdy tylko się da możemy zrobić tak:

def 2d_range(width, height):
    for y in range(height):
        for x in range(width):
            yield x, y

for col, row in 2d_range(width, height):
    if sth_found(col, row):
        break # Cool, works!:)

Na koniec Fibonacci, jak można znaleźć w necie jest kilka sposobów liczenia tej funkcji, od najgorszego:

  • naiwny rekurencyjny, O(n) i pamięć również O(n), więc po chwili rozsadza nam stos;
  • dynamic programming, dalej O(n), ale w stałej pamięci, więc można jako tako używać;
  • matrix exponentiation, lepiej O(logn), stała pamięć, bardzo szybki;
  • double fibonacci identities, ulepszenie poprzedniego, złożność ta sama, ale mniej operacji stałych.

Wszystko to, dobrze opisane, np., tutaj. Właśnie ten przedostatni algorytm poniżej, w zasadzie jest to zapisanie identyczności z tytułu, trzeba tylko pamiętać o exponentiation by squaring, które można też, jak widać, stosować dla macierzy. Może się okazać, że ten tip jest nawet praktyczny (interview), przed kliknięciem zachęcam do samodzielnego zakodowania, link tutaj:
https://nbviewer.jupyter.org/[...]trix_Exponetiation_Fibo.ipynb

#Python #python #algorithms #theory

lion137

Pewnie, że się da, ale czasem mogą być problemy z wysłaniem i mutowaniem zmiennych (nie można sobie wyslać adresu jak w C++ :)), więcej klepania, a tak jest prościej i przejrzyściej.

Afish

No ale jakie wysyłanie zmiennych? Masz domknięcie, nie musisz nic przekazywać przez parametry. Pokaż może przykład, gdzie to nie zadziała.

lion137
2018-06-06 02:27

Różniczkowanie Symboliczne

Ha, ha, dawno już nie było teorii, wszyscy stęsknieni:P; to do roboty: Symbolic differentiation!
Temat jest bliski moim zainteresowaniom, tzn. AI; więc, studiując słynne PAIG plus (obowiązkowo!) SICP natkąłem się na jak w temacie.
Jako że ciężko mi było nadążyć z oryginalnym kodem w Lispie:

(defun deriv-poly (p x)
"Return the derivative, dp/dx, of the polynomial p."
" If P is a number or a polynomial with main-var > x,
" then p is free of x, and the derivative is zero;
" otherwise do real work.
" But first, make sure X is a simple variable,
" of the form #(X 0 1).
(assert (and (typep x 'polynomial) (= (degree x) 1)
(eql (coef x 0) 0) (eql (coef x 1) 1)))
(cond
« numberp p) 0)
«var> (main-var p) (main-var x)) 0)
(var= (main-var p) (main-var x))
d(a + bx + cx2 + dx3)/dx = b + 2cx + 3dx A 2
;; So. shift the sequence p over by 1. then
;; put x back in. and multiply by the exponents
(let «r (subseq pI)))
(setf (main-var r) (main-var x))
(loop for i from 1 to (degree r) do
(setf (coef r i) (poly*poly (+ i 1) (coef r i))))
(normalize-poly r)))
(t ;; Otherwise some coefficient may contain x. Ex:
d(z + 3x + 3zx2 + z2x3)/dz
;; = 1 + 0 + 3x2 + 2zx A 3
;; So copy p. and differentiate the coefficients.
(let «r (copy-poly p)))
(loop for i from 0 to (degree p) do
(setf (coef r i) (deriv-poly (coef r i) x)))
(normalize-poly r)))))

Plus pomoc z SICP (Scheme):

(define (deriv exp var)
(cond ((number? exp) 0)
((variable? exp) (if (same-variable? exp var) 1 0))
((sum? exp) (make-sum (deriv (addend exp) var)
(deriv (augend exp) var)))
((product? exp)
(make-sum
(make-product (multiplier exp)
(deriv (multiplicand exp) var))
(make-product (deriv (multiplier exp) var)
(multiplicand exp))))
(else
(error "unknown expression type: DERIV" exp))))

(Nawiasem pisząc, formatka 4programmers ssie (nie formatuje Lispu)), to wysmażyłem swój w Pythonie (kompilowalny pseudokod :-D). Idea jest prosta: Na najniższym poziomie - pochodna liczby jest zerem, pochodna zmiennej jest jedynką, jeśli to ta sama zmienna, zerem jeśli inna (jak stała). Problemem jest, jak to odwzorować w notacji infiksowej (jak Python i inne języki głównonurtowe). Okazuje się, że można wykorzystać algorytm parsowania formuły do (binarnego na razie!!) drzewa składniowego.
Idea: Pochodna sumy jest sumą pochodnych, a pochodna iloczynu, jest sumą iloczynu pierwszego członu razy pochodna drugiego plus pochodna pierwszego razy drugi; doskonale wpisuje się to w drzewo składniowe:

def evaluate(tree):
    opers = {'+': sym_add, '-': op.sub, '*': sym_mul, '/': op.truediv}

    leftT = tree.getLeftChild()
    rightT = tree.getRightChild()

    if leftT and rightT:
        fn = opers[tree.getRootVal()]
        return fn(evaluate(leftT), evaluate(rightT))
    else:
        return tree.getRootVal()

def evaluate_parse_tree(tree, var):
    opers = {'+': sym_add, '-': op.sub, '*': sym_mul, '/': op.truediv}

    leftT = tree.getLeftChild()
    rightT = tree.getRightChild()

    if leftT and rightT:
        fn = opers[tree.getRootVal()]
        # changes start:
        if fn is sym_mul:
            fn = sym_add
            return fn(sym_mul(evaluate(leftT), evaluate_parse_tree(rightT, var)), sym_mul(evaluate_parse_tree(leftT, var),
                    evaluate(rightT)))
        else:
            return fn(evaluate_parse_tree(leftT, var), evaluate_parse_tree(rightT, var))
    else:
        return derrive(tree.getRootVal(), var)

def derrive(exp, var):
    if is_number(exp): return 0
    elif is_variable(exp):
        return 1 if same_variable(exp, var) else 0

Dwie funkcje evaluate_parse_tree i evaluate działają zgodnie z powyższym schematem; chociaż widać tu dużo braków: Trzeba upraszczać końcowe wyrażenia (albo poprawić sym-add i sym_mul). Lecz pomysł działa. Kod (już nie GitHubie, bo go nie ma:-D), z funkcjami sym_add i tak dalej (w razie nieścisłości zapraszam do dyskusji) tutaj: https://bitbucket.org/lion137[...];fileviewer=file-view-default
Działa?

 form = "((1 + (2 * x)) * (1 + x))"
   p_tree = BinaryTree('')
   p_tree = build_parse_tree(form)
   print(evaluate_parse_tree(p_tree, "x"))# -> 1(1 + 2x) + (0 + 2 + 0x)*(1 + x)
   form = "(x + 1)"
   p_tree = BinaryTree('')
   p_tree = build_parse_tree(form)
   print(evaluate_parse_tree(p_tree, "x")) # -> 1

Jak widać tak (chociaż): Wyrażenie musi być poprawnie znawiasowane i forma wyniku przedstawia sporo do życzenia.
Ale droga jest właściwa! :-D)
#theory #Python

#theory #Python

lion137
2018-01-30 13:34

Faktoryzacja Do Liczb Pierwszych

Jak już były testy(hashtag Primes), to teraz faktoryzacja, artykuł na Wikipedii szczegółowo opisuje dostępne metody. Wiele z nich (tych do dużych liczb) jest skomplikowana i raczej dla nas (jeśli nie pracujemy w kryptografii czy coś w tym stylu) mało użyteczna.

Ale w matematyce rekreacyjnej, jak adventofcode, hackerrank, projecteuler czy inne, wystarczy, w miarę szybko faktoryzować do, powiedzmy 10^16. A to można łatwo zrobić mając obliczoną tablicę liczb pierwszych do pierwiastka z n.

def factorize_trial(n, primes):
    """Factorizing n"""
    if miller_rabin(n):
        return n
    result = []
    for p in primes:
        while n % p == 0:
            result.append(p)
            n //= p
        if miller_rabin(n):
            result.append(n)
            break
    return result

primes można wygenerować używając, na przykład tego: simple_sieve.
Użyłem też, trochę zoptymalizowanego miller_rabin:

def miller_rabin(n):
    d = n - 1
    s = 0
    while d % 2 == 0:
        d >>= 1
        s += 1
    for repeat in range(20):
        a = 0
        while a == 0:
            a = random.randrange(n)
        if not miller_rabin_pass(a, s, d, n):
            return False
    return True

def miller_rabin_pass(a, s, d, n):
    a_to_power = pow(a, d, n)
    if a_to_power == 1:
        return True
    for i in range(s - 1):
        if a_to_power == n - 1:
            return True
        a_to_power = (a_to_power * a_to_power) % n
    return a_to_power == n - 1

Algorytm, jak widać iteruje po primes, próbując dzielić i dołączać do tablicy wynikowej, na końcu pętli jest sprawdzenie pierwszości, żeby ewentualnie zakończyć, gdy po podzieleniu n stanie się pierwsze.
Złożoność będzie primes(sqrt(n)) - liczba liczb pierwszych do pierwiastka z n, czyli będzie działał najwolniej, gdy liczba będzie mieć dwa duże czynniki pierwsze, ale dla wiekszości jest w miarę szybki. Kod. #Primes #Theory #Python #Numerical

lion137
2018-01-29 09:20

Lucas - Lehmer Test

Żeby zakończyć przegląd efektywnych Prime Tests, część druga, część pierwsza dodam test Lucasa - Lehmera, sprawdza on tylko liczby pierwsze Mersenne'a - Zastosowania w kryptografii - jest szybszy niż Miller - Rabin. Artykuł na Wikipedii wyjaśnia dokładnie ideę i implementację. Dodam tylko kod w Pythonie:

def lucas_lehmer(p):
    s = 4
    m = exp(2, p) - 1
    for i in range(p - 2):
        s = ((s * s) - 2) % m
    return True if s == 0 else False

W sumie jest w tej samej klasie złożoności co Miller Rabin ale wykonuje mniej operacji. Program sprawdził 27 liczbę Mersenne'a w 79 sekund, gdzie Miller Rabin wyłączyłem po kilkunastu minutach.
Złożoność algorytmu O(p^2 log p log log p) - zakładając liniowy algorytm mnożenia (jak np. w GMP)

Czyli generalne program do sprawdzania dowolnych liczb: najlepszy wybór Miller - Rabin.
Do liczb Mersenne'a: Lucas - Lehmer.

Z probablistycznych jeszcze warto wspomnieć o teście Frobeniusa, jego jedna runda jest dłuższa, ale osiąga szybsze ograniczenie prawdopodobieństwa.
Są jeszcze algorytmy deterministyczne, jak na przykład AKS, ale jest niewydajny O((logn)^6).
Poza tym, AKS jest deterministyczny, ale w praktyce, jeśli Miller Rabin I AKS zwrócą różne wyniki na tej samej liczbie, to jest większa szansa, że to promieniowanie kosmiczne(!) przypadkiem zmmieniło bit i wynik w komputerze, na którym był zapuszczony AKS (gdyż liczył dłużej!). #Primes #Python #Theory #Numerical

Aryman1983

Takie wpisy to ja rozumiem, a nie płacze nad swoją niedolą :-) Dzięki!

lion137
2018-01-26 23:38

Miller Rabin Test

Obiecany Miler Rabin, część pierwsza. Test Fermata, był prawie dobry, tylko dla pewnej klasy liczb (niewielkiej, ale jednak) zawsze się mylił (były złożone, ale przechodziły). Okazuje się, że można to usunąć, (Miller Rabin).

Najpierw Twierdzenie Fermata trochę inaczej zapisane:
Jeśli n jest liczbą pierwszą, a a dowolną liczbą dodatnią mniejszą od n, to a podniesione do potęgi n - 1przystaje do 1 modulo n.

Czyli(opis algorytmu):
Bierzemy losową liczbę 0 < a < n i podnosimy ją do potęgi n - 1 modulo n;
Podnosząc do potęgi, w każdym momencie podnoszenia do kwadratu szukamy "nietrywialnego pierwiastka kwadratowego jedynki modulo n", czyli sprawdzamy czy kolejna liczba, nie równa jeden lub n - 1 i podniesiona do kwadratu nie jest równa jeden modulo n (jeśli taką znajdziemy to liczba jest złożona). Dla liczby nieparzystej i złożonej, dla co najmniej połowy liczb a < n -1, sprawdzanie tym sposobem odkryje nietrywialny pierwiastek. Dlatego procedura Miller Rabin nie może być oszukana i można w granicy zredukowac błąd do zera gdy liczba prób dąży do nieskończoności. Pseudokod:

test_nontrivial(a, n):
    # t and u, t >=1, u odd and 2^tu = n - 1
    x[0] = modular-exponenation(a, u, n);
    for i -> 1 to t:
        x[i] = x[i - 1]^2 mod n;
        if x[i] == 1 and x[i - 1] != 1 x[i - 1] != n - 1:
            return True (composite)
        if x[t] != 1:
            return True
    return False

Dzięki sztuczce z dobraniem t i u, na początku i ustaleniu x[0], potem, aby podnieść a do n -1 modulo n, wystarczy podnieść x do kwadratu t razy, co umożliwia sprawdzanie warunku. Po iteracji, gdy nie mamy jedynki if x[t] != 1 to zwracamy True i finalnie False (nie znaleźliśmy nietrywialnych pierwiastków).
Teraz wystarczy tylko wykonać odpowiednia ilość testów w pętli, aby uzyskać zadaną dokładność, pseudokod:

miller-rabin(n, s):
    for j <- 1 to s:
        a = random(1, n - 1);
        if nontrivial_root(a, n):
            return False (composite)
        else:
            return True

I kod w Pythonie (nie wiele się różniący:)):

# 2. Miler - Rabin

# procedure nontrivial_root, which looking for a non trivial square root of one mod n

def nontrivial_root(a, n):
    """checking Fermat Theorem, and non trivial square root"""

    # find t and u, such that u odd, t >= 1 and n - 1 = 2^tu:
    t, u = 0, n - 1
    while u % 2 == 0:
        t += 1
        u //= 2

    x0 = modularexponenation(a, u, n)
    for i in range(1, t + 1):
        x1 = x0 ** 2 % n
        if x1 == 1 and x0 != 1 and x0 != n - 1:
            return True
    if x1 != 1:
        return True
    return False

def miller_rabin(n, s):
    """Miler Rabin Test"""

    # Return True if n = 2
    if n == 2:
        return True
    # Return False if even
    if n % 2 == 0:
        return False

    for i in range(s):
        a = randint(1, n - 1)
        if nontrivial_root(a, n):
            return False
    return True

Dokładność, można udowodnić, że błąd (gdy miller_rabin zwróci True) wynosi najwyżej 1 / 2^s, czyli, zgodnie z tą odpowiedzią na stacoverflow: https://stackoverflow.com/a/6330138, s = 40 wystarcza w zupełności. I na koniec złożoność obliczeniowa: O(k * (logn)^3) - jak daleko odeszliśmy od wykładniczej (co prawda w bitach:)). Kod na githubie. Mam nadzieję, że wywód był jasny, to tyle, dziękuję za cierpliwość. #Theory #Primes #Python #python #Numerical

enedil

Tbh, pominąłeś wyjaśnienie najbardziej kluczowej kwestii, dlaczego właściwie można się ograniczyć do sprawdzania tych właśnie liczb, jako kandydatów na nietrywialne pierwiastki.

lion137

@enedil: Zobacz też w części pierwszej, nic nie pominałęm, tak jest definicja liczb Carmichael'a, i procedura nontrivial_root te liczby wykrywa, co powoduje, że miller-rabin jest ulepszonym testem Fermata - Przeprowadza go dla s liczb, przy okazji sprawdzając czy poza tym, że spełniają one warunki twierdzenia Fermata, nie są relatywnie pierwsze do n (czyli Carmichael numbers). Finalnie wracamy do warunków z pierwszego testu - mamy liczbę pierwszą z prawdobodobieństwem.

lion137
2018-01-26 01:49

Primality Tests, Test Fermata

Jak już powstał wpis o algorytmach numerycznych, to w temacie nie można nie wspomnieć o liczbach pierwszych, a zacząć wypada od testów pierwszości. Wypada również pominąć naiwne testy, oparte na szukaniu dzielników, bezużyteczne dla dużych liczb (w teorii liczb bilion to mała liczba) - mają złożoność wykładniczą w bitach.

Jako pierwszy, chyba najprostszy, opierający się na Małym Twierdzeniu Fermata:
Jeśli n jest liczbą pierwszą, a jest jakąkolwiek dodatnią liczbą mniejszą od n, to a podniesione do n - tej potęgi przystaje do a modulo n
Dwie liczby przystają do siebie modulo n, jeśli mają tą samą reszte przy dzieleniu przez n.

Idea, Jeśli n jest złożone, to większość liczb nie będzie spełniać relacji z twierdzenia, co prowadzi do algorytmu:
Mając dane n, weź liczbę 0 < a < n i oblicz a do potęgi n modulo n, jeśli wynik nie równa się a, to liczba jest (na 100%) złożona, jeśli się równa, to są spore szanse na to, że n jest pierwsza (skoro w przypadku złożonej większość liczb a nie spełnia zależności); teraz weź kolejne losowe 0 < a < n i jeśli spełnia, to szanse na pierwszość n rosną jeszcze bardziej, i każde kolejne sprawdzenie zwieksza naszą pewność. Taki algorytm, to Test Fermata.

Potrzebne do algorytmu potęgowanie modulo, znajduje się, na przykład tutaj. Kod (Python):

# Primality Tests
# Copyleft 2018 lion137
# 1. Fermat Test

from Euler import modularexponenation
import random

def fermat_test_single(n):
    """performs single Fermat test"""
    def test(a):
        return modularexponenation(a, n, n) == a
    return test(random.randint(1, n - 1))

Dość oczywista zgodność kodu ze specyfikacją. Następna funkcja wykonuje test zadaną ilośc razy, dając zaskakująco skuteczny test:

def fermat_extended(n, cnt):
    """performing fermat tests in a loop,
    for #cnt of numbers (0, n) """
    for i in range(cnt):
        if not fermat_test_single(n):
            return False
    return True

Można ją potestować, ale najpierw o ograniczeniach. Po pierwsze, jako algorytm zrandomizowany nie zawsze jest poprawny, tzn., istnieje bardzo duża szansa, że daje dobry wynik, ale tylko szansa, precyzyjniej, jeśli test nie przejdzie, to liczba, na pewno jest złożona, ale jeśli przejdzie, to jest pierwsza z prawdopodobieństwem. Ale i to nie jest prawdą, są liczby, które są złożone, ale przechodzą test (mają własność, że dla każdego a, a do n przystaje do a modulo n), są to Carmichael Numbers. Liczby te są niezwykle rzadkie - poniżej stu milionów jest ich tylko 255, ale jednak. Pora na testy, tutaj liczby:

    print(fermat_extended(2 ** 2203 - 1, 5)) # -> True (a ** b w Pythonie to a do potęgi b)
    print(fermat_extended(2 ** 4423 - 1, 5)) # -> True
    print(fermat_extended(137, 6))           # -> True
    print(fermat_extended(6601, 6))          # -> True, złożona (Carmichael)
    print(fermat_extended(2821, 6))          # -> True, złożona (Carmichael)
    print(fermat_extended(999999999999999999, 6)) # -> False (oczywiście złożona)

Jesli chcemy przetestować losową liczbę pierwszą, to ta procedura, z praktycznych względów będzie najlepsza, prawdopodobieństwo błędu na losowym m bitowym numerze dąży do zera, gdy m dąży do nieskończoności (szansa błędu na losowej 1024 bitowej liczbie jest jeden do 10^41 całkiem nieźle:)). Ale gdy liczby nie są losowe to już tak nie jest - na liczbach Carmichael'a, błędu wyeliminować się nie da. To zrobi dopiero (ale później:) ) algorytm, bazujący na teście Fermata, Miller - Rabin. Jeszcze tylko złożoność :O(k (logn)^2 log log n * log log log n),całkiem przyzwoita - dowód raczej do pominięcia:-D. To tyle.
Żródła: Wikipedia, CLRS. #Theory #Primes #python #Python #Numerical

lion137
2018-01-23 14:06

Parę Algorytmów Numerycznych

Przeglądając sobie stare notatki, znalazłem jak męczyłem jakieś algosy w C (pewnie z Knutha albo innej CLRS), a jak łatwo i miło się to pisze w Pythonie:)
Na rozgrzewkę wszytkie podzbiory zer i jedynek danej długości (czyli powerset: dla n = 2 -> 00, 01, 11, 10); idea jest prosta: Zaczynając od zera (liczby binarnej, jako string), dodajemy jeden, aż do zadanej długości, prefiksując, jak potrzeba zerami, kod w Pythonie, gładko to odwzorowujący, wydajnie, jako iterator(pozdrawljaju, pisać to w C:) - chociaż, jeśli interpreter Pythona jest napisany w C, to wszystko co w Pythonie, musi się dać napisać i w C):

def bin_pow_set(n):
    """returns power set of the all bit sequences of
    length n as an iterator"""
    def append_zeroes(x, n):
        """helper function to append zeroes"""
        x = bin(x)
        x = x[2:]
        return '0' * (n - len(x)) + x
    i = 0
    k = 0
    yield '0' * n
    while k < 2 ** n - 1:
        k += 1
        i = append_zeroes(k, n)
        yield i

Zdefiniowana funkcja pomocnicza prefuksuje zerami, a w pętli inkrementujemy k, jako i jest już zwrócony odpowiedni ciąg binarny, i yield i kończy działanie iteratora.:

it_bin = bin_pow_set(3)
print(list(it_bin)) # -> ['000', '001', '010', '011', '100', '101', '110', '111']

Tylko troszkę bardziej skomplikowane będzie uogólnienie tego na dowolne ciagi n wyrazowe z m elementów (wszystkie podzbiory o długości pięć z czterech elementów). Robimy to tak, że zamiast liczb binarnych bierzemy dowolna bazę i robimy praktycznie to samo, co wyżej, kod (równiez iterator):

def any_pow_set(n, m):
    """returns the all length n sequences of
    m elements, as an iterator"""
    def append_zeroes(x, n, base):
        """helper function to append zeroes"""
        x = dec_to_any(x, base)
        return '0' * (n - len(x)) + x
    i = 0
    k = 0
    yield '0' * n
    while k < m ** n - 1:
        k += 1
        i = append_zeroes(k, n, m)
        yield i

Warunek się zmienia, w binarnym było 2^n - 1, tu jest m^n, i jeszcze funkcja pomocnicza konwertująca liczbę z notacji dziesietnej na dowolną (tu, dla prostoty, ograniczona do szesnatkowego):

def dec_to_any(n,base):
    conv_string = "0123456789ABCDEF"
    if n < base:
       return conv_string[n]
    else:
        return dec_to_any(n // base, base) + conv_string[n % base]

Przykład:

it_any = any_pow_set(3, 3)
print(list(it_any)) # -> ['000', '001', '002', '010', '011', '012', ... , '220', '221', '222']
Gray Code

Można iść dalej, gdyby generować podzbiory binarne w takiej kolejności, że w każdym następnym rózni sie tylko jeden bit od poprzedniego to mamy Gray Code, o którego zastosowaniach można poczytać na Wikipedii.
Najprostszy i najmniej wydajny jest sposób rekurencyjny, idea:
Startując od pustego stringa (G[0]), kolejne ciągi G[n+1] tworzą się tak(równanie rekurencyjne):
G[0] = ""
G[n+1] = 0prefG[n], 1pref_rG[n]
gdzie Opref to operacja prefiksujaca zerem każdy ciąg, a 1pref_r - odwraca i prefiksuje jedynką każdy ciąg. I znowu przetłumaczenie tego na kod w Pythonie jest, praktycznie bezpośrednie:

def gray(n, xs=[""]):
    """recursively computes a list of gray codes of
    length n"""
    def zero_pref(xs):
        if not xs:
            return ["0"]
        return ["0" + x for x in xs]

    def one_pref_rev(xs):
        if not xs:
            return ["1"]
        else:
            return ["1" + x for x in [y for y in xs]][::-1]

    if n == 0:
        return ""
    else:
        return zero_pref(gray(n - 1, xs)) + one_pref_rev(gray(n - 1, xs))

Funkcje pomocnicze, przy użyciu list comprehension, robią co powinny, po czym następuje warunek stopu: if n == 0: return "" i, zgodne ze specyfikacją wywołanie rekurencyjne: zero_pref(gray(n - 1, xs)) + one_pref_rev(gray(n - 1, xs))

Wydruk:

print(gray(3)) # -> ['000', '001', '011', '010', '110', '111', '101', '100']

Mozna to troche usprawnić, pozbyć się rekurencji, wersja iteracyjna:

def gray_iter(n):
    """iteratively computes a list of gray codes of
        length n"""
    def zero_pref(xs):
        if not xs:
            return ["0"]
        return ["0" + x for x in xs]

    def one_pref_rev(xs):
        if not xs:
            return ["1"]
        else:
            return ["1" + x for x in [y for y in xs]][::-1]

    xs= [""]
    k = 0
    while k < n:
        xs = zero_pref(xs) + one_pref_rev(xs)
        k += 1
    return xs

To tyle, dziękuję za uwagę, kod na githubie, a challenge na dziś, zmienić Gray Code na strumień (iterator):) #Python #Theory #python #Numerical

lion137
2017-12-13 01:00

Co Mnie Naszło Przeglądając Notatki

Znowu hashtag teoria:) Patrząc na hierarchię Chomskiego:

  • Języki Regularne (wyrażenia regularne);
  • Context Free;
  • Context Sensitive;
  • Języki stopnia zero, nazwijmy je Turing-equivalent systems.

Można jeszcze dołożyć kolejny poziom: Nieobliczalne - rozwiązujące problemy, dla których jeszcze nie powstał język, jak parsowanie Perla:-D, czy problem stopu

Możemy to przedstawić praktycznie:

  • a* - regular
  • a{n}b{n} - context free, rozpoznaje domknięcie nawiasów
  • a{n}b{n}c{n} - context sensitive, domknięcie nawiasów i jeszcze więcej
  • Turing-equivalent systems - generalnie programowanie i obliczenia, jak: czy ta lista zawiera wszystkie liczby pierwsze mniejsze od 99999?
  • Undecidable: czy dana funkcja zwraca stałą???

Rozważmy języki rozpoznawalne przez regexy:

  • 'r' - jakiś regex;
  • 'r' 'r' - konkatenacja dwóch regexów;
  • 'r' | 'r' - regex lub regex;
  • 'r' | 'ε' - jakiś regex lub pusty string;
  • poprzednią regułę mozna już uprościć do tzw. Kleene Star: r* - jakaś ilość regex'a 'r' lub znak pusty; więc: a rozpoznaje: 'a', 'ε', 'aa', 'aaaa', i tak dalej.

Czyli języki rozpoznawane przez języki stopnia 2 (liczymy od zera:)) to języki regularne. A jaką gramatyką rozpoznać język regiularny? Regularną? Hm...

Skonstruujmy gramatykę rozpoznającą regexy:

  • REGEX = REGEX'*' - jeśli regex, to regex 'Kleene Star' to też regex
  • REGEX = '(' REGEX ')' - regex "otoczony" nawiasami jest regex'em
  • REGEX = REGEX '|' REGEX - regex "lub" regex to regex
  • REGEX = REGEX REGEX - konkatenacja.

Patrząc na te reguły widzimy, że druga reguła nie jest regularna, gdyż musi rozpoznać dowolnie zgnieżdżone nawiasy - czyli jest context free. Wniosek: To co jest rozpoznawalne przez języki nie ma zwiazku z tym co one [języki] rozpoznają, jak C - rozpoznaje context free, ale parser C jest context sensitive. Raczej mało intuicjne!

Ale idźmy dalej, spróbujmy skonstruować gramatykę (język) do rozpoznawania języków:

  • GRAMMAR = DEF
  • GRAMMAR = DEF GRAMMAR
  • DEF = SYMBOL '=' RHS '\n'

Nie będę dalej męczył:) to nie jest cała gramatyka, ale okazuje się, że język do rozpoznawania języków jest regularny!!! Czyli teza udwodniona: Nie ma intuicji są tylko dowody!!!

EDIT Patrząc na ten wpis widzę, że może nie być zbyt zrozumiały:), więc w razie wątpliwości zachęcam do pytań w komentarzach

#computation #theory