MATH9520_4.23_CodingHW

400 days ago by aburto2

1. Compute the list of champions and runts for the uncorrected number of Goldbach representations.

N=10000 prime_list = [] for n in srange(N): if is_prime(n): prime_list.append(n)
len(prime_list)
 1229 1229
Goldbach_reps=[0 for n in srange(N)] for p1 in prime_list: for p2 in prime_list: if p1+p2<N: Goldbach_reps[p1+p2]+=1
list_plot(Goldbach_reps)
gr_evens=[] for i in srange(1,N/2): gr_evens.append([2*i,Goldbach_reps[2*i]])
list_plot(gr_evens[:500])
current_max = 0 champions = [] for i in gr_evens: if i[1]>current_max: champions.append(i) current_max=i[1]
list_plot(champions)
current_min=10000 runts=[] for i in reversed(gr_evens): if i[1]<current_min: runts.append(i) current_min=i[1]
list_plot(runts)

2. Compute the vertices of the upper convex hull of the champions, and the lower convex hull of the squared-transformed data runts. How many vertices are there relative to the number of champions/runts respectively?

We will do this using Graham's scan.

First, dealing with the vertices of the upper convex hull of the champions: we should note that all champions are listed in increasing $y$ value, if this was not true then a point in question would not be a champion. Thus the first champion has the lowest $y$-corrdiante and is the leftmost point. The code used to find the convex hull was orginally created by Tom Switzer and can be found at https://gist.github.com/arthur-e/5cf52962341310f438e96c1f3c3398b8.

def convex_hull_graham(points): TURN_LEFT, TURN_RIGHT, TURN_NONE = (1, -1, 0) def cmp(a, b): return (a > b) - (a < b) def turn(p, q, r): return cmp((q[0] - p[0])*(r[1] - p[1]) - (r[0] - p[0])*(q[1] - p[1]), 0) def _keep_left(hull, r): while len(hull) > 1 and turn(hull[-2], hull[-1], r) != TURN_LEFT: hull.pop() if not len(hull) or hull[-1] != r: hull.append(r) return hull points = sorted(points) l = reduce(_keep_left, points, []) u = reduce(_keep_left, reversed(points), []) return l.extend(u[i] for i in range(1, len(u) - 1)) or l
champ_vertices = [] champ_vertices = convex_hull_graham(champions) A = list_plot(champions,color="blue") B = list_plot(champ_vertices,color="red") show(A+B)
len(champ_vertices)
 11 11

By looking at this graph, we can see that in order to close the convex hull, the algorithm picked up 1 extra point. Because the algorithm works in a counter clockwise motion, this point appears at the second point. Thus, we say the number of vertices on the upper convex hull is 10.

champ_vertices
 [[4, 1], [9030, 606], [9240, 658], [6930, 536], [2310, 228], [840, 102], [630, 82], [210, 38], [120, 24], [24, 6], [10, 3]] [[4, 1], [9030, 606], [9240, 658], [6930, 536], [2310, 228], [840, 102], [630, 82], [210, 38], [120, 24], [24, 6], [10, 3]]
champ_vertices.pop(1)
 [9030, 606] [9030, 606]
10/52
 5/26 5/26

We need to transform the set of runts in order to be able to find the upper convex hull.

def f(x): return(x[0],x[1]^2)
SquareRunts=[] SquareRunts = map(f,runts) runt_vertices = [] runt_vertices = convex_hull_graham(SquareRunts) A = list_plot(SquareRunts,color="blue") B = list_plot(runt_vertices,color="red") show(A+B)

We should note that since here we are looking for the lower convex hull, we do not have the extra point problem caused by the algorithm searning in the counter clockwise direction.

len(runt_vertices)
 13 13
