N_target = 44 EVEN_N = True # look only for even n-gons #EVEN_N = False # looking at all n # A different target value of N can be given on the command line. # "even" or "all" can also be given on the command line. """ Lemma 1: Every edge is a primitive vector. Lemma 2: If n is even, then P can be assumed to be centrally symmetric. The correctness of the program (knowing when one can stop) depends crucially on two lemmas that relate area to lattice width. The first lemma is about centrally symmetric polygons (arising when n is even), and the second lemma about arbitrary polygons. Definition: The *lattice width* of a lattice polygon P *in a given lattice direction* (a,b) in Z² is one less than the number of lattice lines orthogonal to (a,b) intersecting P. The *lattice width* of P is the smallest number lattice width over all lattice directions. Lemma 3: A centrally symmetric convex polygon P of lattice width w has area at least w²/2. Proof: [ A weaker bound of w²/4 is given in the following paper (p.175, before remark 2): Imre Bárány and Norihide Tokushige, The minimum area of convex lattice n-gons. Combinatorica, 24 (No. 2, 2004), 171-185. https://www.renyi.hu/~barany/cikkek/94.pdf Our proof uses the same setup but refines the argument in the end. The improvement to w²/2 speeds up the computations by a factor 5-6. ] 1. We assume that the lattice width direction is vertical. Thus P touches two horizontal lines H+: y=b and H-: y=-b, with b=w/2. Let (a,0) be the intersection of P with the positive x-axis. 2. If a>=b, then the area is at least 2b²=w²/2, and we are done. Thus, let us assume a<b. 3. By horizontal shearing, we assume that the tangent at (a,0) (or some supporting line) is not vertical and has slope s at least +1. (Since horizontal lattice-preserving shearing changes the "inverse slope" 1/s in increments of 1, the inverse slope can always be brought into the interval (0,1].) 4. The lattice width in direction (1,0) is b'>=b. Thus P touches two vertical lines V+: x=b' and V-: x=-b' 5. The lattice width in direction (1,-1) is b">=b. Thus P touches two diagonal lines D+: y-x=b" and D-: y-x=-b" (Draw a picture!) 6. If (a,0) had a tangent of slope 1, P could not touch D-, since a<b". Thus the slope s is strictly larger than 1, and P must touch V+ in the upper half-plane, and P must touch D- in the lower half-plane. 7. If the intersection of V+ and D- is below the x-axis, we swap the roles of V+ and D- by a horizontal shearing transformation combined with a vertical reflection: Each point (x,y) is mapped to (x-y,-y). This turns V+ into a line of slope 1 and makes D- vertical. Point on the x-axis are fixed. Thus we can assume that the intersection of V+ and D- is above or on the x-axis, at (b',c) with c>=0. 8. Pick a point where P touches each of the six lines, in a symmetric way. Denote them by B,BR,R (Bottom, Bottom-Right, and Right) and the reflected points by T,TL,L. (Draw a picture!) Replace P by the convex hull of these points (making the area smaller). 9. Holding L,B,R, and T fixed, we can move the points BR and TL on their respective lines D- and D+, keeping them symmetric. The area depends linearly on the movement. If the are would decrease when BR is moved towards R and towards V+, we can stop this movement at the point (b",0) where D- intersects the x-axis. The convex hull of (b",0),(-b",0),T,B has area 2b"b>w²/2. 10. Thus we are left with the case that we can decrease the area by moving BR towards H-, reaching the point B'=(d,-b) with d := b"-1>=0. Then we consider the convex hull of B', R, and their two symmetric points. R has coordinates R=(b',e) with e>=c>=0. The determinant of (d,-b) and (b',e) is bb'+de>=bb'>=b², and we are done. QED. Lemma 4: A (not necessarily symmetric) convex polygon P of lattice width w has area at least w²/3. The inequality is strict except when P is a triangle. """ # Since we use "volume", which is TWICE the area, the following constants # are twice as big as the constants in the lemmas. lattice_width_factor_even = 1 lattice_width_factor_odd = 2/3 + 0.00001 # +0.00001 is because of the strict inequality (except for triangles, but # those are found in the first iteration.) def lower_bound(k,w): """lower bound on vol (=2*area) for k-gon of lattice width w """ if k%2: lattice_width_factor = lattice_width_factor_odd else: lattice_width_factor = lattice_width_factor_even LB = ceil( w**2 * lattice_width_factor ) if k%2 != LB%2: LB +=1 # vol must be odd for odd k and even for even k return LB """ Proof of Lemma 4: The polygon P and the symmetrized version Q = 1/2 * (P-P) (which is half the "difference body" P-P={ a-b | a,b in P }), have the same width in every direction, and Q is centrally symmetric. We can apply Lemma 3 to Q. The area of Q is at most 3/2 times the area of P, and equality holds only when P is a triangle. See any of the following papers. QED Hans Rademacher, Über den Vektorenbereich eines konvexen ebenen Bereiches. Jahresbericht d. Deutschen Math.-Verein. 34 (1926): 64-78. http://eudml.org/doc/145705 Theodor Estermann, Zwei neue Beweise eines Satzes von Blaschke und Rademacher. Jahresbericht d. Deutschen Math.-Verein. 36 (1927): 197-200. http://eudml.org/doc/145763 Theodor Estermann, Über den Vektorenbereich eines konvexen Körpers. Math Z. 28 (1928), 471-475. https://doi.org/10.1007/BF01181177 C. A. Rogers and G. C. Shephard, The difference body of a convex body. Arch. Math. 8 (1957), 220-233. doi:10.1007/BF01899997 I believe the true factor for Lemma 4 might be 3/8, and it is sharp for the triangle (-w/2,-w/2),(w/2,0),(0,w/2), for the case of even w. (By analogy with the ordinary width, not lattice width: "It is well-known that of all convex sets of a given width, the equilaterial triangle has the smallest area.") In the proof, we would consider the 8-gon formed as the intersection of a horizontal slab of width w, a vertical slab of width w1>=1, and two slabs of directions 45° and -45° of *vertical* widths w3>=w and w4>=w, respectively. (Their Euclidean widths are >=w/sqrt(2).) Some sides of the 8-gon may degenerate into a point, in case several lines are concurrent. The body must touch each of these edges (even if they are degenerated to a point). Lemma 5. A) For three successive vertices v,v',v" of P, the triangle vv'v" contains no interior point. B) If n is even, then there cannot even be points on the diagonal vv". Proof: A) Otherwise v' could be replaced by that point, reducing the area. B) Bárány and Tokushige proved that the set E of edge vectors is convex in the following sense: Every vector e in the convex hull of E that is a primitive vector, must belong to E. Let e=v->v' and e'=v'->v". If vv"=e+e' would have gcd d>2 the vector (e+e')/d, which is in the convex hull of e,e',-e,-e', is missing from E, a contradiction. This Lemma is used to shortcut the propagation process. Property B is really strong. It does not hold for the optimal 5-gons, 7-gons, and 11-gons, but by preliminary experiments, it seems that for larger k, it might hold also in the odd case. Property B means the following. Consider a pair (q,f). Then the outgoing edges q->p that need to be considered end on the next lattice line parallel to f. The points q form an arithmetic progression p0, p0+f, p0+2f, .... When considering a range of diameter D (roughly the maximum height), then these are about D/|f| successor points altogether (on all levels py). Under Property A, we have to consider in addition all multiples of the vectors (q-f)->(p0+t*f). I.e., p = (q-f) + s*(p0+t*f-(q-f)) = s*p0 - (s-1)*q + (s*t-s-1)*f This would be D/|f| log(D*|f|) ~ D log D/|f| points. Typically, the average length of the lists (the number of vectors f for a point q) is very small (about 4 for height py=100), and |f| is not always short. (By contrast, currently, every point p gets O(D²) inputs, one from every other point q (looking at the propagation from the incoming side).) Lemma 6: If n is even, then P has two edges parallel to the lattice width direction Proof: Suppose not. Then there is a pair of extreme VERTICES where the strip of minimum lattice width touches P. (as opposed to a pair of extreme EDGES.) Let v be such a vertex, and let v- and v+ its neighbors. As usual, we assume that the lattice width direction is horizontal, as we assume that v=(0,0). Let u- and u+ be the intersection of the edges v -> v- and v -> v+ with the line y=1. The open interval (u-,u+) cannot contain an integer point, by Lemma 5B. Hence, by horizontal shearing we can assume 0 <= u- < u+ <= 1. P is contained in the wedge between v -> v- and v -> v+. Hence the coordinates highest point v'=(x,y) (opposite to v) fulfill the inequality 0 <= x <= y. If x<y, then the width in vertical direction is smaller than the horizontal width (which is y). If x=0, then the width in direction (1,1) is smaller than y. Both cases lead to a contradiction. QED =================================================== Also in the algorithm (as in the proof of Lemma 3), we make the assumption that the lattice width direction is vertical. (The smallest number of parallel lattice lines intersecting P are the horizontal lines.) We assume that the lowest point (or one of the two lowest points) is the origin O=(0,0). The approach is by dynamic programming. It is the algorithm in the following paper: Finding minimum area k-gons. David Eppstein, Mark Overmars, Günter Rote und Gerhard Woeginger. Discrete and Computational Geometry 7 (1992), 45-58. https://doi.org/10.1007/BF02187823 see also Counting convex polygons in planar point sets. Joseph S. B. Mitchell, Günter Rote, Gopalakrishnan Sundaram, and Gerhard Woeginger. Information Processing Letters 56 (1995), 45-49 http://page.mi.fu-berlin.de/rote/Papers/abstract/Counting+k-subsets+and+convex+k-gons+in+the+plane.html ========================================================= For each pair of points p and q and every k, we remember the smallest counterclockwise k-gon (0,0), .... , q, p Restriction: * The successive points move strictly upwards. (Thus, p is the top point, and there is a "long edge" from O to p.) Polygons with one or two horizontal edges will be treated specially. By horizontal lattice-preserving shearing transformations, it suffices to store results only for the range 0 <= px < py (The answer for any point p in the plane can then be figured out by an appropriate shearing transformation.) In the end we put together a counterclockwise k-gon ending in p with a counterclockwise k'-gon starting in O or in (-1,0) and ending in p or in p-(1,0). Main data structure: ==================== min_gon is a dictionary min_gon[(px,py)] is a dictionary min_gon[(px,py)][k] is a list of pairs (f,vol)=((fx,fy),vol), where the incoming edge f=(fx,fy)=q-p is a primitive vector (pointing upwards, fy>0) and vol = 2 * min-area of a k-gon ((0,0), .... ,q,p). In the list, the vectors f are sorted clockwise by direction, and the corresponding values "vol" are strictly increasing. (If pairs don't fit this order, we can eliminate one of them, as DOMINATED.) We fill this dictionary row by row, increasing py. [ If we want to enumerate ALL optimal solutions, we should allow the "vol" values to be weakly increasing. We should then use the weaker lower bound also for even k, in case we are interested in solutions that are not centrally symmetric. Another possibility is COUNTING the number of solutions. (For each height separately. Every solution has up to 2k representations with at least one flat edge, and infinitely many without flat edge.) We should consider all solutions of height up until the value is confirmed. ] [ We should have proceeded from left to right, swapping x and y. Then the notion of slope could have been used more naturally. ] "vol" is TWICE the area. https://oeis.org/A070911 1,2,5,6,13,14,21,28,43,48,65,80,103,118, / 151,174,213,242,289,328,387,420,497 The odd values after the slash were unconfirmed according to OEIS, as of August 2023. POTENTIAL IMPROVEMENTS: In the n-gons computed, the left and right halves have always equal size +-1, for odd n? So k needs to go only to n/2, in practice. (This is definitely warranted for even-n-gons). Can one show that the difference is at most 1? Is the "best"-list for a given point a convex function of k? [ NO: This is very much not the case, not even in the middle of the range: Example: (px,py) = (20,29), k in range(2, 31) Best vols:[0,5,6,13,14,29,64,81,122,171,210,255,334,421,502,641,746,865,1072,1235,1438,1657,1924,2175,2440,2821,3214,3631,4096] Differences are far from monotone: 5 1 7 1 15 35 17 41 49 39 45 79 87 81 139 105 119 207 163 203 219 267 251 265 381 393 417 465 ] The program should have an option to evaluate only the even ones; they are faster anyway because of the better constant factor in the lower bound, AND for the above reason. For computing best paths to the top vertex (without constraints on the incoming slopes, when we take best1[0] in the final combination), one could do a meet-in-the-middle divide-and-conquer approach. Partition into k/2 + k/2, and try all meeting points. -> reduce range of k by 1/2; factor-2 improvement Another big opportunity for improvement would be to restrict the range of outgoing points p that can be meaningfully reached from a given q. """ known_values = {i+3:v for i,v in enumerate(( 1,2,5,6,13,14,21,28,43,48,65,80,103,118,151,174,213,242,289, 328,387,420,497,548,625,690,783,860,967,1046,1177,1264,1409, 1498,1671,1780,1955,2078,2279,2444,2651,2824,3051,3240,3493, 3676,3969,4176,4489,4714,5045,5302,5623,5906,6243,6556,6991, 7224,7731,8040,8479,8878,9365, 9804,10313,10774,11323,11796, 12359,12836,13459,13948,14615,15114,15811,16364,17073,17670, 18431,19024,19833,20436,21321,21968,22901,23518,24527,25270, 26283,27050,28105,28896,29993,30798,31959,32830,34019,34946, 36145,37140 # up to n=102 ))} # optional, just for guidance; not used by the algorithm known_values = {} from math import tan,pi,gcd,sqrt,ceil import itertools import sys if len(sys.argv)>1: N_target = int(sys.argv[1]) for s in sys.argv[1:]: if s.lower() == "even": EVEN_N = True if s.lower() == "all": EVEN_N = False if EVEN_N: b_file = open("b089187.txt","w") a_file = open("a089187.txt","w") output_so_far = 1 else: b_file = open("b070911.txt","w") a_file = open("a070911.txt","w") output_so_far = 2 # the values are not confirmed in order # because the even ones come out earlier. min_gon = dict() record_vol = [n**3 for n in range(N_target+1)] # loose upper bound as a start value min_height = [0 for n in range(N_target+1)] # smallest height of a min-area polygon max_height = [0 for n in range(N_target+1)] # largest height of a min-area polygon with at least one horizontal edge, # AS FOUND BY THE PROGRAM, not definite. The program does not look for the largest heights. k_START = 2 if EVEN_N else 3 # the values below k_start are trivial confirmed = [True]*k_START + [False]*(N_target+1-k_START) # checks which entries are confirmed. how_achieved = [None]*(N_target+1) # Auxiliary procedures to construct the solution, once the optimal # area has been determined def find_polygon(px,py,k,alpha0=0): # find the smallest k-gon ending in (px,py) by backtracing # 0<px<py. The startpoint (0,0) is not produced. alphaT = alpha0 result = [(px+alphaT*py,py)] (fx,fy),vol = min_gon[px,py][k][0] for k in range (k-1,1,-1): qx,qy = px-fx,py-fy vol -= qx*py-qy*px alpha = qx//qy px,py = qx-alpha*qy,qy alphaT += alpha # alpha is accumulated result.append((px+alphaT*py,py)) for (fx,fy),vol2 in min_gon[px,py][k]: if vol2==vol: # this must be the right entry, recognition by vol is easiest break else: error("not found",k,px,py) return result def print_solution(k,how_achieved,vol,for_b_file=False): # parameters k and vol are only for output, have no effect on the procedure. how,px,py,k1,k2 = how_achieved #print(f"{(k,how,px,py,k1,k2)=}") p1 = find_polygon(px,py,k1) p1.reverse() if how=="1 FLAT": top_point = px-1 else: top_point = px ppx = (-top_point)%py alpha0 = (-top_point - ppx) // py # usually, alpha0==-1 leftshift = 1 if how=="2 FLAT" else 0 p2 = [(-x-leftshift,y) for x,y in find_polygon(ppx,py,k2,alpha0)] #print(f"-- {how=} {alpha0=}, {p2[0][0]=},{top_point=},{ppx=},{px=},{py=},{k=}") assert p2[0][0]==top_point-leftshift if how=="DIAGONAL": p2 = p2[1:] #eliminate the common point last = [(-1,0)] if how=="2 FLAT" else [] polygon = [(0,0)]+p1+p2+last # counterclockwise if not for_b_file: print (f"--- a smallest {2*k if EVEN_N else k}-gon,", f"area = {vol}:" if EVEN_N else f"area = {vol}/2:", polygon) else: result = vol # //2 if EVEN_N else vol print(k,result, file=b_file, flush=True) print(k,result, ",".join(f"({x},{y})" for x,y in polygon), file=a_file, flush=True) # compact format without spaces ## Start the computation py_target_start = 1 py_target_end=20 # will be iteratively doubled when not sufficient target_height = 1 while py_target_start<= target_height: collect=tuple([[[] for k in range(N_target+1)] for px in range(py)] for py in range(py_target_end+1) ) print(f"ONE ROUND STARTS: propagate from rows {1}-{py_target_end-1}", f"to rows {py_target_start}-{py_target_end}") for qy in range(py_target_end): # expand from row qy into rows in py_target_start <= py <= py_target_end py_lower = max(py_target_start,qy+1) print(f"propagate from row {qy}", f"to rows {py_lower}-{py_target_end}", flush=True) for qx in range(qy): for k,sides in min_gon[qx,qy].items(): if k>=N_target: continue for (fx,fy),vol in sides: # search first point q-d above qy-fy on the next lattice line # This could be done by the extended Euclidean algorithm. # Linear search is o.k. in this context. for dy in range(fy): dx = fx*dy//fy + 1 # assert dx*fy-fx*dy>0 if dx*fy-fx*dy==1: break #else: assert 0 #print(f"{(dx,dy)=} {(fx,fy)=}") if EVEN_N: px,py = qx+fx-dx,qy+fy-dy while py<=py_target_end: if py>=py_target_start: Delta_x,Delta_y = px-qx,py-qy # gcd(Delta_x,Delta_y)==1 by construction alpha = px//py collect[py][px-alpha*py][k].append( (vol + qx*py-qy*px, Delta_x-alpha*Delta_y,Delta_y)) px,py = px+fx,py+fy else: rx,ry = qx-fx,qy-fy k0x,k0y = qx-dx,qy-dy # shoot from r in direction k0 while k0y<=py_target_end: stepx,stepy = k0x-rx,k0y-ry px,py = k0x,k0y # the inner loop is necessary for odd n, as far as we know. while py<=py_target_end: if py>=py_lower: Delta_x,Delta_y = px-qx,py-qy if gcd(Delta_x,Delta_y)==1: alpha = px//py collect[py][px-alpha*py][k].append( (vol + qx*py-qy*px, Delta_x-alpha*Delta_y,Delta_y)) px,py = px+stepx,py+stepy k0x,k0y = k0x+fx,k0y+fy py = qy+1 # consolidate the lists in collect[py] if py<py_target_start: continue for px in range(py): # 2-gons are the base case. if gcd(px,py)==1: min_gon[px,py]={2:[((px,py),0)]} else: min_gon[px,py]={} work = num_work = 0 # statistics for px,l2 in enumerate(collect[py]): for k,triples in enumerate(l2): if not triples: continue work += len(triples) num_work += 1 # statistics result = [] prev_fx, prev_fy = -1,0 # ensure that test is true at the first iteration prev_vol = -1 for vol,fx,fy in (sorted(triples)): # smallest area first <==> result should be ordered clockwise. if (# this is the first time OR (fx,fy) is clockwise from prev: # fx/fy > prev_fx/prev_fy fx*prev_fy > prev_fx*fy): if vol!=prev_vol: result.append(((fx,fy),vol)) else: # replace the last entry: result[-1] = ((fx,fy),vol) # adjust this treatment when generating ALL optimal solutions: # Then multiple vol-entries should be kept, but SORTED BY SLOPE. prev_fx,prev_fy,prev_vol = fx,fy,vol min_gon[px,py][k+1]=result # print (f"{(px,py)=} {k=} {triples=}"); print(result) # put together solutions # update target_height # compute new record areas improved = [False]*(N_target+1) def check_record(k, vol, how, text=""): if k<k_START or k>N_target: return if vol<record_vol[k]: #if text: print(text) if confirmed[k]: print(k,vol,how,text) assert not confirmed[k] record_vol[k] = vol how_achieved[k] = how min_height[k] = py improved[k] = True if vol==record_vol[k] and how[0] != "DIAGONAL": # DIAGONAL case not counted: Without a flat edge, height is unbounded. max_height[k] = py # Combine "right" k1-gons with "left" k2-gons for px in range(py): if EVEN_N: # even case is easy: insert a horizontal edge at the top: ppx0 = (-px)%py # unnecessarily complicated: # One could just turn the solution by 180 degrees. # It is possible that this produces an asymmetric solution. for k1,best1 in min_gon[px,py].items(): _,vol1 = best1[0] check_record(k1, # a (2k)-gon vol1+py, # actual area, not "volume" ("2 FLAT",px,py,k1,k1)) else: ppx0 = (-px)%py # A) exact match along long edge (0,0)-(px,py) ppx1 = (1-px)%py # B) horizontal edge of length 1 at the top for k1,best1 in min_gon[px,py].items(): _,vol1 = best1[0] for k2,best2 in min_gon[ppx0,py].items(): _,vol2 = best2[0] check_record(k1+k2-2, vol1+vol2, ("DIAGONAL",px,py,k1,k2)) check_record(k1+k2, vol1+vol2+2*py, ("2 FLAT",px,py,k1,k2)) for k2,best2 in min_gon[ppx1,py].items(): _,vol2 = best2[0] check_record(k1+k2-1, vol1+vol2+py, ("1 FLAT",px,py,k1,k2)) target_height = 0 # how high we still have to go for k in range(k_START,N_target+1): if not confirmed[k]: num_vertices = 2*k if EVEN_N else k vol_factor = 2 if EVEN_N else 1 if lower_bound(num_vertices, py+1)>=record_vol[k]*vol_factor: # evaluated for increased py in the next iteration #print("LB",k,py+1,LB) confirmed[k] = True print_solution(k,how_achieved[k],record_vol[k]) else: # find necessary target height to confirm the current record: t_k = py+2 while lower_bound(num_vertices, t_k)<record_vol[k]*vol_factor: t_k += 1 target_height = max(target_height,t_k-1) if confirmed[k] and k==output_so_far+1: print_solution(k,how_achieved[k],record_vol[k], True) # print to a_file and b_file output_so_far += 1 print("N=%d, height %d"%(N_target,py) + (".\nRESULT:" if target_height==0 else "->%d:"%target_height), ",".join(("*" if improved[k] else "")+ str(record_vol[k])+ ("" if confirmed[k] else "?") for k in range(k_START,N_target+1)), # # Meaning of the signs: # ? still unconfirmed # * has just been improved # "work=%d"%work, # work done in this step, total items in "collect" "" if num_work==0 else "avg=%.2f"%(work/num_work), # average list size in "collect" flush=True) if py%10==0: print("smallest heights:",*("%2d"%x for x in min_height[k_START:])) print("largest heights: ",*("%2d"%x for x in max_height[k_START:])) if py%10==0: # average list length rises very slowly: # 3.5 for height py=70, 3.9 for height py=100, # max length, on the other hand is 71 for py=100. # the long lists are mostly associated with steep vectors p, # and they contain many f-vectors of the form (-1,..) or (1,..) print("avg list length", sum(len(best1) for ppy in range (1,py+1) for ppx in range (ppy) for k1,best1 in min_gon[ppx,ppy].items()) / sum(1 for ppy in range (1,py+1) for ppx in range (ppy) for k1,best1 in min_gon[ppx,ppy].items())) print("max list length", max(len(best1) for ppy in range (1,py+1) for ppx in range (ppy) for k1,best1 in min_gon[ppx,ppy].items())) if all(confirmed): break py_target_start = py_target_end+1 py_target_end = py_target_end * 2 # make smaller steps if memory space is an issue. py_target_end = min(py_target_end,target_height) a_file.close() b_file.close() print("k: ",*("%2d"%x for x in range(k_START,N_target+1))) print("smallest heights:",*("%2d"%x for x in min_height[k_START:])) print("largest heights: ",*("%2d"%x for x in max_height[k_START:])) # "largest heights" are preliminary, since the algorithm stops # too early to find the largest heights. # "smallest heights" is definitely the smallest lattice width that # a convex lattice n-gon of minimum area can have. # There is quite a large variation of heights, in line with the # result of Bárány and Tukushige that, in the limit, # the shapes resemble ellipses that are more and more oblong. # But already a single n-gon has necessarily a large variation of # lattice widths, when considering all edge directions. # This is because the lattice width in direction (a,b) with gcd(a,b)=1 # is the Euclidean width times the length of (a,b). # Since the lengths of edge vectors vary a lot (a vector (1,0) # may be adjacent to a vector (L,1) for large L, while the # Euclidean width will hardly change when turning to the next edge, # a large fluctuation is inherent. # The limiting behavior predicted by Bárány and Tukushige # may set in much later than how far the calculations go.