|
Size: 13084
Comment:
|
Size: 30457
Comment:
|
| Deletions are marked like this. | Additions are marked like this. |
| Line 1: | Line 1: |
| = Sage Interactions - Number Theory = goto [:interact:interact main page] [[TableOfContents]] |
<<TableOfContents>> = Integer Factorization = == Divisibility Poset == by William Stein {{{#!sagecell @interact def _(n=(5..100)): Poset(([1..n], lambda x, y: y%x == 0) ).show() }}} {{attachment:divposet.png}} |
| Line 8: | Line 18: |
| {{{ | {{{#!sagecell |
| Line 53: | Line 63: |
| attachment:factortree.png === Continued Fraction Plotter === |
{{attachment:factortree.png}} More complicated demonstration using Mathematica: http://demonstrations.wolfram.com/FactorTrees/ == Factoring an Integer == by Timothy Clemans Sage implementation of the Mathematica demonstration of the same name. http://demonstrations.wolfram.com/FactoringAnInteger/ {{{#!sagecell @interact def _(r=selector(range(0,10000,1000), label='range', buttons=True), n=slider(0,1000,1,2,'n',False)): if not r and n in (0, 1): n = 2 s = '$%d = %s$' % (r + n, factor(r + n)) s = s.replace('*', '\\times') html(s) }}} = Prime Numbers = == Illustrating the prime number theorem == |
| Line 57: | Line 86: |
| {{{ @interact def _(number=e, ymax=selector([None,5,20,..,400],nrows=2), clr=Color('purple'), prec=[500,1000,..,5000]): c = list(continued_fraction(RealField(prec)(number))); print c show(line([(i,z) for i, z in enumerate(c)],rgbcolor=clr),ymax=ymax,figsize=[10,2]) }}} attachment:contfracplot.png === Illustrating the prime number thoerem === by William Stein {{{ |
{{{#!sagecell |
| Line 73: | Line 92: |
| attachment:primes.png === Computing Generalized Bernoulli Numbers === by William Stein (Sage-2.10.3) {{{ @interact def _(m=selector([1..15],nrows=2), n=(7,(3..10))): G = DirichletGroup(m) s = "<h3>First n=%s Bernoulli numbers attached to characters with modulus m=%s</h3>"%(n,m) s += '<table border=1>' s += '<tr bgcolor="#edcc9c"><td align=center>$\\chi$</td><td>Conductor</td>' + \ ''.join('<td>$B_{%s,\chi}$</td>'%k for k in [1..n]) + '</tr>' for eps in G.list(): v = ''.join(['<td align=center bgcolor="#efe5cd">$%s$</td>'%latex(eps.bernoulli(k)) for k in [1..n]]) s += '<tr><td bgcolor="#edcc9c">%s</td><td bgcolor="#efe5cd" align=center>%s</td>%s</tr>\n'%( eps, eps.conductor(), v) s += '</table>' html(s) }}} attachment:bernoulli.png === Fundamental Domains of SL_2(ZZ) === by Robert Miller {{{ L = [[-0.5, 2.0^(x/100.0) - 1 + sqrt(3.0)/2] for x in xrange(1000, -1, -1)] R = [[0.5, 2.0^(x/100.0) - 1 + sqrt(3.0)/2] for x in xrange(1000)] xes = [x/1000.0 for x in xrange(-500,501,1)] M = [[x,abs(sqrt(x^2-1))] for x in xes] fundamental_domain = L+M+R fundamental_domain = [[x-1,y] for x,y in fundamental_domain] @interact def _(gen = selector(['t+1', 't-1', '-1/t'], nrows=1)): global fundamental_domain if gen == 't+1': fundamental_domain = [[x+1,y] for x,y in fundamental_domain] elif gen == 't-1': fundamental_domain = [[x-1,y] for x,y in fundamental_domain] elif gen == '-1/t': new_dom = [] for x,y in fundamental_domain: sq_mod = x^2 + y^2 new_dom.append([(-1)*x/sq_mod, y/sq_mod]) fundamental_domain = new_dom P = polygon(fundamental_domain) P.ymax(1.2); P.ymin(-0.1) P.show() }}} attachment:fund_domain.png === Computing modular forms === by William Stein {{{ j = 0 @interact def _(N=[1..100], k=selector([2,4,..,12],nrows=1), prec=(3..40), group=[(Gamma0, 'Gamma0'), (Gamma1, 'Gamma1')]): M = CuspForms(group(N),k) print j; global j; j += 1 print M; print '\n'*3 print "Computing basis...\n\n" if M.dimension() == 0: print "Space has dimension 0" else: prec = max(prec, M.dimension()+1) for f in M.basis(): view(f.q_expansion(prec)) print "\n\n\nDone computing basis." }}} attachment:modformbasis.png === Computing the cuspidal subgroup === by William Stein {{{ html('<h1>Cuspidal Subgroups of Modular Jacobians J0(N)</h1>') @interact def _(N=selector([1..8*13], ncols=8, width=10, default=10)): A = J0(N) print A.cuspidal_subgroup() }}} attachment:cuspgroup.png === A Charpoly and Hecke Operator Graph === by William Stein {{{ # Note -- in Sage-2.10.3; multiedges are missing in plots; loops are missing in 3d plots @interact def f(N = prime_range(11,400), p = selector(prime_range(2,12),nrows=1), three_d = ("Three Dimensional", False)): S = SupersingularModule(N) T = S.hecke_matrix(p) G = Graph(T, multiedges=True, loops=not three_d) html("<h1>Charpoly and Hecke Graph: Level %s, T_%s</h1>"%(N,p)) show(T.charpoly().factor()) if three_d: show(G.plot3d(), aspect_ratio=[1,1,1]) else: show(G.plot(),figsize=7) }}} attachment:heckegraph.png === Demonstrating the Diffie-Hellman Key Exchange Protocol === by Timothy Clemans (refereed by William Stein) {{{ @interact def diffie_hellman(button=selector(["New example"],label='',buttons=True), bits=("Number of bits of prime", (8,12,..512))): maxp = 2^bits p = random_prime(maxp) k = GF(p) if bits>100: g = k(2) else: g = k.multiplicative_generator() a = ZZ.random_element(10, maxp) b = ZZ.random_element(10, maxp) print """ <html> <style> .gamodp { background:yellow } .gbmodp { background:orange } .dhsame { color:green; font-weight:bold } </style> <h2>%s-Bit Diffie-Hellman Key Exchange</h2> <ol style="color:#000;font:12px Arial, Helvetica, sans-serif"> <li>Alice and Bob agree to use the prime number p=%s and base g=%s.</li> <li>Alice chooses the secret integer a=%s, then sends Bob (<span class="gamodp">g<sup>a</sup> mod p</span>):<br/>%s<sup>%s</sup> mod %s = <span class="gamodp">%s</span>.</li> <li>Bob chooses the secret integer b=%s, then sends Alice (<span class="gbmodp">g<sup>b</sup> mod p</span>):<br/>%s<sup>%s</sup> mod %s = <span class="gbmodp">%s</span>.</li> <li>Alice computes (<span class="gbmodp">g<sup>b</sup> mod p</span>)<sup>a</sup> mod p:<br/>%s<sup>%s</sup> mod %s = <span class="dhsame">%s</span>.</li> <li>Bob computes (<span class="gamodp">g<sup>a</sup> mod p</span>)<sup>b</sup> mod p:<br/>%s<sup>%s</sup> mod %s = <span class="dhsame">%s</span>.</li> </ol></html> """ % (bits, p, g, a, g, a, p, (g^a), b, g, b, p, (g^b), (g^b), a, p, (g^ b)^a, g^a, b, p, (g^a)^b) }}} attachment:dh.png === Plotting an elliptic curve over a finite field === {{{ E = EllipticCurve('37a') @interact def _(p=slider(prime_range(1000), default=389)): show(E) print "p = %s"%p show(E.change_ring(GF(p)).plot(),xmin=0,ymin=0) }}} attachment:ellffplot.png === Prime Spiral - Square === |
{{attachment:primes.png}} == Prime Spiral - Square FIXME == |
| Line 240: | Line 96: |
| {{{ | {{{#!sagecell |
| Line 372: | Line 228: |
{{attachment:SquareSpiral.PNG}} == Prime Spiral - Polar == by David Runde {{{#!sagecell @interact def polar_prime_spiral(start=1, end=2000, show_factors = false, highlight_primes = false, show_curves=true, n = 0): #For more information about the factors in the spiral, visit http://www.dcs.gla.ac.uk/~jhw/spirals/index.html by John Williamson. if start < 1 or end <=start: print "invalid start or end value" if n > end: print "WARNING: n is greater than end value" def f(n): return (sqrt(n)*cos(2*pi*sqrt(n)), sqrt(n)*sin(2*pi*sqrt(n))) list =[] list2=[] if show_factors == false: for i in [start..end]: if i.is_pseudoprime(): list.append(f(i-start+1)) #Primes list else: list2.append(f(i-start+1)) #Composites list P = points(list) R = points(list2, alpha = .1) #Faded Composites else: for i in [start..end]: list.append(disk((f(i-start+1)),0.05*pow(2,len(factor(i))-1), (0,2*pi))) #resizes each of the dots depending of the number of factors of each number if i.is_pseudoprime() and highlight_primes: list2.append(f(i-start+1)) P = plot(list) p_size = 5 #the orange dot size of the prime markers if not highlight_primes: list2 = [(f(n-start+1))] R=points(list2, hue = .1, pointsize = p_size) if n > 0: print 'n =', factor(n) p = 1 #The X which marks the given n W1 = disk((f(n-start+1)), p, (pi/6, 2*pi/6)) W2 = disk((f(n-start+1)), p, (4*pi/6, 5*pi/6)) W3 = disk((f(n-start+1)), p, (7*pi/6, 8*pi/6)) W4 = disk((f(n-start+1)), p, (10*pi/6, 11*pi/6)) Q = plot(W1+W2+W3+W4, alpha = .1) n=n-start+1 #offsets the n for different start values to ensure accurate plotting if show_curves: begin_curve = 0 t = var('t') a=1 b=0 if n > (floor(sqrt(n)))^2 and n <= (floor(sqrt(n)))^2 + floor(sqrt(n)): c = -((floor(sqrt(n)))^2 - n) c2= -((floor(sqrt(n)))^2 + floor(sqrt(n)) - n) else: c = -((ceil(sqrt(n)))^2 - n) c2= -((floor(sqrt(n)))^2 + floor(sqrt(n)) - n) print 'Pink Curve: n^2 +', c print 'Green Curve: n^2 + n +', c2 def g(m): return (a*m^2+b*m+c); def r(m) : return sqrt(g(m)) def theta(m) : return r(m)- m*sqrt(a) S1 = parametric_plot(((r(t))*cos(2*pi*(theta(t))),(r(t))*sin(2*pi*(theta(t)))), begin_curve, ceil(sqrt(end-start)), rgbcolor=hue(0.8), thickness = .2) #Pink Line b=1 c= c2; S2 = parametric_plot(((r(t))*cos(2*pi*(theta(t))),(r(t))*sin(2*pi*(theta(t)))), begin_curve, ceil(sqrt(end-start)), rgbcolor=hue(0.6), thickness = .2) #Green Line show(R+P+S1+S2+Q, aspect_ratio = 1, axes = false) else: show(R+P+Q, aspect_ratio = 1, axes = false) else: show(R+P, aspect_ratio = 1, axes = false) }}} {{attachment:PolarSpiral.PNG}} = Modular Forms = == Computing modular forms FIXME == by William Stein {{{#!sagecell j = 0 @interact def _(N=[1..100], k=selector([2,4,..,12],nrows=1), prec=(3..40), group=[(Gamma0, 'Gamma0'), (Gamma1, 'Gamma1')]): M = CuspForms(group(N),k) print j; global j; j += 1 print M; print '\n'*3 print "Computing basis...\n\n" if M.dimension() == 0: print "Space has dimension 0" else: prec = max(prec, M.dimension()+1) for f in M.basis(): view(f.q_expansion(prec)) print "\n\n\nDone computing basis." }}} {{attachment:modformbasis.png}} == Computing the cuspidal subgroup == by William Stein {{{#!sagecell html('<h1>Cuspidal Subgroups of Modular Jacobians J0(N)</h1>') @interact def _(N=selector([1..8*13], ncols=8, width=10, default=10)): A = J0(N) print A.cuspidal_subgroup() }}} {{attachment:cuspgroup.png}} == A Charpoly and Hecke Operator Graph FIXME == by William Stein {{{#!sagecell # Note -- in Sage-2.10.3; multiedges are missing in plots; loops are missing in 3d plots @interact def f(N = prime_range(11,400), p = selector(prime_range(2,12),nrows=1), three_d = ("Three Dimensional", False)): S = SupersingularModule(N) T = S.hecke_matrix(p) G = Graph(T, multiedges=True, loops=not three_d) html("<h1>Charpoly and Hecke Graph: Level %s, T_%s</h1>"%(N,p)) show(T.charpoly().factor()) if three_d: show(G.plot3d(), aspect_ratio=[1,1,1]) else: show(G.plot(),figsize=7) }}} {{attachment:heckegraph.png}} = Modular Arithmetic = == Quadratic Residue Table FIXME == by Emily Kirkman {{{#!sagecell from numpy import array as narray @interact def quad_res_plot(first_n_odd_primes = (20,200),display_size=[7..15]): # Compute numpy matrix of legendre symbols r = int(first_n_odd_primes) np = [nth_prime(i+2) for i in range(r)] leg = [[legendre_symbol(np[i], np[j]) for i in range(r)] for j in range(r)] na = narray(leg) for i in range(r): for j in range(r): if na[i][j] == 1 and Mod((np[i]-1)*(np[j]-1)//4,2) == 0: na[i][j] = 2 m = matrix(na) # Define plot structure MP = matrix_plot(m, cmap='Oranges') for i in range(r): if np[-1] < 100: MP += text('%d'%nth_prime(i+2),(-.75,r-i-.5), rgbcolor='black') MP += text('%d'%nth_prime(i+2), (i+.5,r+.5), rgbcolor='black') if len(np) < 75: MP += line([(0,i),(r,i)], rgbcolor='black') MP += line([(i,0),(i,r)], rgbcolor='black') if np[-1] < 100: for i in range(r): # rows for j in range(r): # cols if m[j][i] == 0: MP += text('0',(i+.5,r-j-.5),rgbcolor='black') elif m[j][i] == -1: MP += text('N',(i+.5,r-j-.5),rgbcolor='black') elif m[j][i] == 1: MP += text('A',(i+.5,r-j-.5),rgbcolor='black') elif m[j][i] == 2: MP += text('S',(i+.5,r-j-.5),rgbcolor='black') MP += line([(0,r),(r,r)], rgbcolor='black') MP += line([(r,0),(r,r)], rgbcolor='black') MP += line([(0,0),(r,0)], rgbcolor='black') MP += line([(0,0),(0,r)], rgbcolor='black') if len(np) < 75: MP += text('q',(r/2,r+2), rgbcolor='black', fontsize=15) MP += text('p',(-2.5,r/2), rgbcolor='black', fontsize=15) MP.show(axes=False, ymax=r, figsize=[display_size,display_size]) html('Symmetry of Prime Quadratic Residues mod the first %d odd primes.'%r) }}} {{attachment:quadres.png}} {{attachment:quadresbig.png}} == Cubic Residue Table FIXME == by Emily Kirkman {{{#!sagecell def power_residue_symbol(alpha, p, m): if p.divides(alpha): return 0 if not p.is_prime(): return prod(power_residue_symbol(alpha,ell,m)^e for ell, e in p.factor()) F = p.residue_field() N = p.norm() r = F(alpha)^((N-1)/m) k = p.number_field() for kr in k.roots_of_unity(): if r == F(kr): return kr def cubic_is_primary(n): g = n.gens_reduced()[0] a,b = g.polynomial().coefficients() return Mod(a,3)!=0 and Mod(b,3)==0 from numpy import array as narray @interact def cubic_sym(n=(10..35),display_size=[7..15]): # Compute numpy matrix of primary cubic residue symbols r = n m=3 D.<w> = NumberField(x^2+x+1) it = D.primes_of_degree_one_iter() pp = [] while len(pp) < r: k = it.next() if cubic_is_primary(k): pp.append(k) n = narray([ [ power_residue_symbol(pp[i].gens_reduced()[0], pp[j], m) \ for i in range(r) ] for j in range(r) ]) # Convert to integer matrix for gradient colors for i in range(r): for j in range(r): if n[i][j] == w: n[i][j] = int(-1) elif n[i][j] == w^2: n[i][j] = int(-2) elif n[i][j] == 1: n[i][j] = int(1) m = matrix(n) # Define plot structure MP = matrix_plot(m,cmap="Blues") for i in range(r): MP += line([(0,i),(r,i)], rgbcolor='black') MP += line([(i,0),(i,r)], rgbcolor='black') for j in range(r): if m[i][j] == -2: MP += text('$\omega^2$',(i+.5,r-j-.5),rgbcolor='black') if m[i][j] == -1: MP += text('$\omega $',(i+.5,r-j-.5),rgbcolor='black') if m[i][j] == 0: MP += text('0',(i+.5,r-j-.5),rgbcolor='black') if m[i][j] == 1: MP += text('R',(i+.5,r-j-.5),rgbcolor='white') MP += line([(0,r),(r,r)], rgbcolor='black') MP += line([(r,0),(r,r)], rgbcolor='black') MP += line([(0,0),(r,0)], rgbcolor='black') MP += line([(0,0),(0,r)], rgbcolor='black') MP += text('$ \pi_1$',(r/2,r+2), rgbcolor='black', fontsize=25) MP += text('$ \pi_2$',(-2.5,r/2), rgbcolor='black', fontsize=25) html('Symmetry of Primary Cubic Residues mod ' \ + '%d primary primes in $ \mathbf Z[\omega]$.'%r) MP.show(axes=False, figsize=[display_size,display_size]) }}} {{attachment:cubres.png}} = Cyclotomic Fields = == Gauss and Jacobi Sums in Complex Plane == by Emily Kirkman {{{#!sagecell def jacobi_sum(e,f): # If they are both trivial, return p if e.is_trivial() and f.is_trivial(): return (e.parent()).order() + 1 # If they are inverses of each other, return -e(-1) g = e*f if g.is_trivial(): return -e(-1) # If both are nontrivial, apply mult. formula: elif not e.is_trivial() and not f.is_trivial(): return e.gauss_sum()*f.gauss_sum()/g.gauss_sum() # If exactly one is trivial, return 0 else: return 0 def latex2(e): return latex(list(e.values_on_gens())) def jacobi_plot(p, e_index, f_index, with_text=True): # Set values G = DirichletGroup(p) e = G[e_index] f = G[f_index] ef = e*f js = jacobi_sum(e,f) e_gs = e.gauss_sum() f_gs = f.gauss_sum() ef_gs = (e*f).gauss_sum() # Compute complex coordinates f_pt = list(f.values_on_gens()[0].complex_embedding()) e_pt = list(e.values_on_gens()[0].complex_embedding()) ef_pt = list(ef.values_on_gens()[0].complex_embedding()) f_gs_pt = list(f_gs.complex_embedding()) e_gs_pt = list(e_gs.complex_embedding()) ef_gs_pt = list(ef_gs.complex_embedding()) try: js = int(js) js_pt = [CC(js)] except: js_pt = list(js.complex_embedding()) # Define plot structure S = circle((0,0),1,rgbcolor='yellow') S += line([e_pt,e_gs_pt], rgbcolor='red', thickness=4) S += line([f_pt,f_gs_pt], rgbcolor='blue', thickness=3) S += line([ef_pt,ef_gs_pt], rgbcolor='purple',thickness=2) S += point(e_pt,pointsize=50, rgbcolor='red') S += point(f_pt,pointsize=50, rgbcolor='blue') S += point(ef_pt,pointsize=50,rgbcolor='purple') S += point(f_gs_pt,pointsize=75, rgbcolor='black') S += point(e_gs_pt,pointsize=75, rgbcolor='black') S += point(ef_gs_pt,pointsize=75, rgbcolor='black') S += point(js_pt,pointsize=100,rgbcolor='green') if with_text: S += text('$J(%s,%s) = %s$'%(latex2(e),latex2(f),latex(js)), (3,2.5),fontsize=15, rgbcolor='black') else: html('$$J(%s,%s) = %s$$'%(latex2(e),latex2(f),latex(js))) return S @interact def single_jacobi_plot(p=prime_range(3,100), e_range=(0..100), f_range=(0..100)): e_index = floor((p-2)*e_range/100) f_index = floor((p-2)*f_range/100) S = jacobi_plot(p,e_index,f_index,with_text=False) S.show(aspect_ratio=1) }}} {{attachment:jacobising.png}} == Exhaustive Jacobi Plotter == by Emily Kirkman {{{#!sagecell def jacobi_sum(e,f): # If they are both trivial, return p if e.is_trivial() and f.is_trivial(): return (e.parent()).order() + 1 # If they are inverses of each other, return -e(-1) g = e*f if g.is_trivial(): return -e(-1) # If both are nontrivial, apply mult. formula: elif not e.is_trivial() and not f.is_trivial(): return e.gauss_sum()*f.gauss_sum()/g.gauss_sum() # If exactly one is trivial, return 0 else: return 0 def latex2(e): return latex(list(e.values_on_gens())) def jacobi_plot(p, e_index, f_index, with_text=True): # Set values G = DirichletGroup(p) e = G[e_index] f = G[f_index] ef = e*f js = jacobi_sum(e,f) e_gs = e.gauss_sum() f_gs = f.gauss_sum() ef_gs = (e*f).gauss_sum() # Compute complex coordinates f_pt = list(f.values_on_gens()[0].complex_embedding()) e_pt = list(e.values_on_gens()[0].complex_embedding()) ef_pt = list(ef.values_on_gens()[0].complex_embedding()) f_gs_pt = list(f_gs.complex_embedding()) e_gs_pt = list(e_gs.complex_embedding()) ef_gs_pt = list(ef_gs.complex_embedding()) try: js = int(js) js_pt = [CC(js)] except: js_pt = list(js.complex_embedding()) # Define plot structure S = circle((0,0),1,rgbcolor='yellow') S += line([e_pt,e_gs_pt], rgbcolor='red', thickness=4) S += line([f_pt,f_gs_pt], rgbcolor='blue', thickness=3) S += line([ef_pt,ef_gs_pt], rgbcolor='purple',thickness=2) S += point(e_pt,pointsize=50, rgbcolor='red') S += point(f_pt,pointsize=50, rgbcolor='blue') S += point(ef_pt,pointsize=50,rgbcolor='purple') S += point(f_gs_pt,pointsize=75, rgbcolor='black') S += point(e_gs_pt,pointsize=75, rgbcolor='black') S += point(ef_gs_pt,pointsize=75, rgbcolor='black') S += point(js_pt,pointsize=100,rgbcolor='green') if with_text: S += text('$J(%s,%s) = %s$'%(latex2(e),latex2(f),latex(js)), (3,2.5),fontsize=15, rgbcolor='black') else: html('$$J(%s,%s) = %s$$'%(latex2(e),latex2(f),latex(js))) return S @interact def exhaustive_jacobi_plot(p=prime_range(3,8)): ga = [jacobi_plot(p,i,j) for i in range(p-1) for j in range(p-1)[i:]] for i in range(len(ga)): ga[i].save('j%d.png'%i,figsize=4,aspect_ratio=1, xmin=-2.5,xmax=5, ymin=-2.5, ymax=2.5) # Since p is odd, will have n = p-1 even. So 1+2+...+n = (n/2)*(n+1). # We divide this by rows of 2. rows = ceil(p*(p-1)/4) s='<table bgcolor=lightgrey cellpadding=2>' for i in range(rows): s+='<tr><td align="center"><img src="cell://j%d.png"></td>'%(2*i) s+='<td align="center"><img src="cell://j%d.png"></td></tr>'%(2*i+1) s+='</table>' html(s)}}} {{attachment:jacobiexh.png}} = Elliptic Curves = == Adding points on an elliptic curve == by David Møller Hansen {{{#!sagecell def point_txt(P,name,rgbcolor): if (P.xy()[1]) < 0: r = text(name,[float(P.xy()[0]),float(P.xy()[1])-1],rgbcolor=rgbcolor) elif P.xy()[1] == 0: r = text(name,[float(P.xy()[0]),float(P.xy()[1])+1],rgbcolor=rgbcolor) else: r = text(name,[float(P.xy()[0]),float(P.xy()[1])+1],rgbcolor=rgbcolor) return r E = EllipticCurve('37a') list_of_points = E.integral_points() html("Graphical addition of two points $P$ and $Q$ on the curve $ E: %s $"%latex(E)) @interact def _(P=selector(list_of_points,label='Point P'),Q=selector(list_of_points,label='Point Q'), marked_points = checkbox(default=True,label = 'Points'), Lines = selector([0..2],nrows=1), Axes=True): curve = E.plot(rgbcolor = (0,0,1),xmin=25,xmax=25,plot_points=300) R = P + Q Rneg = -R l1 = line_from_curve_points(E,P,Q) l2 = line_from_curve_points(E,R,Rneg,style='--') p1 = plot(P,rgbcolor=(1,0,0),pointsize=40) p2 = plot(Q,rgbcolor=(1,0,0),pointsize=40) p3 = plot(R,rgbcolor=(1,0,0),pointsize=40) p4 = plot(Rneg,rgbcolor=(1,0,0),pointsize=40) textp1 = point_txt(P,"$P$",rgbcolor=(0,0,0)) textp2 = point_txt(Q,"$Q$",rgbcolor=(0,0,0)) textp3 = point_txt(R,"$P+Q$",rgbcolor=(0,0,0)) if Lines==0: g=curve elif Lines ==1: g=curve+l1 elif Lines == 2: g=curve+l1+l2 if marked_points: g=g+p1+p2+p3+p4 if P != Q: g=g+textp1+textp2+textp3 else: g=g+textp1+textp3 g.axes_range(xmin=-5,xmax=5,ymin=-13,ymax=13) show(g,axes = Axes) def line_from_curve_points(E,P,Q,style='-',rgb=(1,0,0),length=25): """ P,Q two points on an elliptic curve. Output is a graphic representation of the straight line intersecting with P,Q. """ # The function tangent to P=Q on E if P == Q: if P[2]==0: return line([(1,-length),(1,length)],linestyle=style,rgbcolor=rgb) else: # Compute slope of the curve E in P l=-(3*P[0]^2 + 2*E.a2()*P[0] + E.a4() - E.a1()*P[1])/((-2)*P[1] - E.a1()*P[0] - E.a3()) f(x) = l * (x - P[0]) + P[1] return plot(f(x),-length,length,linestyle=style,rgbcolor=rgb) # Trivial case of P != R where P=O or R=O then we get the vertical line from the other point elif P[2] == 0: return line([(Q[0],-length),(Q[0],length)],linestyle=style,rgbcolor=rgb) elif Q[2] == 0: return line([(P[0],-length),(P[0],length)],linestyle=style,rgbcolor=rgb) # Non trivial case where P != R else: # Case where x_1 = x_2 return vertical line evaluated in Q if P[0] == Q[0]: return line([(P[0],-length),(P[0],length)],linestyle=style,rgbcolor=rgb) #Case where x_1 != x_2 return line trough P,R evaluated in Q" l=(Q[1]-P[1])/(Q[0]-P[0]) f(x) = l * (x - P[0]) + P[1] return plot(f(x),-length,length,linestyle=style,rgbcolor=rgb) }}} {{attachment:PointAddEllipticCurve.png}} == Plotting an elliptic curve over a finite field == {{{#!sagecell E = EllipticCurve('37a') @interact def _(p=slider(prime_range(1000), default=389)): show(E) print "p = %s"%p show(E.change_ring(GF(p)).plot(),xmin=0,ymin=0) }}} {{attachment:ellffplot.png}} = Cryptography = == The Diffie-Hellman Key Exchange Protocol == by Timothy Clemans and William Stein {{{#!sagecell @interact def diffie_hellman(bits=slider(8, 513, 4, 8, 'Number of bits', False), button=selector(["Show new example"],label='',buttons=True)): maxp = 2 ^ bits p = random_prime(maxp) k = GF(p) if bits > 100: g = k(2) else: g = k.multiplicative_generator() a = ZZ.random_element(10, maxp) b = ZZ.random_element(10, maxp) html(""" <style> .gamodp, .gbmodp { color:#000; padding:5px } .gamodp { background:#846FD8 } .gbmodp { background:#FFFC73 } .dhsame { color:#000; font-weight:bold } </style> <h2 style="color:#000;font-family:Arial, Helvetica, sans-serif">%s-Bit Diffie-Hellman Key Exchange</h2> <ol style="color:#000;font-family:Arial, Helvetica, sans-serif"> <li>Alice and Bob agree to use the prime number p = %s and base g = %s.</li> <li>Alice chooses the secret integer a = %s, then sends Bob (<span class="gamodp">g<sup>a</sup> mod p</span>):<br/>%s<sup>%s</sup> mod %s = <span class="gamodp">%s</span>.</li> <li>Bob chooses the secret integer b=%s, then sends Alice (<span class="gbmodp">g<sup>b</sup> mod p</span>):<br/>%s<sup>%s</sup> mod %s = <span class="gbmodp">%s</span>.</li> <li>Alice computes (<span class="gbmodp">g<sup>b</sup> mod p</span>)<sup>a</sup> mod p:<br/>%s<sup>%s</sup> mod %s = <span class="dhsame">%s</span>.</li> <li>Bob computes (<span class="gamodp">g<sup>a</sup> mod p</span>)<sup>b</sup> mod p:<br/>%s<sup>%s</sup> mod %s = <span class="dhsame">%s</span>.</li> </ol> """ % (bits, p, g, a, g, a, p, (g^a), b, g, b, p, (g^b), (g^b), a, p, (g^ b)^a, g^a, b, p, (g^a)^b)) }}} {{attachment:dh.png}} = Other = == Continued Fraction Plotter == by William Stein {{{#!sagecell @interact def _(number=e, ymax=selector([None,5,20,..,400],nrows=2), clr=Color('purple'), prec=[500,1000,..,5000]): c = list(continued_fraction(RealField(prec)(number))); print c show(line([(i,z) for i, z in enumerate(c)],rgbcolor=clr),ymax=ymax,figsize=[10,2]) }}} {{attachment:contfracplot.png}} == Computing Generalized Bernoulli Numbers == by William Stein (Sage-2.10.3) {{{#!sagecell @interact def _(m=selector([1..15],nrows=2), n=(7,(3..10))): G = DirichletGroup(m) s = "<h3>First n=%s Bernoulli numbers attached to characters with modulus m=%s</h3>"%(n,m) s += '<table border=1>' s += '<tr bgcolor="#edcc9c"><td align=center>$\\chi$</td><td>Conductor</td>' + \ ''.join('<td>$B_{%s,\chi}$</td>'%k for k in [1..n]) + '</tr>' for eps in G.list(): v = ''.join(['<td align=center bgcolor="#efe5cd">$%s$</td>'%latex(eps.bernoulli(k)) for k in [1..n]]) s += '<tr><td bgcolor="#edcc9c">%s</td><td bgcolor="#efe5cd" align=center>%s</td>%s</tr>\n'%( eps, eps.conductor(), v) s += '</table>' html(s) }}} {{attachment:bernoulli.png}} == Fundamental Domains of SL_2(ZZ) == by Robert Miller {{{#!sagecell L = [[-0.5, 2.0^(x/100.0) - 1 + sqrt(3.0)/2] for x in xrange(1000, -1, -1)] R = [[0.5, 2.0^(x/100.0) - 1 + sqrt(3.0)/2] for x in xrange(1000)] xes = [x/1000.0 for x in xrange(-500,501,1)] M = [[x,abs(sqrt(x^2-1))] for x in xes] fundamental_domain = L+M+R fundamental_domain = [[x-1,y] for x,y in fundamental_domain] @interact def _(gen = selector(['t+1', 't-1', '-1/t'], nrows=1)): global fundamental_domain if gen == 't+1': fundamental_domain = [[x+1,y] for x,y in fundamental_domain] elif gen == 't-1': fundamental_domain = [[x-1,y] for x,y in fundamental_domain] elif gen == '-1/t': new_dom = [] for x,y in fundamental_domain: sq_mod = x^2 + y^2 new_dom.append([(-1)*x/sq_mod, y/sq_mod]) fundamental_domain = new_dom P = polygon(fundamental_domain) P.ymax(1.2); P.ymin(-0.1) P.show() }}} {{attachment:fund_domain.png}} |
Contents
Integer Factorization
Divisibility Poset
by William Stein
Factor Trees
by William Stein
More complicated demonstration using Mathematica: http://demonstrations.wolfram.com/FactorTrees/
Factoring an Integer
by Timothy Clemans
Sage implementation of the Mathematica demonstration of the same name. http://demonstrations.wolfram.com/FactoringAnInteger/
Prime Numbers
Illustrating the prime number theorem
by William Stein
Prime Spiral - Square FIXME
by David Runde
Prime Spiral - Polar
by David Runde
Modular Forms
Computing modular forms FIXME
by William Stein
Computing the cuspidal subgroup
by William Stein
A Charpoly and Hecke Operator Graph FIXME
by William Stein
Modular Arithmetic
Quadratic Residue Table FIXME
by Emily Kirkman
Cubic Residue Table FIXME
by Emily Kirkman
Cyclotomic Fields
Gauss and Jacobi Sums in Complex Plane
by Emily Kirkman
Exhaustive Jacobi Plotter
by Emily Kirkman
Elliptic Curves
Adding points on an elliptic curve
by David Møller Hansen
Plotting an elliptic curve over a finite field
Cryptography
The Diffie-Hellman Key Exchange Protocol
by Timothy Clemans and William Stein
Other
Continued Fraction Plotter
by William Stein
Computing Generalized Bernoulli Numbers
by William Stein (Sage-2.10.3)
Fundamental Domains of SL_2(ZZ)
by Robert Miller