runt_vertices
 [(2, 0), (68, 16), (128, 36), (398, 169), (992, 676), (1718, 1681), (2672, 3136), (2936, 3844), (5948, 12100), (9602, 23409), (9866, 29241), (9974, 34225), (9998, 38809)] [(2, 0), (68, 16), (128, 36), (398, 169), (992, 676), (1718, 1681), (2672, 3136), (2936, 3844), (5948, 12100), (9602, 23409), (9866, 29241), (9974, 34225), (9998, 38809)]
13/65
 1/5 1/5

3. Consider the factorizations of the n-values of the vertices of the convex hulls. Do the champions get progressively more small primes? Do the vertices of the runts show any pattern you can discern?

champ_vertices
 [[4, 1], [9240, 658], [6930, 536], [2310, 228], [840, 102], [630, 82], [210, 38], [120, 24], [24, 6], [10, 3]] [[4, 1], [9240, 658], [6930, 536], [2310, 228], [840, 102], [630, 82], [210, 38], [120, 24], [24, 6], [10, 3]]
runt_vertices
 [(2, 0), (68, 16), (128, 36), (398, 169), (992, 676), (1718, 1681), (2672, 3136), (2936, 3844), (5948, 12100), (9602, 23409), (9866, 29241), (9974, 34225), (9998, 38809)] [(2, 0), (68, 16), (128, 36), (398, 169), (992, 676), (1718, 1681), (2672, 3136), (2936, 3844), (5948, 12100), (9602, 23409), (9866, 29241), (9974, 34225), (9998, 38809)]
Champ_N=[] for i in champ_vertices: Champ_N.append(i[0]) Runt_N=[] for i in runt_vertices: Runt_N.append(i[0])
Champ_N.sort() Champ_N
 [4, 10, 24, 120, 210, 630, 840, 2310, 6930, 9240] [4, 10, 24, 120, 210, 630, 840, 2310, 6930, 9240]
Runt_N.sort() Runt_N
 [2, 68, 128, 398, 992, 1718, 2672, 2936, 5948, 9602, 9866, 9974, 9998] [2, 68, 128, 398, 992, 1718, 2672, 2936, 5948, 9602, 9866, 9974, 9998]

4. Do all the same things for the corrected data.

def HL(n): # want this to work for odd numbers too, but how to know whether to skip the first factor if it is 2? g=list(factor(2*n)) # this forces 2 to be a prime factor, which makes it easy to skip cf=1 #correction factor for i in g[1:]: # this skips the prime 2, and goes through all the other prime,power pairs p=i[0] # here p is the current odd prime cf*=(p-1)/(p-2) return(cf)
HL_values=[] for i in srange(1,20000): HL_values.append([2*i,HL(2*i)])
list_plot(HL_values)
corrected_gr=[] for i in gr_evens: corrected_gr.append([i[0],i[1]/HL(i[0])])
list_plot(corrected_gr)
current_max = 0 corrected_champ = [] for i in corrected_gr: if i[1]>current_max: corrected_champ.append(i) current_max=i[1]
A = list_plot(corrected_champ,color="red") B = list_plot(champions,color="blue") show(A+B)
cor_champ_vertices = [] cor_champ_vertices = convex_hull_graham(corrected_champ) A = list_plot(corrected_champ,color="blue") B = list_plot(cor_champ_vertices,color="red") show(A+B)
len(cor_champ_vertices)
 13 13

Again we see that the algorithm as it stands picks up an extra point to enclose the set. Thus we conclude that the number of corrected champions on the upper convex hull is 12.

cor_champ_vertices
 [[4, 1], [9052, 35003/180], [9946, 1118475/4972], [6304, 7995/49], [2974, 144045/1486], [1618, 50841/808], [1108, 6875/138], [694, 12765/346], [274, 2835/136], [142, 207/14], [64, 10], [34, 105/16], [16, 4]] [[4, 1], [9052, 35003/180], [9946, 1118475/4972], [6304, 7995/49], [2974, 144045/1486], [1618, 50841/808], [1108, 6875/138], [694, 12765/346], [274, 2835/136], [142, 207/14], [64, 10], [34, 105/16], [16, 4]]
