Teorija brojeva u SAGE-u
verzija: SageMath 8.4
In [1]:
%display ascii_art
Skupovi brojeva
In [2]:
[NN,ZZ,QQ,RR,CC]
Out[2]:
In [3]:
[24 in NN, 24 in ZZ, -32 in NN, -32 in ZZ]
Out[3]:
In [4]:
[2/3 in ZZ, 2/3 in QQ, 2.3 in QQ, 2.3 in RR, sqrt(2) in QQ, sqrt(2) in RR]
Out[4]:
In [5]:
[sqrt(-2) in RR, sqrt(-2) in CC, I in CC]
Out[5]:
In [6]:
type(NN)
Out[6]:
In [7]:
type(CC)
Out[7]:
Cjelobrojne funkcije "pod" i "strop"
In [8]:
floor(5.23)
Out[8]:
In [9]:
ceil(5.23)
Out[9]:
In [10]:
floor(pi)
Out[10]:
In [11]:
ceil(pi)
Out[11]:
In [12]:
floor(-5.23)
Out[12]:
In [13]:
ceil(-5.23)
Out[13]:
In [14]:
plot(floor(x),(x,-3,3),thickness=3,figsize=[4,4])
Out[14]:
In [15]:
plot(floor(x),(x,-3,3),thickness=3,exclude=[-3..3],figsize=[4,4])
Out[15]:
In [16]:
plot(ceil(x),(x,-3,3),thickness=3,figsize=[4,4])
Out[16]:
In [17]:
plot(ceil(x),(x,-3,3),thickness=3,exclude=[-3..3],figsize=[4,4])
Out[17]:
In [18]:
slika1=plot(floor(x),(x,-3,3),color="blue",thickness=3,exclude=[-3..3])
slika2=plot(ceil(x),(x,-3,3),color="red",thickness=3,exclude=[-3..3])
show(slika1+slika2,figsize=[4,4])
Teorem o dijeljenju s ostatkom
In [19]:
500//17
Out[19]:
In [20]:
500%17
Out[20]:
In [21]:
-500//17
Out[21]:
In [22]:
-500%17
Out[22]:
Neočekivani rezultati
In [23]:
500//-17
Out[23]:
In [24]:
500%-17
Out[24]:
Definiramo svoje funkcije koje će raditi prema našim definicijama
In [25]:
def kvocijent(a,b):
if b==0:
return "Error: dijeljenje s nulom"
elif b>0:
return floor(a/b)
else:
return ceil(a/b)
In [26]:
def ostatak(a,b):
if b==0:
return "Error: dijeljenje s nulom"
elif b>0:
return a-floor(a/b)*b
else:
return a-ceil(a/b)*b
In [27]:
kvocijent(500,17)
Out[27]:
In [28]:
ostatak(500,17)
Out[28]:
In [29]:
kvocijent(500,-17)
Out[29]:
In [30]:
ostatak(500,-17)
Out[30]:
In [31]:
ostatak(-50,0)
Out[31]:
Najveća zajednička mjera
In [32]:
gcd(30,18)
Out[32]:
In [33]:
gcd(-30,18)
Out[33]:
In [34]:
gcd(-15,-17)
Out[34]:
In [35]:
gcd(252,200)
Out[35]:
In [36]:
gcd([15,30,48])
Out[36]:
Jedno cjelobrojno rješenje jednadžbe 252x+200y=M(252,200)
In [37]:
xgcd(252,200)
Out[37]:
Najmanji zajednički višekratnik
In [38]:
lcm(8,12)
Out[38]:
In [39]:
lcm(-8,12)
Out[39]:
In [40]:
lcm([12,5,14])
Out[40]:
In [41]:
lcm(252,200)
Out[41]:
Prosti brojevi
Ispitivanje da li je broj prost ili složen
In [42]:
is_prime(1023)
Out[42]:
In [43]:
is_prime(17)
Out[43]:
In [44]:
map(lambda x: is_prime(x),[1023,17])
Out[44]:
Da li je prirodni broj potencija prostog broja
In [45]:
is_prime_power(64)
Out[45]:
In [46]:
is_prime_power(111934216)
Out[46]:
Rastav broja na proste faktore
In [47]:
%display latex
In [48]:
factor(40)
Out[48]:
In [49]:
factor(4063)
Out[49]:
In [50]:
factor(10468)
Out[50]:
In [51]:
factor(10^30+1)
Out[51]:
Djeljitelji prirodnog broja
In [52]:
%display ascii_art
In [53]:
divisors(40)
Out[53]:
In [54]:
divisors(4063)
Out[54]:
In [55]:
divisors(10468)
Out[55]:
Prosti djeljitelji prirodnog broja
In [56]:
[faktor for faktor in divisors(40) if is_prime(faktor)]
Out[56]:
In [57]:
[faktor for faktor in divisors(4063) if is_prime(faktor)]
Out[57]:
In [58]:
[faktor for faktor in divisors(10468) if is_prime(faktor)]
Out[58]:
In [59]:
[faktor for faktor in divisors(10^30+1) if is_prime(faktor)]
Out[59]:
Skup prostih brojeva
In [60]:
Primes()
Out[60]:
In [61]:
23 in Primes()
Out[61]:
In [62]:
1 in Primes()
Out[62]:
Prvih 30 prostih brojeva
In [63]:
primes_first_n(30)
Out[63]:
In [64]:
list_plot(primes_first_n(30),figsize=[5,5])
Out[64]:
In [65]:
list_plot(primes_first_n(30),plotjoined=True,figsize=[5,5])
Out[65]:
In [66]:
list_plot(primes_first_n(100),color="red",figsize=[5,5])
Out[66]:
In [67]:
list_plot(primes_first_n(100),plotjoined=True,color="red",figsize=[5,5])
Out[67]:
Deseti po redu prosti broj
In [68]:
primes_first_n(10)[-1]
Out[68]:
Skup prostih brojeva većih ili jednakih od 41 i strogo manjih od 200
In [69]:
prime_range(41,200)
Out[69]:
Skup svih brojeva između 41 i 150 koji su potencije prostih brojeva
In [70]:
prime_powers(41,150)
Out[70]:
Najmanji prosti broj veći od n
In [71]:
next_prime(100)
Out[71]:
In [72]:
next_prime(31)
Out[72]:
In [73]:
next_prime(44)
Out[73]:
In [74]:
next_prime(2^100)
Out[74]:
Najveći prosti broj manji od n
In [75]:
previous_prime(100)
Out[75]:
In [76]:
previous_prime(31)
Out[76]:
In [77]:
previous_prime(44)
Out[77]:
In [78]:
previous_prime(2^100)
Out[78]:
Funkcija π(x)
In [79]:
prime_pi(100)
Out[79]:
In [80]:
prime_pi(10000)
Out[80]:
In [81]:
prime_pi(34.56)
Out[81]:
In [82]:
prime_pi(pi^100)
In [83]:
floor(pi^100).ndigits(10)
Out[83]:
In [84]:
n(pi^100)
Out[84]:
In [85]:
prime_pi(pi^20)
Out[85]:
In [86]:
[prime_pi(n) for n in xrange(1,50)]
Out[86]:
In [87]:
[(n,prime_pi(n)) for n in xrange(1,50)]
Out[87]:
Teorem o prostim brojevima
In [88]:
prvi=list_plot([prime_pi(n) for n in xrange(1,101)])
show(prvi,figsize=[7,5])
In [89]:
drugi=plot(x/log(x),(x,1,100),color="red")
show(drugi,figsize=[7,5])
In [90]:
show(prvi+drugi,figsize=[7,5])
In [91]:
def teorem_prosti_brojevi(n,xy=[7,5]):
slika1=list_plot([prime_pi(k) for k in xrange(1,n)])
slika2=plot(x/log(x),(x,1,n),color="red")
slika3=show(slika1+slika2,figsize=xy)
return slika3
In [92]:
teorem_prosti_brojevi(100)
In [93]:
teorem_prosti_brojevi(1000)
In [94]:
teorem_prosti_brojevi(10000)
In [95]:
teorem_prosti_brojevi(100000)
In [96]:
[2e4,8e4,1e5]
Out[96]:
Slučajni prosti brojevi
slučajni prosti broj manji od 200
In [97]:
random_prime(200)
Out[97]:
In [98]:
random_prime(200)
Out[98]:
lista od 10 slučajnih prostih brojeva manjih od 200
In [99]:
[random_prime(200) for n in xrange(10)]
Out[99]:
slučajni prosti broj između 200 i 500
In [100]:
random_prime(500,lbound=200)
Out[100]:
slučajni prosti broj sa 12 znamenaka
In [101]:
random_prime(10^12,lbound=10^11)
Out[101]:
3x4 matrica slučajnih prostih brojeva manjih od 1000
In [102]:
[[random_prime(1000) for n in xrange(3)] for m in xrange(4)]
Out[102]:
In [103]:
matrix([[random_prime(1000) for n in xrange(3)] for m in xrange(4)])
Out[103]:
dva slučajna prosta broja za RSA kriptosustav
In [104]:
(p,q)=(random_prime(2^640,lbound=2^639),random_prime(2^640,lbound=2^639))
In [105]:
p
Out[105]:
In [106]:
q
Out[106]:
brojevi p i q u bazi 10 imaju 193 znamenke, a u bazi 2 imaju 640 znamenaka (bitova)
In [107]:
p.ndigits(10)
Out[107]:
In [108]:
p.ndigits(2)
Out[108]:
In [109]:
[q.ndigits(10),q.ndigits(2)]
Out[109]:
lista znamenaka broja p u bazi 2 (pazite, znamenke su navedene u obrnutom redoslijedu)
In [110]:
p.digits(2)
Out[110]:
In [111]:
broj=4
(broj.digits(2),list(reversed(broj.digits(2))))
Out[111]:
opcija "proof=False" daje pseudoprosti broj
In [112]:
time random_prime(2^640,lbound=2^639)
Out[112]:
In [113]:
time random_prime(2^640,lbound=2^639,proof=False)
Out[113]:
In [114]:
time random_prime(2^640,lbound=2^639,proof=True)
Out[114]:
Linearne kongruencije
Multiplikativni inverz modulo n
$15x\equiv 1\!\pmod{341}$ rješenje: $x=91$
In [115]:
inverse_mod(15,341)
Out[115]:
$10x\equiv 1\!\pmod{12}$ nema rješenja
In [116]:
inverse_mod(10,12)
za rješavanje bilo koje linearne kongruencije možemo lagano napisati svoju funkciju
In [117]:
def lin_kong(a,b,n):
d=gcd(a,n)
if b%d!=0:
return False
a1=a//d
b1=b//d
n1=n//d
x=(inverse_mod(a1,n1)*b1)%n1
rjesenja=[x+k*n1 for k in xrange(d)]
return rjesenja
$79x\equiv 15\!\pmod{722}$ rješenje: $119$
In [118]:
lin_kong(79,15,722)
Out[118]:
$155x\equiv 1185\!\pmod{1405}$ rješenja: $198, 479, 760, 1041, 1322$
In [119]:
lin_kong(155,1185,1405)
Out[119]:
$15x\equiv 1\!\pmod{341}$ rješenje: $91$
In [120]:
lin_kong(15,1,341)
Out[120]:
$10x\equiv 1\!\pmod{12}$ nema rješenja
In [121]:
lin_kong(10,1,12)
Out[121]:
$65x\equiv 27\!\pmod{169}$ nema rješenja
In [122]:
lin_kong(65,27,169)
Out[122]:
Funkcija koja daje rješenja linearne kongruencije $ax\equiv b\!\pmod{n}$ zajedno s tablicom kakvu dobivamo ručnim rješavanjem.
In [123]:
def lin_kong_rucno(a,b,n):
qi=[]
yi=[0,1]
n1,q,r1,d=a,n//a,n%a,a
while r1!=0:
qi.append(q)
yi.append(yi[-2]-q*yi[-1])
q=n1//r1
n1,r1=r1,n1%r1
if r1!=0: d=r1
if b%d!=0: return False
a1,b1,n1=a//d,b//d,n//d
x1=(yi[-1]*b1)%n1
rjesenja=[x1+k*n1 for k in xrange(d)]
tablica=table([["$i$"]+range(-1,len(yi)-1)+["rjesenja"],["$q_i$"," "," "]+qi+[rjesenja],["$y_i$"]+yi+
["$(a_1,b_1,n_1)=(%d,%d,%d)$"%(a1,b1,n1)]])
return tablica
$79x\equiv 15\!\pmod{722}$ rješenje: $119$
In [124]:
lin_kong_rucno(79,15,722)
Out[124]:
$155x\equiv 1185\!\pmod{1405}$ rješenja: $198, 479, 760, 1041, 1322$
In [125]:
lin_kong_rucno(155,1185,1405)
Out[125]:
$10x\equiv 1\!\pmod{12}$ nema rješenja
In [126]:
lin_kong_rucno(10,1,12)
Out[126]:
Kineski teorem o ostacima
$x\equiv 1\!\pmod{4},\quad x\equiv 8\!\pmod{9},\quad x\equiv 10\!\pmod{25}$ rješenje: $x=485$
In [127]:
crt([1,8,10],[4,9,25])
Out[127]:
$x\equiv 4\!\pmod{5},\quad x\equiv 1\!\pmod{12},\quad x\equiv 7\!\pmod{14}$ rješenje: $x=49$
In [128]:
crt([4,1,7],[5,12,14])
Out[128]:
$x\equiv 2\!\pmod{3},\quad x\equiv 1\!\pmod{12},\quad x\equiv 7\!\pmod{14}$ nema rješenja
In [129]:
crt([2,1,7],[3,12,14])
Funkcija koja daje rješenje sustava kongruencija pomoću kineskog teorema o ostacima na način kojeg provodimo ručnim rješavanjem. Na izlazu se vraća rječnik u kojemu su dane vrijednosti za $n, k_i, x_i, x_0, x_{0,mod}$.
In [130]:
def KTO_rucno(A,N):
"""Daje rjesenje sustava kongruencija x=a_i(mod n_i) pomocu kineskog teorema o ostacima
na nacin kojeg provodimo rucnim rjesavanjem.
Moduli n_i moraju biti u parovima relativno prosti. Na izlazu se vraca rjecnik u kojemu su
dane vrijednosti za n, k_i, x_i, x_0 i x0_mod.
A je lista a_i-ova, a N je lista n_i-ova."""
for t in Combinations(N,2).list():
if gcd(t[0],t[1])!=1:
return "Error: moduli nisu u parovima relativno prosti."
n=prod(N)
ki=map(lambda x: n/x,N)
xi=[]
for i in xrange(len(N)):
xi.append(int(solve_mod([ki[i]*x==A[i]],N[i])[0][0]))
x0=sum(map(lambda x: x[0]*x[1],zip(ki,xi)))
x0mod=x0%n
izlaz=dict(zip(['n','k_i','x_i','x0','x0_mod'],[n,ki,xi,x0,x0mod]))
return izlaz
$x\equiv 1\!\pmod{4},\quad x\equiv 8\!\pmod{9},\quad x\equiv 10\!\pmod{25}$ rješenje: $x=485$
In [131]:
KTO_rucno([1,8,10],[4,9,25])
Out[131]:
$x\equiv 4\!\pmod{5},\quad x\equiv 1\!\pmod{12},\quad x\equiv 7\!\pmod{14}$ rješenje: $x=49$
In [132]:
KTO_rucno([4,1,7],[5,12,14])
Out[132]:
Eulerova funkcija
In [133]:
euler_phi(17)
Out[133]:
In [134]:
euler_phi(28)
Out[134]:
In [135]:
[euler_phi(n) for n in xrange(1,51)]
Out[135]:
In [136]:
list_plot([euler_phi(n) for n in xrange(1,51)],xmax=55,figsize=[5,4])
Out[136]:
Kako dobiti listu svih pozitivnih brojeva manjih od 100 koji su relativno prosti sa 100?
In [137]:
filter(lambda n: gcd(100,n)==1,range(100))
Out[137]:
In [138]:
def reducirani_sustav_ostataka(m):
ostaci=filter(lambda n: gcd(m,n)==1,range(m))
return ostaci
In [139]:
reducirani_sustav_ostataka(100)
Out[139]:
In [140]:
reducirani_sustav_ostataka(1000)
Out[140]:
Posljednje dvije znamenke broja $3^{502}\cdot 7^{201}$
In [141]:
3^502*7^201%100
Out[141]:
Za velike potencije naredba power_mod puno brže računa ostatak
In [142]:
time 45^33333333%342
Out[142]:
In [143]:
time power_mod(45,33333333,342)
Out[143]:
Sustavi kongruencija
In [144]:
var('x,y')
Out[144]:
$x+5y\equiv 3\!\pmod{8},\quad 4x+5y\equiv 1\!\pmod{8}$ rješenje: $(2,5)$
In [145]:
solve_mod([x+5*y==3,4*x+5*y==1],8)
Out[145]:
$3x+5y\equiv 14\!\pmod{28},\quad 5x+9y\equiv 6\!\pmod{28}$ rješenja: $(20,2),\,\,(6,16)$
In [146]:
solve_mod([3*x+5*y==14,5*x+9*y==6],28)
Out[146]:
$2x+3y\equiv 6\!\pmod{7},\quad 3x+y\equiv 2\!\pmod{7}$ rješenja: $(0,2),\,\,(1,6),\,\,(2,3),\,\,(3,0),\,\,(4,4),\,\,(5,1),\,\,(6,5)$
In [147]:
solve_mod([2*x+3*y==6,3*x+y==2],7)
Out[147]:
$2x+3y\equiv 6\!\pmod{5\cdot6\cdot7\cdot11},\quad 3x+y\equiv 2\!\pmod{5\cdot6\cdot7\cdot11}$
In [148]:
solve_mod([2*x+3*y==6,3*x+y==2],5*6*7*11)
Out[148]:
$3x+5y\equiv 14\!\pmod{100},\quad 5x+9y\equiv 6\!\pmod{100}$ rješenja: $(48,74),\,\,(98,24)$
In [149]:
solve_mod([3*x+5*y==14,5*x+9*y==6],100)
Out[149]:
$79x\equiv 15\!\pmod{722}$ rješenje: $119$
In [150]:
solve_mod([79*x==15],722)
Out[150]:
Red od a modulo n
red od 2 modulo 7 je jednak 3
In [151]:
Mod(2,7).multiplicative_order()
Out[151]:
In [152]:
type(Mod(2,7))
Out[152]:
In [153]:
type(mod(2,7))
Out[153]:
red od 16 modulo 29 je jednak 7
In [154]:
Mod(16,29).multiplicative_order()
Out[154]:
red od 12 modulo 28 ne postoji
In [155]:
Mod(12,28).multiplicative_order()
In [156]:
Mod(2,7) in IntegerModRing(7)
Out[156]:
In [157]:
mod(2,7) in IntegerModRing(7)
Out[157]:
Primitivni korijen modulo n
primitivni korijen modulo 7
In [158]:
primitive_root(7)
Out[158]:
svi primitivni korijeni modulo 7
In [159]:
filter(lambda x: Mod(x,7).multiplicative_order()==euler_phi(7), range(1,7))
Out[159]:
svi primitivni korijeni modulo 41
In [160]:
filter(lambda x: Mod(x,41).multiplicative_order()==euler_phi(41), range(1,41))
Out[160]:
ne postoji primitivni korijen modulo 12
In [161]:
filter(lambda x: gcd(x,12)==1 and Mod(x,12).multiplicative_order()==euler_phi(12), range(1,12))
Out[161]:
In [162]:
primitive_root(12)
primitivni korijen modulo n generira reducirani sustav ostataka modulo n
$n=41$
In [163]:
reducirani_sustav_ostataka(41)
Out[163]:
In [164]:
[mod(6^k,41) for k in xrange(euler_phi(41))]
Out[164]:
In [165]:
sorted([mod(6^k,41) for k in xrange(euler_phi(41))])
Out[165]:
$n=50$
In [166]:
primitive_root(50)
Out[166]:
In [167]:
reducirani_sustav_ostataka(50)
Out[167]:
In [168]:
sorted([mod(27^k,50) for k in xrange(euler_phi(50))])
Out[168]:
Nelinearne kongruencije
$x^5\equiv 2\!\pmod{7}$
In [169]:
solve_mod([x^5==2],7)
Out[169]:
$x^{12}\equiv 37\!\pmod{41}$
In [170]:
solve_mod([x^12==37],41)
Out[170]:
$7^x\equiv 6\!\pmod{17}$
In [171]:
filter(lambda x: mod(7^x,17)==6, range(17))
Out[171]:
$3^x\equiv 2\!\pmod{23}$
In [172]:
filter(lambda x: mod(3^x,23)==2, range(23))
Out[172]:
$2x^8\equiv 5\!\pmod{13}$
In [173]:
solve_mod([2*x^8==5],13)
Out[173]:
$x^5\equiv 2\!\pmod{12}$ nema rješenja
In [174]:
solve_mod([x^5==2],12)
Out[174]:
$x^5\equiv 4\!\pmod{12}$
In [175]:
solve_mod([x^5==4],12)
Out[175]:
In [176]:
filter(lambda x: mod(x^5,12)==4, range(12))
Out[176]:
$x^7\equiv 17\!\pmod{322122548}$
In [177]:
time solve_mod([x^7==17],322122548)
Out[177]:
In [178]:
%display latex
In [179]:
factor(322122548)
Out[179]:
Još jedna napomena
Ovdje želimo naglasiti razliku između naredbi mod (ili Mod) i naredbe % (za dobivanje ostatka pri dijeljenju dva cijela broja). Na prvi pogled se čini da te naredbe rade istu stvar, ali nije baš tako. Doduše, na izlazu daju kao isti rezultate, ali te rezultate SAGE potpuno drukčije doživljava. Pogledajmo na jednom primjeru tu razliku.
In [180]:
%display ascii_art
In [181]:
37%14
Out[181]:
In [182]:
type(37%14)
Out[182]:
Kako na 37%14 SAGE gleda kao na element iz prstena $\mathbb{Z}$, zbrajanje i množenje se odvija u prstenu $\mathbb{Z}$ i dobivamo očekivani rezultat.
In [183]:
37%14+5*14
Out[183]:
Da li je tako i s naredbom mod?
In [184]:
mod(37,14)
Out[184]:
In [185]:
type(mod(37,14))
Out[185]:
Kako na mod(37,14) SAGE gleda kao na element iz prstena $\mathbb{Z}_{14}$, tada se zbrajanje i množenje odvijaju u prstenu $\mathbb{Z}_{14}$. Stoga rezultat kojeg je vratio SAGE nije bug, već je to stvarno dobar rezultat. Naime, u prstenu $\mathbb{Z}_{14}$ svi višekratnici broja 14 su zapravo jednaki 0, pa je tako i $5\cdot14$ također jednako 0 jer se operacije zbrajanja i množenja u prstenu $\mathbb{Z}_{14}$ rade modulo 14.
In [186]:
mod(37,14)+5*14
Out[186]:
In [ ]: