CC BY
64Archaeoastronomy and Ancient Technologies 2016, 4(1), 1-18; http://aaatec.org/art/a_mg1
Archaeoastronomy and Ancient Technologies
www.aaatec.org ISSN 2310-2144
Theoretical Analysis of the Lunar Epicyclic Pin and Aperture System in the Antikythera Mechanism
Maria Giannopoulou
Independent researcher in collaboration with the University of Athens National and Kapodistrian University of Athens, Zografou, Attica, Greece;
E-Mail: [email protected]
Abstract
In this paper I develop a theoretical analysis of a two gear eccentric system in the Antikythera Mechanism, the oldest known analog computer. One gear carries a pin and as it rotates it transmits the motion to the other gear, which has an aperture. The angular velocity of the second gear depends from the shape of the aperture. Until now it was believed that the aperture was an orthogonal parallelogram. Prof. Xenophon Moussas expressed the idea that its shape could also be an ellipse or ellipse like. In this work the shapes Slot, Ellipse and Parabola are examined. The resulting motion of the gears is determined by a model in the framework of Lagrangian Mechanics. The theoretical results of the three shapes are compared with ephemeris data in order to specify which of them fits better with the data.
Keywords: Antikythera mechanism, ancient astronomy, ancient physics, prehistoric science, ancient technologies.
Introduction
The Antikythera Mechanism was discovered by Greek sponge divers in 1900 at the sea region of Antikythera island in Greece. Named after its place of discovery in a Roman shipwreck, it is technically more complex than any known device for at least a millennium afterwards [1]. It is the oldest computer and was produced probably in the second half of the 2n century B.C. (between 150 and 100). It could make mathematical calculations by using gears made of bronze and estimate the positions of Sun and Moon relative with the stars and represents them to cyclic scales which show the zodiac, that is to say the sky with the stars (Fig.1) [2]. It also gives the date and probably the hour too. It displays the lunar phases and also keeps various calendars like the tropical which is still used for the agrarian works and lunar-sun calendars that the Greeks were using for their festivals, like the Olympic or the Nemea Games and nowadays for Easter. It can also predict the Sun and Lunar eclipses. Its dimensions were about 32x22x7 cm [2].
It is notable that the Antikythera Mechanism gives lunar motion by using a mechanism which emulate Kepler's 2nd law. That is to say the Moon moves faster when it is at the perigee and slower at the apogee. The exact definition of Luna's location in the sky is important because it can be used to calculate the geographical length [2].
Archaeoastronomy and Ancient Technologies 2016, 4(1), 1-18 In Figure 2a is represented a schematic sectional diagram of the gearing.
Figure 1. Reconstruction of the Antikythera mechanism in the National Archaeological Museum, Athens1.
a
1 https://archeocomputing.files.wordpress.com/2011/03/antikythera-mechanism2.jpg
b
Figure 2. Antikythera mechanism: a - schematic sectional diagram of the gearing [3, Fig. 6], [4]; b - detailed view of Pin and Slot Hipparchos' lunar mechanism (gears k1, k2, e5, e6).
Figure 2b shows the gears that emulate Luna's motion. They are e5, e6, k1 and k2. As e5 rotates, it forces k1 which transmits the motion to k2 via a pin. Then k2 forces e6 to rotate. The gear k2 has an aperture and it's center is lightly offset from k1's center. The last two gears reproduce the differential lunar motion as seen from Earth. The result is displayed via e1 and b3 gears. Figure 3 shows the diagram of the lunar anomaly mechanism.
Figure 3. Diagram of the Lunar anomaly mechanism [3, Fig. 8], [5].
Prof. X. Moussas expressed the idea that the shape of the aperture might not be a parallelogram (slot- as drawn in Figure 3) but could also be an ellipse or ellipse like -maybe parabola, as shown in Figure 4.
a b
Figure 4. Lunar anomaly mechanism: a - CT slide, b - detailed view of ellipse by X. Moussas.
[both pictures offered by X. Moussas ]
In the following theoretical analysis of the two-gear system k1 and k2, three different aperture shapes Slot, Ellipse and Parabola are examined in order to find out which of the three fits better to the ephemeris data.
Slot - Geometrical analysis
Let us consider two circles with the same radius ro (Fig. 5) . Let d be the distance between their centers O and O' which is constant. So the length of the slot is 2d. The black circle represents the gear that carries the pin K, which moves in a circular path around O. Because the pin has dimensions we will from now on consider that point K is the center of the pin. (So ro is the
distance between the center of gear kl and the center of the pin K.) Its angular velocity coQ is
constant. The vector OK = fo has also constant length rQ. The blue circle represents the gear that carries the slot. The center of the slot C, circuits around O'. It is obvious that the distance O'C is
also rQ. We have 00' = d so the radius r =0'K is f = r -d.
Figure 5. A schematic illustration of the gears with the relative offset of their centers.
2 ©2013 Univercity of Athens
3 Figures 5 - 22 are the author's copyright
So the length of r is r{6) = Jif + d2 -2rad ■ cos9a and it depends from angle 6>0. We can estimate the angular velocity of K as seen from O'. We have:
ro ==> ro = ro®odo, 0o = c»j f = rr => r = rr + rcoO, 0 = CO
r =r0-d => r = r0 => r0 ■ f0 = r-r = f2 + 7-V =>
/ (r <°nd ■ sin® t )2 / 9 ,9 , \ 9
O) = 2[y 1®o) t +(r2 + d - 2r0d ■ cos60 )0 r + d -2rd■ cos® t v '
CO
O = ■
(l - z cos coot) d
2
l + z -2zcos® t
z = — <^: 1 r„
The blue gear carries the slot, which is parallel to O'C. So the gear rotates with rn, too. If we expand in Taylor series, we get:
CO — COo+ cooz cos COot + COoZ2 COS2COot + .... (l)
A body that is being moved by a central force, generally follows a conic section. We are interested for an elliptical one. Let the focus (c, 0) be the tractive point and rQ the length of the semi-major axis. (Figure 6). In this case the body has a constant angular momentum pe = mr20 .
The angular velocity is then 6 = pe / mr2.
Figure 6. Motion in elliptical orbit.
w u l -e J l l
We have r =-ra and —; = —; +
l + e cos 6 r ro
where e is the orbits eccentricity.
2 cos 6
e +
e=0
2+cos26
e2 +.
-le=0
P
For e = 0, we take coa = —^-r = 6 and 6 = caj, so we can write:
mr
9(t,s)~co0 + elco0cosw0t-e2co0{l+cos2a/j+... (2)
2
r
u
The eccentricity of Moon's orbit is 0.0549 [6]. From Mechanism's measurements [7] we have:
1.1 [mm]
_d__ r0 9.6 [mm]
0.114583 = 2x0.0572917
That is to say that the relative offset of the two centers is estimated to be almost double than the eccentricity of Moon's orbit. Equations (1) and (2) coincide up to 1st order expansion.
Theoretical analysis of the two-gear system
Figure 7 represents an elliptical aperture and its positions as the gears rotate. The center of the Pin K moves along the half periphery of this ellipse, around C. Clockwise for 0<6<n and counterclockwise for n<Q<2n. Mind that it's r0 = 9.6[mm] and d = 1.l[mm] [7].
Figure 7. The positions of an elliptical aperture as the gears rotate. We are referring to the Oxy Cartesian Coordinate System (CCS). We consider
KA = 00 ' = d =
(d\ v0,
d
d
= constant vector. It is then 0'A = 0K=>—0'A = —OK. So O'A
dt dt
forms an angle 0o with the Ox semi axis and it moves with an coQ angular velocity, as OK. The angle between O'C and Ox semi axis, which we are interested for, is 0 = 0 + SO. It is obvious that SO > 0 at the upper semicircle and SO < 0 at the lower one. We also consider the rotating O'xy CCS. As can be seen in Figure 7 the pin always moves on the half-plane with positive y .
Since O'C is fixed on the blue gear, its angular velocity 0 = co0 +SO is the one we are looking
for. Note that for 0 = 0 it must be SO > 0 and for 0 = n it is SO < 0. So we need to describe the motion as seen from the rotating O'xy CCS. The Lagrangian of the rotating pin K at the Oxy CCS (we denote as m it's mass) is:
L = ^m(x2+ y2)-\mco20(x2+ y2)
We are going to rewrite the Lagrangian to a more "compact" form. We set q1 = x, q2 = y . We denote as the elements of the Euclidean metric, were i,j =1,2. We can use the dimensionless
variables: qk ~^qk / ro, d -^d !r0, so becomes \rQ I = 1.
1
1
We have then: L = —mniiq'qJ -—mca2niiq'qJ 2 j ^
We are going to rewrite the Lagrangian in terms of 5c, y coordinates of the rotating system. The transformation is:
f cos 6 - sin6^| sin# cos#
The new form of the Lagrangian is:
q'=R1(6)q1+d' where R1 (0) =
L = ^ mñ,/fc7 + ^mÖ2fiMqkq1 + mOR í ^ j cf'q1mco2ßklqkq' - ^me/ (d f -mco2R'j (0)q1di
Where (c/)2 = d-d . The term -ma2 (c/)2 /2 may be omitted because it is constant and does not
affect the equations of motion. We can also divide by the mass m, which is constant, too. We finally get:
L = \njkq' +\02nklqkq' + 0r[^] cfq' ~^co2onklqkq' -co2R[- (d)q'd,
d ( dL^ The equations of motion are —
dt
dL
. This leads us to:
J
qr + 0R ^ ] qk +2OR ^ ] qk = 02qr -co2oqr -cd2R(-Oj& d
k
v 2 J ¿ V 2 j k
7Z
We set d = 0.114583, mo = 2n [rad/sec], and 0 = ®ot + 56 [rad]. The period becomes T = 1 [sec].
Figure 8. Coordinates of K (5c,y) in the rotating O'xy CCS.
The equations of motion are then:
^t)-y(t)S0 -(¿25664 + 2S0)^(t)~x(t)(\2.5664 + 80)S0 + 4.52351 cos(27it + S0) = 0 ^ (í) + je (í) (l2.5664+ 2^) JE (í)-^(í)(l2.5664 + ^9) ¿0-4.52357 sin (2^ + á9) = 0
For a given shape of the aperture with known parametric equations i(/),>■(/) (Figure 8), we finally take two differential equations of the unknown function 06(f). We then solve the second equation for 36 and replace it at the first one, and we take a second order equation of SO .
-862 [x2 (t) + y2 (t)] + 80 {-12.5664 [i2 (t) + y2 (t)] + 25c (t) y (t) - 2x (t) (t)}
+4.523 57 [x (t) cos (2m + 8e)-y(t)ún(27it + 8e)\ + l2.566A(k(t)y(t)-x(t)^(t)) +x(t)'5c(t) + y(t)'Í;(t) = 0
The algebraic solutions are two non-linear differential equations of first order 8012=f12 (*, x, x, y, y, y, 80, t). One of them is going to be solved numerically.
For Ellipse the parametric equations verify the
(i-l^+J^.
where e is the aperture's eccentricity. It is also x2 + y2 =r2 = 1 + d2 — 2d• cosco0t.
From these two equations the accepted solution for x(t) is:
jc(i) = --—+ — + d2s4 ~ 2s2d■ cos(0ot
The solutions for y(i) are y(t) = + d2 - 2d ■ cos o0t - x2 (t), but since we want y > 0 the corresponding equation is:
y (t) = + d2 - 2d ■ cos coat - x2 (t)
As an estimate, based on Figure 4, the eccentricity of the elliptical aperture is about 0.98. We set semi-major axis d = 0.114583, e = 0.98, rn0 = 2n [rad/sec] and we get the following graphs of Figure 9:
Figure 9. The parametric equations i (/),>> (7) and the time derivative y(t}.
Equation x(t) varies between jc(0) = 0.885417 (or 8.5 mm) and jc(0.5) = 1.11458 (or 10.7 mm). Equation y(t) varies between ,y(0) = 0 and j(0.5) = 0.0228018 (or -0.22 mm). The time derivative y(t) shows a gap at t = kT/2. This occurs when 0 = кж, к<=Ъ. From the parametric equations for initial conditions tQ =0.00001, 80(to) = 0.00001, we get 80, (ta) - -13.551 and 802(to) - 0.641597.
Figure 10. The modulation of the gear's angular velocity for the cases y > 0 and y < 0.
At 0 = 0 <50 (t) is positive, so we are going to solve the second differential equation the
S02=f2( x,x,x,y,y,y,S62,t}. We find a numerical solution for 80(t) and we plot the 80(t) at
Figure 10. In the same Figure has been also plotted the case for which it would always be y < 0. The graph also shows a gap at 0 = kn, so Ellipse may not be a suitable shape for the aperture because it produces an intermitted movement of the gear.
Since we wanted to restrain the movement at the positive half axis, the Lagrangian should had been written as:
1 ~ , 1 />2„~ ~k~J , nn\ Я I 7.1 1 ~k~l 2~ ~k
L = ~nHcfq' +-0\!qk'q1 \ Я Я ~<nHqk R(~e)'sds-cQ(-cf)
There also should be a potential -c0(-v) where c is a suitable positive constant. The corresponding equation of motion should then be:
1 = 4.52357 sin (2^ +^6») + y (l2.5664 + 50)50- (l 2.5664 + 2SO)£- xSO + cS(y)
There exists a force c8(y) that causes the reflection at point v = 0. We tacitly assumed the existence of a potential -c©(-v) when we selected the positive solution for y . The solution of the differential equation 562 = f2(x,x,x,y,y,y,S62,t) also shows a discontinuity.
Slot is a limit case of Ellipse with e = 1 and parametric equations:
I x(t) = y]l + d2 - 2d ■ cos
m„ t
For initial conditions t = 0, 80(0) = 0, it is S01 (0) - -13.3774 and S02 (0) - 0.811045 . So we will solve again the Sd2 = f2{x,x,x,y,y,y,ô02,t). The solution is drawn in Figure 11 (blue) along with the Ellipse's curves. In case of Slot, y(t) = 0 is a continuous function. So the solution of 5û( t ) is continuous as well.
Figure 11. The modulation of the gear's angular velocity for Ellipse ( y > 0, y < 0 )
and for Slot ( v = 0 ).
Let us examine now Parabola. We consider the same length of axes as in Ellipse. The parametric equations are:
X = l-£/-COS(Z>(V)
y = d ^Jl-s2 sin2 (p(t) where p( t) is a nonlinear function of time. They lead us to the relation
(£-if
d2
y
d-y/l
-s
d2
(3)
We have:
x2 + y2 = 1 + d1 - 2d ■ cos coct
y2=d2(\-s2)
d2
ll
x~ +
d2(\-s2)
d2
= \ + d2 - 2d-cos coj (4)
The above equation (4) has four solutions. We call the one of interest xs (/). It shows a time delay of T/2 from the expected one. So our final result is i(/) = i,s (/ + 0.5) and is drawn below. This also varies between 0.885417 (or 8.5 mm) and 1.11458 (or 10.7 mm) and the corresponding y(t) oscillates between 0 and 0.0228018 (or 0.22 mm) as we expected.
Figure 12. Parabola's parametric equation i(/). The graph of the equation (3) is printed in Figure13.
2
2
Figure 13. Parabola''s parametric equation >'(/)
Again, we have to solve the equation 862 = f2 {x,x,x,y,y,y,Sd2,t}. Finally, we get a numerical solution for the angular velocity SO(7) which is shown below (Figure 14):
Figure 14. The modulation of the gear's angular velocity for Parabola.
Figure 15 shows a closed trajectory in the phase space. This means that we have a periodic curve, as we expected.
Figure 15. Closed trajectory in the phase space for Parabola. The time derivative _y(7) is continuous (Figure 16).
Figure 16. Time derivative >'(/) for Parabola.
This is also confirmed by calculating the limits at points t = 0 and t = 0.5:
p (0) = -3.9633 • 1 (T17, limjp(O) = -3.9633-1 (T17, lim p (0) = -3.9633 • 1 (T17 ^(0.5) = -6.29641 -lCT17, lim ^(0.5) = -6.2964MCT17, lim j5(0.5) =-6.2964MCT
/—>0.5
Since the trajectory in phase space is closed, _y(7) is continuous for every 0 = kn, k e Z. The step function c©(-j>) is not needed to be present in Lagrangian's formula any more and this yields a continuous solution. The two shapes of the aperture are printed below (Figure 17):
^ Parabol?
Ellipse
f f 0.01Í \\
/ / 0.01 \\
(/ 0.005 1_Li_1_._._._._1_._._._._ _________1_________1_■'■" ■
0.S 0.S5 1.05 1.1 j(J-)
Figure 17. The shapes of an elliptical and a paraboloid aperture.
Note that the difference in the final results is due to the form of the shape, as can be seen in Figure 18 which shows the three cases synoptically.
Figure 18. The modulation of the gear's angular velocity for Slot (green), for Ellipse ( y > 0, violet)
and for Parabola (blue).
Convert units to deg/day
The anomalistic period of the Moon is Tan — 27.55455 [days]4 and the corresponding frequency
360°
1S ®an =
2n
-13.065 [deg/d] or coan= — - 0.228027 [rad / d].
rpT L *—' J u/l y-T
an an
We said before that 0(l) = coQ + S6(coJ) where we had set coa = In [radsec]. The expression of S0(coj) includes terms of coscoj where the argument is expressed in rad. We have to replace the angular velocity coa = In [radsec] by the anomalistic velocity 13.065 [tfeg7<s/] and the period T= 1 [sec] in SO(oot) by Tan. So we have:
ê(t) = 27r + ôê(27rt) [rad / sec] —» 6
deg
day
= 13.065 + -
13.065
2n rad
_ sec _
m 2n-
27.55455
The curves below can now be compared with the data:
4 http://eclipse.gsfc.nasa.gov/SEsaros/SEsaros.html
Figure 19. The modulation of the gear's angular velocity for Slot (green), for Ellipse (y > 0, violet) for Parabola (blue) and the anomalistic frequency (red).
Comparison with astronomical ephemeris data
The data cover the period from 1/1/2009 to 31/12/2013. The calculations were performed using an astronomical computer program HORIZONS System5. For each day we have the right ascension RA [h, m, s] and the declination S [deg, m, s]. In order to compute the angular velocity we need to consider the celestial sphere (Fig. 20).
RA\h\
We also consider a spherical coordinate system where 6 = ^ 15 [deg] and ( = 90 -8
RA\h\ n\rad \ n n
[deg] or 6 = ^V 15° —--± [rad] and (p = n-8°[rad].
1\ h\ 180° 2 180°
The position vectors of the moon in Figure 20, are:
rx=R (sin (px cos 0X sin (px sin 0X cos (px)
f2= R{sm(p2cosd2 sm(p2smd2 cos<z»2) where R is the radius of the celestial sphere.
Let u be the angle between these two vectors. It is then:
fx-f2 = R2 cosu = R2 [cos(px cos<p2 + cos(Qx -02)sin<px sin<p2]
u = arccos (cos ( cos (p2 + cos (6X - 6 ) sin ( sin ( ) From the sourced data we can now calculate the angular velocity in deg/day or in rad/day.
Y i+2
We then compute the moving average of five terms: u = _ ^ u .
5 j=i-2 '
The results are represented in dots along with Parabola's curve (blue) and Slot's curve (green) (Fig. 21). The computed mean values of minima and maxima of the ephemeris data are ^ = 11.8883 [deg/d] and = 14.689 [deg/d]. They are shown in orange in Figures 21.
5 http://ssd.jpl.nasa.gov/horizons.cgi#results
Celestial North Pole
Celestial South Pole
Figure 20. A schematic illustration of moon's trajectory in the Celestial Sphere.
b
e
Figure 21. Comparison of the theoretical curves of Slot (green) and Parabola (blue) with the ephemeris data (dots), the anomalistic frequency (red) and the mean values of minima and maxima (orange). a - year 2009, b - year 2010, c - year 2011, d - year 2012, e - year 2013.
Now we are going to calculate the distribution of the deviation of the theoretical curves of Slot and Parabola from the computed moving averages:
deVSlot =-
(t )
100, devPar =
{t, )
100
ut ui
The statistics are:
Slot: Mean = 0.970374, Median = 1.3069 , Standard Deviation = 2.0113, Skewness = -0.664744 Parabola: Mean = 0.973558, Median = 1.13315, Standard Deviation = 2.54583, Skewness = -0.428225
Parabola's distribution shows a slightly larger standard deviation. The Histograms are given in below Figures 22a, b. In the above analysis of Parabola we assumed that the ratio of the minor
semi axis (a) to the major semi axis (b) is \jl-s2 — 0.199 which corresponds to £• = 0.98. If a = 1.1mm then b = 0.22mm. The question is which was the lower level of that age's micro-technology.
Histogram for Parabola
M
a u Q
O 1H
CL
mean s deviation %
b
Figure 22. Histogram of the percent deviation of the theoretical curve from data:
a - for Slot, b - for Parabola.
If we take into account that the pin has an unknown diameter then the major semi axis remains the same, since it is the distance between the center of gear k1 and the center of the pin K (as discussed earlier). But then the minor semi axis must be smaller. This means that the corresponding s must be larger than 0.98. Furthermore, according to the geometrical analysis, since the eccentricity of moon's orbit is ~0.0549 [6], it should be d = 2x0.0549 = 0.1098.
a
For example, if we have s = 0.99 and d = 0.1098, the statistics are: Slot: Mean = 0.9442, Median = 1.0994, Standard Deviation = 1.8914, Skewness = -0.4178 Parabola: Mean = 0.9463, Median = 1.1076, Standard Deviation = 2.16, Skewness = -0.3664 The curves now seem to fit better with the data. The Mean and the Standard Deviation are smaller than the previous case. Slot still has a smaller standard deviation than Parabola.
The mean distances at perigee and apogee are 363300 [Km] and 405500 [Km] and the mean maximum and minimum orbital velocities are 1.076 [Km/sec] and 0.964 [Km/sec]6 respectively. We can compute the mean angular velocities at these points. They are <n = 11.7685 \deg / day\
and <x = 14.6617 \deg / day\. Remind that the computed minima and maxima of the data were
11.8883° / day and = 14.689° / day.
In the first case (where we had set s = 0.98 and d = 0.114583) the minima and maxima of Slot and Parabola were
<0.1 = 14.7558° /day, co^ = 11.7219° /day, = 14.8627° /day, <£1 = 11.6229° / day
In the second case (s = 0.99 and d = 0.1098) the minima and maxima of Slot and Parabola were <2 = 14.6765° / day, <°n 2 = 11.7724° /day, 2 = 14.7307° / day, <£2 = 11.7207° / day
The ratios min/max are:
— „slot par „slot „par
CO CO — , — , — o — T
rmL = 0.803, —= 0.809, = 0.794, ^^ = 0.782, = 0.802, ^^ = 0.796
c — —do' —pa\ —slot7 —par7
max max max, 1 max, 1 max, 2 max, 2
In both cases Slot's value is closer to 0.803.
1 - e2
In the geometrical analysis for an elliptical orbit we had found: r =-ro. If rmm, 6m..
1 + e cos 6
are the distance and angular velocity at the perigee (where 6 = 0) and rmax, 0mm are the distance and angular velocity at the apogee (6 = n), due to angular momentum's conservation law it is
/¿„¿L = rlJ_ ^ 0. / ¿L =rlJrL
So: ^ = — ^
r„„„ 1 + e
f„ VA „A2
r
mm
r
max V max J
For e = 0.0549 it is
f V
r .
mm
r
V max J
1 - e 1 + e
- 0.803.
= 1 -4e + 8e - 12eJ +.... (5)
CO (l - z cos 0 ) 0 1 - Z
The equation of Slot was co =-^-— . So at points 6a = 0 and 6a = it is =-
1 + z - 2z cos 0„ 1 + z
max
For z = 2e, = 1 -4e + 8e2-16e3 +.... (6) 1 + 2e
For e = 0.0549, it is 6>mn/6^ - 0.802 .
,nd
Equations (5) and (6) coinside up to 2 order expansion. Slot's results are more accurate.
6 http://nssdc.gsfc.nasa.gov/planetary/factsheet/moonfact.html
Conclusions
The above theoretical analysis illustrates the dynamics of the pin-aperture system. It is important to have an accurate measure of the eccentricity of Luna's orbit. The shape of the aperture introduces a potential which forces the gear to a differential rotation that simulates Lunar motion. The differences between the shapes are important especially at the minima and maxima The advantage of Slot is that it gives more accurate results and is easy to be manufactured. Other shapes with known parametric equations could also be examined in order to compare the produced motion with astronomical ephemeris data.
Acknowledgements
I wish to thank Prof. Xenophon Moussas for fruitful discussions and the provision of bibliography and photos.
References
1. Freeth, T.; Bitsakis, Y.; Moussas, X.; Seiradakis, J.H.; Tselikas, A.; Magkou, E.; Zafeiropoulou, M.; Hadland, R.; Bate, D.; Ramsey, A.; Allen, M.; Crawley, A.; Hockley, P.; Malzbender, T.; Gelb, D.; Ambrisco, W.; Edmunds, M.G. Decoding the ancient Greek astronomical calculator known as the Antikythera Mechanism. Nature 2006, 444, 587-591.
2. Moussas, X. Antikythera Mechanism, PINAX (TABLET) the first computer and mechanical Cosmos (in Greek), Ed. Hellenic Physical Union, Athens, Greece, 2011 and 2012 2nd edition.
3. Freeth, T.; Jones, A. The Cosmos in the Antikythera Mechanism, ISAWPapers 4, 2012. http://dlib.nyu.edu/awdl/isaw/isaw-papers/47 (Access 05.02.2016)
4. Freeth, T.; Jones, A.; Steele, J. M.; Bitsakis, Y. Calendars with Olympiad display and eclipse prediction on the Antikythera Mechanism. Nature 2008, 454, 614-617 (partially available on http://www.nature.com/nature/journal/v454/n7204/extref/nature07130-s1.pdf).
5. Wright, M.T. The Antikythera Mechanism and the Early History of the Moon-Phase Display. Antiquarian Horology 2006, 29, 319-329.
6. Meeus, J. Mathematical astronomy morsels. Richmond, VA: Willmann-Bell, 1997; 11-12.
7. Gourtsoyannis, E. Hipparchus vs. Ptolemy and the Antikythera Mechanism: Pin-Slot device models lunar motions. Advances in Space research 2010, 46, 540-544.
© This article is an open access article distributed under the terms and conditions of the Creative Commons by Attribution (CC-BY) license (http://creativecommons.org/licenses/by/4.0/).