cor_champ_vertices.pop(1)
 [9052, 35003/180] [9052, 35003/180]
len(cor_champ_vertices)
 12 12
len(corrected_champ)
 103 103
12/103
 12/103 12/103
current_min=10000 corrected_runts=[] for i in reversed(corrected_gr): if i[1]<current_min: corrected_runts.append(i) current_min=i[1]
C = list_plot(corrected_runts,color="green") D = list_plot(runts,color="black") show(D+C)
SquareCorRunts=[] SquareCorRunts = map(f,corrected_runts) cor_runt_vertices = [] cor_runt_vertices = convex_hull_graham(SquareCorRunts) A = list_plot(SquareCorRunts,color="blue") B = list_plot(cor_runt_vertices,color="red") show(A+B)
len(cor_runt_vertices)
 17 17
len(corrected_runts)
 92 92
17/95
 17/95 17/95
cor_champ_vertices
 [[4, 1], [9946, 1118475/4972], [6304, 7995/49], [2974, 144045/1486], [1618, 50841/808], [1108, 6875/138], [694, 12765/346], [274, 2835/136], [142, 207/14], [64, 10], [34, 105/16], [16, 4]] [[4, 1], [9946, 1118475/4972], [6304, 7995/49], [2974, 144045/1486], [1618, 50841/808], [1108, 6875/138], [694, 12765/346], [274, 2835/136], [142, 207/14], [64, 10], [34, 105/16], [16, 4]]
cor_runt_vertices
 [(2, 0), (6, 1/4), (12, 1), (30, 81/16), (68, 225/16), (98, 25), (128, 36), (398, 6558721/39204), (992, 142129/225), (1718, 1234608769/736164), (2672, 21344400/6889), (4472, 9966649/1296), (5948, 6670805625/552049), (9602, 59902073001/2560000), (9866, 8777628721/300304), (9994, 583560649/17161), (9998, 969061079281/24980004)] [(2, 0), (6, 1/4), (12, 1), (30, 81/16), (68, 225/16), (98, 25), (128, 36), (398, 6558721/39204), (992, 142129/225), (1718, 1234608769/736164), (2672, 21344400/6889), (4472, 9966649/1296), (5948, 6670805625/552049), (9602, 59902073001/2560000), (9866, 8777628721/300304), (9994, 583560649/17161), (9998, 969061079281/24980004)]
Cor_Champ_N=[] for i in cor_champ_vertices: Cor_Champ_N.append(i[0]) Cor_Runt_N=[] for i in cor_runt_vertices: Cor_Runt_N.append(i[0])
Cor_Champ_N.sort() Cor_Champ_N
 [4, 16, 34, 64, 142, 274, 694, 1108, 1618, 2974, 6304, 9946] [4, 16, 34, 64, 142, 274, 694, 1108, 1618, 2974, 6304, 9946]
Cor_Runt_N.sort() Cor_Runt_N
 [2, 6, 12, 30, 68, 98, 128, 398, 992, 1718, 2672, 4472, 5948, 9602, 9866, 9994, 9998] [2, 6, 12, 30, 68, 98, 128, 398, 992, 1718, 2672, 4472, 5948, 9602, 9866, 9994, 9998]

5. Can we use the vertices of the convex hulls for the corrected data to produce envelopes which enclose the entire set of corrected data, but meet it at interesting points. Can we use this to measure the spread of data? Compare to random models?

def g(x): return([x[0],sqrt(x[1])])
Flatrunt_vertices=[] Flatrunt_vertices=map(g,runt_vertices)
A=list_plot(champ_vertices,color="red") B=list_plot(Flatrunt_vertices,color="green") C=list_plot(gr_evens,color="blue") show(C+A+B)
Flatcor_runt_vertices=[] Flatcor_runt_vertices=map(g,cor_runt_vertices)
A=list_plot(cor_champ_vertices,color="red") B=list_plot(Flatcor_runt_vertices,color="green") C=list_plot(corrected_gr,color="blue") show(C+A+B)