Creator of the Symmetry454 CalendarThe Length of the Lunar Cycle

Bookmark or cite this page as <http://www.sym454.org/lunar/>

Analysis by Dr. Irv Bromberg, University of Toronto, Canada email icon

[Click here to go back to the Symmetry454 / Kalendis home page]

Menu of Topics:

  1. Overview
  2. Glossary of Terms
  3. Celestial Mechanics by Numerical Integration
  4. Reference Material
  5. Short-Term Periodic Variations of the Lunar Cycle
  6. A Close Look at the Short-Term Periodic Variations of the Lunar Cycle
  7. Fixed Arithmetic Lunar Calendars
  8. Constant Interval New Moon Estimate
  9. Periodic Lunar Cycle Variations Relative to a Fixed Lunar Calendar
  10. Which Point in the Lunar Cycle has the Greatest Amount of Periodic Variations? new section
  11. Lunar Half-Cycles: Triple the Periodic Variations!
  12. New Moon minus Estimate
  13. The Mean New Moon Moment
  14. Estimation of Fixed Lunar Cycle Calendar Drift
  15. Moment to Lunation Number Conversion
  16. The Mean Synodic Month
  17. Mean Lunar Angular Motion (Mean Sidereal Month)
  18. The Rate of Change of the Mean Synodic Month
  19. Is There a Long-Term Effect of the Mean Earth Orbital Eccentricity?
  20. The Role of Earth Axial Tilt (Obliquity)
  21. Future Research

Overview

This page examines the periodic variations in the length of the lunar cycle, and its long-term trend over a 19000-year period.

Expressions will be presented for computing the moment of the mean lunar conjunction, the mean synodic month, and the rate of change of the mean synodic month.

All of the charts are in GIF format but can be viewed in more detail and higher resolution by clicking the chart, which will open the corresponding Adobe Acrobat PDF (Portable Document Format) file in a new window. To obtain the freeware Acrobat Reader, click here.

Glossary of Terms

Terrestrial Time (TT)
Uniform time base running at the same rate as International Atomic Time at the Prime Meridian (longitude 0°) on the surface of the rotating geoid (surface of Earth, at sea level). Astronomical calculations are typically carried out in terms of TT in order to avoid uncertainties due to irregular, periodic and long-term variations of the Earth rotation rate.
Universal Time (UT)
Mean Solar Time at the Prime Meridian, running at a rate such that each mean solar day has exactly 86400 mean solar seconds.
Delta T = ΔT
ΔT = TT – UT, UT = TT – ΔT, and (approximately, but adequate for the next 10000 years) TT = UT + ΔT.
Typically, astronomical calculations are carried out in terms of TT, then converted to UT by subtracting ΔT.
In order for Delta T to work optimally for the purposes suggested herein, its approximation should be a continuous function of the time elapsed relative to J2000.0, without any annual or monthly stepwise granularity. Suitable expressions (published January 2007) can be found at the NASA Eclipses web site at <http://eclipse.gsfc.nasa.gov/SEhelp/deltatpoly2004.html>, except that to avoid monthly granularity in the Delta T approximation use the following expression when calculating y (the fractional year number):

y = 2000 + ( TTmoment – J2000.0 ) / MARY

where TTmoment is the Terrestrial Time moment and J2000.0 is as defined below, both in terms of the number of days and fraction of a day elapsed relative to a specified ordinal day numbering epoch, and MARY (Mean Atomic Revolution Year) = 365+31/128 atomic days, as explained on "The Lengths of the Seasons" at <http://www.sym454.org/seasons/>. Over the entire range from 500 BC to 2050 AD, this modification never causes more than 4/5 second of difference compared to the unmodified arithmetic of the NASA algorithm.
The studies documented herein were carried out before NASA published the new Delta T approximation.
Click here for Delta T comparison charts of Dershowitz & Reingold (used herein) vs. Meeus Espenak (new at NASA) 42KB.
Click below to access Robert Harry van Gent's web page explaining Delta T and the history of Delta T approximations:
<http://www.staff.science.uu.nl/~gent0113/deltat/deltat.htm>.
J2000.0
January 1, 2000 AD at Noon, Terrestrial Time = Julian Day 2451545.0 = Rata Die 730120.5 (days since Gregorian epoch).
Ecliptic
Path of Earth's orbit around Sun, projected onto the celestial sphere.
Earth Axial Tilt = Obliquity
The angle between the Earth's equator and the ecliptic, or the angle between Earth's axis of rotation and a line perpendicular to the ecliptic.
Lunar Conjunction = New Moon = Dark Moon
Moment when the celestial longitudes of Moon and of Sun are equal when projected onto the ecliptic. It is only rarely that at a lunar conjunction the centers of Sun, Moon, and Earth, in that order, are aligned in a straight line. Such an alignment only occurs during the maximum of a total solar eclipse at the moment that the center of the lunar umbra (darkest part of the lunar shadow) crosses the terrestrial latitude that is equal to the solar declination. The lunar orbital plane is tilted about 5° 9' relative to the ecliptic (mean value), therefore most of the time the lunar latitude is either north or south of the ecliptic. If the lunar latitude is >32 arcminutes (32' = the solar angular diameter ≈ lunar angular diameter) north or south of the ecliptic at the moment of a lunar conjunction then Moon passes north or south of Sun, respectively, without any solar eclipse. Nevertheless, one could properly describe the lunar conjunction as the moment when the centers of Sun, Moon, and Earth, in that order, are all in the plane that is perpendicular to the plane of Earth's orbit (the ecliptic plane).
Mean Lunar Conjunction
Average lunar conjunction moment, assuming that Moon's orbit is a circle with Earth at its center, and that Earth's orbit is a circle with Sun at its center.
Length of the Lunar Cycle = (Duration of the) Lunation = Synodic Month = Lunar Month
Time elapsed from one lunar conjunction until the next.
Mean Lunar Cycle = Mean Synodic Month (MSM)
Average time elapsed from one mean lunar conjunction until the next. Its long-term variations will be examined herein.
Lunation Number
A count of lunar months before or after a specified epoch. There are approximately 12368 lunations per millennium.
The study presented here counts mean lunations relative to the mean lunar conjunction nearest to J2000.0, which was on January 6th, 2000 AD, taken as lunation number zero.
Lunations prior to that epoch are negative, and lunations after that epoch are positive lunation numbers.
Fractional values of 0, 0.25, 0.5, and 0.75 refer to the Mean New Moon, Mean 1st Quarter, Mean Full Moon, and Mean 3rd Quarter, respectively.
Add 71233 to convert J2000 lunation numbers to traditional Hebrew calendar lunation numbers (count of months elapsed since the Hebrew calendar epoch), such as are used to calculate the traditional molad.
Add 24724 to convert J2000 lunation numbers to nth-new-moon lunation numbers, as used in the algorithms of the book Calendrical Calculations, by Reingold & Dershowitz.
Add 953 to convert J2000 lunation numbers to Brown lunation numbers, for example as quoted by the Royal Astronomical Society of Canada (RASC) in their annual Observer's Handbook, which are relative to January 1923 AD, based on the series described by Ernest W. Brown in Planetary Theory (1933).
Add 17038 to convert J2000 lunation numbers to the Islamic lunation number since the epoch of the Observational Islamic calendar. This conversion is approximate, because the J2000 lunation number refers to the mean lunar conjunction at the Prime Meridian whereas the Islamic lunation number refers to the first visible new lunar crescent, whose timing is always later and also varies from locale to locale.
Add 37105 to convert J2000 lunation numbers to Goldstine lunation numbers, as used by Herman Goldstine in his 1973 book New and Full Moons 1001 B. C. to A. D. 1651.
Lunar Apogee and Perigee
Apogee is the point in Moon's elliptical orbit where it is furthest from Earth.
Perigee is the point in Moon's elliptical orbit where it is closest to Earth.
Both points advance eastward around Earth (in the same direction as Moon revolves), completing a full cycle in about 412 days or almost 14 lunar months. NASA has posted a picture showing the apparent lunar diameter at apogee and perigee, see: <http://antwrp.gsfc.nasa.gov/apod/ap071025.html>.
Aphelion and Perihelion
Aphelion is the point in Earth's elliptical orbit where it is furthest from Sun, currently near July 3rd.
Perihelion is the point in Earth's elliptical orbit where it is closest to Sun, currently near January 3rd.
Both points advance through the entire calendar year, completing a full cycle in about 21000 years (faster as orbital eccentricity declines).
The season that contains aphelion is longer because Earth moves slower near aphelion. When aphelion is near an equinox or solstice, the two adjacent seasons are both relatively long.
The season that contains perihelion is shorter because Earth moves faster near perihelion. When perihelion is near an equinox or solstice, the two adjacent seasons are both relatively short.
NASA has posted a picture showing the apparent solar diameter at perihelion and aphelion, see: <http://antwrp.gsfc.nasa.gov/apod/ap070709.html>.

Celestial Mechanics by Numerical Integration

To obtain the accurate moments of the actual astronomical lunar conjunctions, as the basis for evaluating the lengths of the lunar cycle and how it changes over time, and to prepare the collection of charts offered below, I used numerical integration, which is arguably the "gold standard" for celestial mechanics, and which is easy to do using SOLEX version 9.1β (or later), a "postcard-ware" computer program written by Professor Aldo Vitagliano of the Chemistry Department at the University of Naples, Italy. Version 9.1β was the first version which could automatically find lunar conjunction moments, logging those moments to a "MINDISH.OUT" text file. The final version 9.1 was released in January 2007 and is available at the SOLEX web site.

The SOLEX integration was carried out in terms of Terrestrial Time (usually abbreviated TT but indicated as TDT within SOLEX), with Delta T switched off and the geographic locale set to the Equator at the Prime Meridian, starting from date January 1, 2000. SOLEX integrated forward at 1-day intervals to beyond the year 12000 AD, and then starting again from year 2000 SOLEX integrated backward at 1-day intervals to before the year 7000 BC.

According to the SOLEX documentation, its numerical integration takes into account:

SOLEX Limitations:

For more information about SOLEX and to download the program please see its web page at <http://main.chemistry.unina.it/~alvitagl/solex/>.

For those charts that are in terms of mean solar days, I assumed steady tidal slowing of the Earth rotation rate such that mean solar days get longer by 1.75 milliseconds per century. Of course the actual Earth rotation rate does not slow down at such a perfectly steady rate, but goes through short-term fluctuations and long-term periodic cycles, many of which are unpredictable with our present state of knowledge. Tidal slowing was probably greater in the past when the polar ice caps were more massive with lower sea levels and axial tilt was greater than in the present era, and tidal slowing will probably diminish over the coming several millennia due to global warming (reduction of polar ice mass, rising sea levels) and due to declining axial tilt.

Reference Material

Jean Meeus published and excellent chapter entitled "The Duration of the Lunation" on pages 19-31 (chapter 4) in "More Mathematical Astronomy Morsels" by Jean Meeus, published in 2002 by Willmann-Bell, Inc., Richmond, Virginia. Herein, when I mention Jean Meeus without further qualification, it refers to what he wrote in that chapter.

Another relevant chapter by Jean Meeus is chapter 33, "Long-period variations of the orbit of the Earth" in More Mathematical Astronomy Morsels, pages 201-205, published in 2002 by Willmann-Bell, Richmond VA, in which he presented his adaptation of the 18-term astronomical algorithm of Pierre Bretagnon (1984) for determining the mean Earth orbital eccentricity over a ±million year range. For more information and charts of the Earth orbital eccentricity, based on Meeus' adaptation of the Bretagnon algorithm, please see my "Lengths of the Seasons" web page at <http://www.sym454.org/seasons/>.

Short-Term Periodic Variations of the Lunar Cycle

In the present era the median length of the lunar cycle is about 29d 12h 30m, the average (MSM) is slightly more than 29d 12h 44m, the shortest lunations are about 29d 6h 30m, and the longest are about 29d 20h. Thus the length of the synodic month varies over a range spanning about 13h 30m. These variations were greater in the past and will diminish in the future:

Centile trends, per group of 4657 lunar months, based on SOLEX 9.1β numerical integration

I will explain later why I chose the "magic" number 4657 for averaging the lengths of lunar cycles.

The chart above happens to be in terms of Terrestrial Time (TT), but there would not be any discernible difference if the data had been calculated in terms of Universal Time (UT), because of the very wide range of the Y-axis scale.

A Close Look at the Short-Term Periodic Variations of the Lunar Cycle

The following chart shows the variations obtained when the actual length of each lunar cycle is subtracted from the actual length of the mean synodic month, as calculated by SOLEX in terms of Terrestrial Time. Each plotted point is a lunar conjunction, with more than 18000 plotted in total:

The next chart takes a closer view at the periodic variations by zooming into ±333 lunations relative to the J2000 epoch:

The next chart takes an even closer view at the periodic variations by zooming into ±70 lunations relative to the J2000 epoch:

Fixed Arithmetic Lunar Calendars

A variety of arithmetic lunar calendars are in use that are based on a repeating fixed lunar cycle duration. Some examples, sorted in descending order of the assumed lunar cycle length, are listed below, along with a good contemporary estimate of the mean synodic month:

Description Assumed Lunar Cycle Length
(exact days)
Decimal Month (days)
overscored groups repeat
Time in excess
of 29 days
Orthodox Easter computus (365+1/4) × 19/235 = 29+499/940 ≈ 29.530851 12:44:25+25/47s
Hebrew calendar molad 29+12/24+44/1440+1/25920 = 29+13753/25920 29.530594135802469... 12:44:03+1/3s
Yerm Lunar calendar 25101/850 = 29+451/850 29.530588235294117647... 12:44:02+14/17s
Hindu Surya calendar (modern) 29+7087771/13358334 ≈ 29.530587946 ≈ 12:44:02.8
Mean Synodic Month
(2000 AD, mean solar days)
29+82517/155520 ≈ 29.53058770576 12:44:02+7/9s
Tibetan Phugpa calendar 29+3001/5656 29.530586987270155... 12:44:02+506/707s
Gregorian Easter computus 2081882250/70499183 = 29+37405943/70499183 ≈ 29.53058690056 ≈ 12:44:02.7
Cycle of 49 yerms 25101/801 = 29+425/801 29.530588235294117647... 12:44:02+62/89s
Cassidy-Dee Easter computus 48091470/1628531 = 29+864071/1628531 ≈ 29.530583083773 ≈ 12:44:02.4
25 Saros cycle 29+2958/5575 ≈ 29.5305829596412556 12:44:02+82/223s
Hindu Arya calendar (old) 1577917500/53433336 = 29+2362563/4452778 ≈ 29.53058180758 ≈ 12:44:02.3
Fixed Islamic calendar (30×6×59+11)/(30×12) = 29+191/360 29.5305... 12:44:00

Those that are listed above the year 2000 AD mean synodic month have assumed lunar cycle lengths that are too long and are drifting late, and those that appear below it are too short and are drifting early, relative to the present era mean lunar cycle. For optimizing the long-term drift of a fixed interval lunar calendar, however, it is better to choose a mean month that is slightly too short, because the mean lunar cycle is getting progressively shorter in terms of the mean solar days that are relevant to calendars.

Regardless of the assumed length of the lunar cycle, the short-term periodic variations of actual lunar conjunction relative to mean lunar conjunctions calculated using any of the above fixed intervals will be about double the periodic variations of the length of the actual lunar cycle, as will be shown next.

Constant Interval New Moon Estimate

To evaluate the variations of the lunar cycle relative to a fixed lunar calendar cycle, we will examine those variations relative to estimated New Moon moments that are uniformly spaced at constant time intervals of 29 days 12 hours 44 minutes and 2+7/8 seconds (atomic time) per lunation, counting the lunations relative to zero near the J2000,0 epoch on Gregorian January 6, 2000 AD at 14:20:44 TT. (The actual lunar conjunction was at 18:14:42 TT, according to SOLEX, but our estimate is intended to relate to mean lunar conjunction moments.)

The choice of J2000.0 as the lunation epoch is arbitrary, but is convenient because today most modern astronomical algorithms are calculated relative to J2000.0, and because the use of a near-present-era epoch optimizes the floating point arithmetic accuracy by maximizing the number of significant digits to the right of the decimal point.

To compute the fixed day number of the New Moon estimate, relative to J2000.0, in terms of Terrestrial Time:

NewMoonEstimateTT( Lunation ) = MeanNewMoonTT_J2000 + MSM_TT_J2000 × Lunation + J2000

where Lunation is the lunation number relative to zero = January 6, 2000 AD, and the following are constants:

MeanNewMoonTT_J2000 = 5 – 1/2 + 14/24 + 20/1440 + 44/86400
(The 1/2 day is only deducted so that the time components can refer to midnight instead of noon.)
When this value is added to J2000 as above it represents the epoch mean lunar conjunction moment on January 6, 2000 at 14:20:44 TT.

MSM_TT_J2000 = 29 + 12/24 + 44/1440 + (2+7/8) / 86400
This is the approximate MSM at J2000.0, but it doesn't need to be exact for our purposes because any error will be cancelled out later.

Fractional Lunation values of 0, 0.25, 0.5, and 0.75 refer to the Mean New Moon, Mean 1st Quarter, Mean Full Moon, and Mean 3rd Quarter, respectively.

Periodic Lunar Cycle Variations Relative to a Fixed Lunar Calendar

The following chart shows the variations obtained when the constant interval New Moon estimates are subtracted from the corresponding the actual lunar conjunctions, as calculated by SOLEX in terms of Terrestrial Time. Each plotted point is a lunar conjunction, with more than 18000 plotted in total:

The next chart takes a closer view at the periodic variations by zooming into ±333 lunations relative to the J2000 epoch:

The next chart takes an even closer view at the periodic variations by zooming into ±70 lunations relative to the J2000 epoch. This chart shows the difference between the actual lunar conjunction and our constant interval estimate (fixed lunar calendar, black line with "x" symbols), and for comparison also shows the difference between the actual length of each lunar cycle and the mean synodic month (blue line with solid dots):

The rather large periodic variations of the length of the lunar cycle make it impossible to accurately calculate the mean length of the lunar cycle based on the separation between any two well-established lunar conjunctions, such as total solar eclipses, even when that separation spans many centuries.

Which Point in the Lunar Cycle has the Greatest Amount of Periodic Variations?

As impressive as the periodic variations at the lunar conjunction are, that is actually the point in the lunar cycle that has the least periodic variations! The most variations are actually found at the first quarter Moon (lunar quadrature, lunar phase = 90°).

To investigate this, I used the newer version 11 of SOLEX, which has the ability to automatically report the moments of each lunar quarter (actually the moments of conjunctions, quadrature, and opposition of any object). I ran SOLEX 11 in DE421 mode with extended 80-bit precision, 18th order integration, and forced solar mass loss, with automatic searching for lunar conjunctions, quadratures, and oppositions, integrating backward and forward from the present era over the date range from November 30, 1815 through September 17, 3288 AD with Delta T switched off (atomic time), a total of 18216 lunar months, during which time the lunar orbital nodes retrogressed four times 360° with respect to the Earth orbital perihelion. Then I calculated the duration of the lunation measured from conjunction-to-conjunction, quadrature-to-quadrature, and opposition-to-opposition, and found the maximum, minimum, and range of variation in each case:

Lunar Quarter and Lunar PhaseMaximum
>29 days
Minimum
>29 days
Range
(max-min)
Comment
(degrees of ecliptic elongation from Moon to Sun)(HH:MM:SS)(HH:MM:SS)(HH:MM:SS) 
Conjunction (New Moon = 0°)19:54:5206:33:4013:21:12least variations, about 13+1/3 hours
First Quadrature (1st Quarter = 90°)22:12:0904:13:1117:58:58almost 18 hours of variations
Opposition (Full Moon = 180°)19:57:4806:34:1913:23:29a tad more variable than conjunctions
Last Quadrature (3rd Quarter = 270°)22:12:5404:14:0417:58:51almost 18 hours of variations

The variations at oppositions are only 00:02:17 greater than at conjunctions.

The variations at first and last quadrature are essentially equal (within 7 seconds) but are almost 4+2/3 hours greater than at conjunctions.

These findings are in agreement with those reported by Jean Meeus, although he evaluated the lunation variations for only the 20th and 21st centuries (see chapter 2, "About the extreme durations of the lunation" in Mathematical Astronomy Morsels V, published by Willmann-Bell, Inc., 2009).

A sampling of these variations are graphically depicted below (click on the chart title or the chart image to view a high quality PDF version):

Lunar Cycle Duration at All Quarters

Lunar Cycle Variations at All Quarters

In the chart above, the height of the tallest peaks is about an hour greater than the depth of the deepest valleys.

When peaks are tallest the opposing valleys are shallowest. Conversely, when valleys are deepest the opposing peaks are shortest.

Lunar Half-Cycles: Triple the Periodic Variations!

In connection with the intercalation of the Hebrew calendar, some traditional commentators in the Babylonian Talmud stated that the year must be intercalated by insertion of a winter leap month if otherwise the spring equinox will fall later than the waxing or "renewal" half of the lunar cycle, in other words later than the lunar opposition or full moon moment. To examine the astronomical practicality of such an intercalation criterion, I investigated the periodic variations of the lunar half-cycles, that is the duration of the waxing half-cycle from lunar conjunction to opposition (new moon to full moon) and the duration of the waning half-cycle from lunar opposition until the next conjunction (full moon to new moon).

For this evaluation I used the lunar phase algorithms that were published by Jean Meeus in "Astronomical Algorithms", second edition, published in 1998 by Willmann-Bell, Richmond, Virginia, USA and the NASA Espenak-Meeus Delta T polynomials as referenced above, to find sequentially each lunar conjunction and opposition for ±1000 lunations relative to January 2000 AD (that is March 2, 1919 to November 11, 2080 AD) and then calculated the length in days of the waxing and waning half-cycles and full cycles.

The following chart shows the distribution of half-cycle and full cycle lengths obtained, expressed in terms of length in days versus the cumulative centiles. Click here or on the chart to open a higher-resolution PDF version (which also includes the next chart). The left and right y-axes both span a 2-day variation range so that it is easy to compare half-cycles to full cycles:

The following chart, focussing on only ±98 lunations relative to January 2000 AD, shows why this is so: as each half-cycle gets longer the other half gets shorter by a similar about, cancelling about 2/3 of the full-cycle variations. Click here or on the chart to open a higher-resolution PDF version (which also includes the previous chart). Here the left and right y-axes both span a 3-day variation range, again so that it is easy to compare half-cycles to full cycles. The points for the lengths of the waxing half-cycles and full cycles are plotted against whole lunation numbers, but the points for the waxing half-cycles are plotted against half-lunation numbers:

In addition, Jean Meeus has shown that on average the first and last lunar quarters are each about 10 minutes longer than each of the second and third lunar quarters, in other words Moon spends on average 20 more minutes over the day side of Earth than it does over the night side of Earth, per lunar cycle (see pages 12-13 in chapter 1, Jean Meeus"About the phases of the Moon" in "Mathematical Astronomy Morsels IV", published in 2007 by Willmann-Bell, Inc., Richmond, Virginia). This means that although the use of fractional mean lunation numbers as suggested in several places on this page are a good estimate of the mean lunar conjunction (fraction = 0) and opposition (fraction = 0.5), the mean waxing quadrature (fraction = 0.25) tends to be about 10 minutes premature and the mean waning quadrature (fraction = 0.75) tends to be about 10 minutes late.

New Moon minus Estimate

To calculate mean lunar conjunction moments we need to average out (cancel) the short-term periodic variations. I tried averaging in groups corresponding to the strongest beat frequencies of 14 × 111 = 1554 lunations (this is also the least common multiple of 14 and 111), but the results were not as smooth as desired. Eventually I found that tripling the group sizes yielded excellent smoothing, optimal when the lunations were in groups of 4657 lunar months, which is also very close to twice the approximately 2277-lunation period of the weak periodicity.

The following chart shows the difference in days between the actual astronomical lunar conjunction (SOLEX New Moon) minus the Constant Interval New Moon Estimate as was defined above, in Terrestrial Time, for dates from 7000 BC to 12000 AD, a total of 220,000 lunar cycles, averaged in 49 groups of 4657 lunations, along with the minimum and maximum differences in each group. Each average is plotted as a "+" symbol at the middle lunation number of that group:

The polynomial near top center gives the least squares statistical regression line to the averages. Although it appears to be parabolic, over this range of dates a quartic (4th order) polynomial provides a superior fit. Although it is not very obvious, the span between Earliest to Latest gradually tapers from the past to the future, due to the decline of Earth's mean orbital eccentricity over the evaluated date range.

The Mean New Moon Moment

The polynomial above can be used to define an accurate function returning the Mean New Moon Moment in Terrestrial Time:

NewMoonAdjustTT( L ) = 3.5962433E-22 L4 – 7.799103E-17 L3 + 1.005115E-10 L2 + 2.867010E-08 L + 8.945687E-05

MeanNewMoonTT( L ) = NewMoonEstimate( L ) + NewMoonAdjustTT( L )

where L is the lunation number relative to J2000 (as above, fractional L values of 0, 0.25, 0.5, and 0.75 refer to the Mean New Moon, Mean 1st Quarter, Mean Full Moon, and Mean 3rd Quarter, respectively). There may be a performance and numerical stability advantage in rearranging the NewMoonAdjustTT polynomial, in fact any of the polynomials presented herein, into monomial form according to Horner's Rule, if the programming language will allow it, by factoring out the powers of the lunation number, replacing exponentiation with nested multiplication:

NewMoonAdjustTT( L ) = { [ (3.5962433E-22 L – 7.799103E-17) L + 1.005115E-10] L + 2.867010E-08} L + 8.945687E-05

One could combine the arithmetic of NewMoonEstimate( ) and NewMoonAdjustTT( ) into a single polynomial expression, but only the constant and linear terms are combinable, and for clarity I felt is best to retain separate functions.

This polynomial and the others presented on this web page are valid over the range of lunation number included in this study, that is -100500 to 123500. They can probably be pushed an extra 20000 lunations further to the past or future, but don't rely on them beyond that.

To convert the TT New Moon moment to UT for use with a calendar, simply subtract Delta T:

MeanNewMoonUT( L ) = ThisMeanNewMoonDeltaT( ThisMeanNewMoon )

where ThisMeanNewMoon = MeanNewMoonTT( L )

Optionally one can convert the UT New Moon moment to any desired time zone or reference meridian by adding the offset as the appropriate fraction of a day relative to the Prime Meridian.

The following chart confirms that the MeanNewMoonTT function generates accurate mean lunar conjunction moments across the full range of lunations evaluated in this study. It is quite apparent that over many millennia the heights of the Earliest and Latest extreme peaks gradually converge toward intermediate heights (due to declining Earth orbital eccentricity):

Of course the MeanNewMoonUT function performs equivalently, provided that the same Delta T function is used for the astronomical lunar conjunction moments.

Estimation of Fixed Lunar Cycle Calendar Drift

It is easy to use the MeanNewMoonUT function to estimate the drift of a fixed arithmetic lunar or lunisolar calendar with respect to the mean lunar cycle over any range of lunations:

MeanLunarCycleDrift( MSM, FromLunation, ToLunation )

where MSM is the assumed constant lunation interval, and the lunation numbers are relative to zero = January 6, 2000 AD.

StartMeanConjunction = MeanNewMoonUT( FromLunation )

EndMeanConjunction = MeanNewMoonUT( ToLunation )

ElapsedMonths = ToLunationFromLunation

RETURN MSM × ElapsedMonthsEndMeanConjunction + StartMeanConjunction

The result is returned in terms of the integrated total number of days and fraction of a day drift accumulated from the starting to the ending lunation numbers. This drift is over and above any drift that had already accumulated up to the moment of the starting lunation, and it is valid regardless of the selected time zone or reference meridian of the lunar conjunction moments.

For example, we can use this function to calculate the lunar cycle drift inherent in the molad (mean new moon estimate) of the traditional fixed arithmetic Hebrew calendar from its inception (era of Hillel ben Yehudah, in Hebrew year 4119 = 358 AD) until the present era (say Hebrew year 5768 = 2007 AD). First we define a constant offset for converting traditional Hebrew elapsed month numbers to J2000 lunation numbers, based on the Hebrew lunation number for Shevat 5760, which was the Hebrew month that started shortly after the January 6, 2000 lunar conjunction:

HebrewLunationAtJ2000 = 71233

Next, set the starting and ending lunations to Tishrei 4119 and Tishrei 5768, respectively, converted to J2000 lunation numbers:

FromLunation = 50933 – HebrewLunationAtJ2000 = –20300

ToLunation = 71328 – HebrewLunationAtJ2000 = 95

For the fixed interval MSM use the traditional molad interval of the Hebrew calendar, which is 29 days, 12 hours, and 44+1/18 minutes per month:

MoladInterval = 29 + 12/24 + ( 44+1/18 ) / 1440 = 765433/25920 = 29.530594135802469...

The MoladInterval can't be represented as an exact decimal number because the above overscored digits repeat forever, but its double precision floating point representation is certainly adequate to calculate the drift to an accuracy of a second:

Drift = MeanLunarCycleDrift( MoladInterval, FromLunation, ToLunation )

In this case ElapsedMonths = 20395 and the Drift = about 0.0682385 days or about 1 hour 38 minutes and 16 seconds. This implies that the molad reference meridian has been drifting eastward. Each 4 minutes of time drift corresponds to 1° of longitude drift, but it is easiest to multiply our fractional day Drift by the 360° of Earth rotation per full 24-hour period = 24.57° or almost 24° 34' of eastward drift relative to wherever the molad reference meridian of longitude was in the era of Hillel ben Yehudah. Each 360°/24 hours per day = 15° corresponds to one standard time zone, so the molad has drifted by nearly two time zones!

By choosing several ending lunation numbers at appropriate intervals, one can chart how the drift rate changes over time, from which it will be evident that the molad drift is accelerating quadratically. For further information, see my full analysis of the molad of the Hebrew calendar at <http://www.sym454.org/hebrew/molad.htm>.

A convenient feature to add to the MeanLunarCycleDrift function is to allow the user to pass the MSM in terms of the number of seconds in excess of 29 days, 12 hours and 44 minutes. If so, then for the molad example presented above the MSM could be passed as 60/18 or reduced to 10/3 seconds.

I used the arithmetic explained in this section to generate a Microsoft Excel spreadsheet that dynamically depicts lunar calendar drift rates for a wide array of fixed arithmetic lunar cycles from 8000 BC to 12000 AD. The spreadsheet employs a VBA (Visual Basic for Applications) macro that calculates the lunar calendar drift relative to a user-specified zero reference year, displaying the results graphically, optionally allowing the user to alter a Delta T multiplier. In order for the macro to run, you must have the full version of Excel, it won't run in the Excel Viewer environment. In addition, your Excel security settings must allow the macro to run, with or without your confirmation, as you prefer. To enable macros, use the Excel "Tools" menu, choose "Macro", slide over to "Security...", then choose the desired macro security level.

Click here to download the "Lunar Calendar Drift" spreadsheet 1.1MB

The shorter cycles in that spreadsheet are not particularly useful today, but in the future they could each take turns for a few centuries in sequence as the mean lunation interval gets progressively shorter, as can be seen in this progressive yerm era analysis for lunar conjunctions at the Prime Meridian, hence the reference years for the shorter cycles are preset to the distant future.

Moment to Lunation Number Conversion

Given any moment in time, it can be useful to find the corresponding J2000-relative lunation number.

I plotted every 1000 lunations from -100000 to +100000 against the mean lunation moment expressed in terms of J2000-relative Mean Atomic Revolution Years (MARY = 365+31/128 atomic days per year, for more information see "The Length of the Seasons" at <http://www.sym454.org/seasons/>), and then used least squares statistical regression to generate a quadratic polynomial for converting any given TT moment to the corresponding J2000-relative lunation number:

TTmomentToLunation( TTmoment ) = -5.367946E-10 MARYs 2 + 12.3682665 MARYs – 0.172522

where MARYs = ( TTmoment – J2000.0 ) / MARY
and MARY = 365+31/128 atomic days

The result is the lunation number relative to zero = January 6, 2000 AD. Any fractional component indicates the portion of that lunation number elapsed up to the given moment (the mean lunar phase), but one should not rely on more than 2 or 3 decimal points. Values near whole numbers are returned at each TT mean lunation moment. Actually the quadratic coefficient is very small, so one could omit it and adjust the constant coefficient slightly to make a linear expression instead, especially if the intention is to round the result to 2 or 3 decimal points anyway:

TTmomentToLunation( TTmoment ) = 12.3682665 MARYs – 0.184336

To find the lunation number given a UT moment, simply add Delta T and then pass the sum to the TTmomentToLunation function:

UTmomentToLunation( UTmoment ) = TTmomentToLunation( UTmoment + DeltaT( UTmoment ) )

The Mean Synodic Month

Given methods to compute mean New Moon moments, the average length of the lunar cycle, or mean synodic month (MSM) at any given mean lunar conjunction can be calculated as the average of the current and previous mean lunation interval:

MSM_Atomic( L ) = ( MeanNewMoonTT( L + 1 ) – MeanNewMoonTT( L – 1 ) ) / 2

MSM_Solar( L ) = ( MeanNewMoonUT( L + 1 ) – MeanNewMoonUT( L – 1 ) ) / 2

where L is the lunation number relative to J2000, and for the MSM_Solar function to work properly the underlying Delta T approximation must return "non-granular" Delta T values that vary as a continuous function of the time elapsed relative to some specified epoch.

Instead of using these MSM functions, to parallel the 4657-lunation averaging groups, I used the elapsed time between the starting mean lunar conjunction of one group and the start of the next group, divided by 4657, thus computing the MSM corresponding to the middle lunation number of each group. I then subtracted 29 days 12 hours and 44 minutes from each MSM value, and finally multiplied the residual fraction by 86400 (the number of seconds in a day), thus obtaining only the few excess fractional seconds for plotting.

The following chart shows the MSM as the number of seconds in excess of 29 days 12 hours 44 minutes in terms of both atomic time (getting progressively longer due to tidal forces accelerating Moon further away from Earth) as well as mean solar time (getting progressively shorter due to tidal forces slowing Earth's rotation more than the slowing of the lunar revolutions, a natural consequence of the gravitational transfer of angular momentum from Earth to Moon):

The cubic (3rd order) polynomials shown for each MSM regression line give the MSM as the number of seconds in excess of 29 days 12 hours and 44 minutes for any desired lunation number within the evaluated range.

Excess = LunationToMSM_Atomic( L ) = 1.242862E-16 L3 – 2.021546E-11 L2 + 1.7369075E-05 L + 2.877432

Excess = LunationToMSM_Solar( L ) = 1.2434254E-16 L3 – 2.021679E-11 L2 – 2.51203947E-05 L + 2.777861

where L is the lunation number relative to J2000.

To convert the Excess to days, simply divide it by the number of seconds in a day and add 29 days 12 hours and 44 minutes as follows:

MeanSynodicMonth = 29 + 12/24 + 44/1440 + Excess/86400

My Delta T function simplistically assumes an essentially steady rate of tidal slowing such that the length of the mean solar day will get longer by 1.75 milliseconds per century, but the actual remote past tidal slowing was probably substantially greater than that due to the historically greater mass of the polar ice caps with lower global sea levels and due to the greater Earth axial tilt that existed several millennia ago. The two MSM curves cross in the year 1810 AD, which was a few decades prior to the conventional minimum Delta T value, because the mean solar time curve is dominated by my remote past and distant future parabolic Delta T approximations.

If you prefer a different Delta T approximation then use it to derive your own mean solar time MSM polynomial instead of relying on the above.

Although the two lines seem to be reciprocals of each other, they aren't quite. In terms of atomic time, Earth's rotation rate is slowing down more than Moon's, as should be obvious by the fact that the red mean solar time trend line has a steeper slope, going all the way from the top left to the bottom right corner of the graph whereas the atomic time trend rises by only about 3/4 of the same scale. The tidal transfer of angular momentum varies with the arrangement of the continents, the depths of the oceans, the oblateness of Earth's not-quite spherical shape, the tilt of Earth's axis, the tilt of Moon's axis, and the inclination and eccentricity of Moon's orbit, as well as other factors.

The inverse function, returning the lunation number that has a given MSM, can be obtained by mathematical rearrangement of the above polynomials, but that yields a rather unwieldy equation, so instead it is much simpler to derive a direct polynomial by swapping the axes of the above chart and then similarly obtaining the inverse polynomial by least squares statistical regression. The user can pass either the desired MSM or the number of seconds in excess of 29 days 12 hours and 44 minutes (as written, the function considers any passed MSM value <29 to be just the excess seconds):

MSM_AtomicToLunation( MSM ) = -849.9472 Excess3 + 10174.37 Excess2 + 20156.62 Excess – 121610.9

MSM_SolarToLunation( MSM ) = -413.7623 Excess3 + 1822 Excess2 – 40469.25 Excess – 107460.75

where, in both cases, IF MSM > 29 THEN Excess = ( MSM – 29 – 12/2444/1440 ) / 86400 ELSE Excess = MSM

As above, the mean solar time polynomial is of uncertain reliability beyond the valid date range of the underlying Delta T function (typically from 1600 AD to the present era).

Mean Lunar Angular Motion (Mean Sidereal Month)

The length of the mean synodic month (MSM) divided into 360° yields the mean daily change in lunar phase, which is the angular difference between the lunar and solar ecliptic longitudes. This indicates the mean rate at which the selenographic colongitude of the lunar sunrise terminator progresses across the lunar surface, or the sunset terminator, which equals the selenographic colongitude +180°.

For example, based on the above arithmetic the MSM at Lunation zero on January 6, 2000 AD was about 29.5305877 mean solar days. Dividing that into 360° yields a mean change in lunar phase of about 12.19° or 12° 11' 27" per mean solar day.

With respect to the distant stars, however, for example the stars of the zodiac constellations, Moon's angular motion, or sidereal angular motion, is slightly greater. During the course of a solar year it amounts to 360° for each elapsed lunation, plus 360° to account for the motion of Sun through the entire zodiac:

MeanLunationsInYear = SolarYearLength / MeanSynodicMonth

MeanSiderealMotionInYear = ( MeanLunationsInYear × 360° ) + 360°

MeanSiderealAngularMotionPerLunation = MeanSiderealMotionInYear / MeanLunationsInYear

MeanSiderealAngularMotionPerDay = MeanSiderealAngularMotionPerLunation / MeanSynodicMonth

For example, using the Lunation zero MSM and dividing it into the present era mean northward equinoctial year length of 365 days 5 hours and 49 minutes there were about 12.36827268 lunations per year, so the total annual lunar mean sidereal angular motion was about 4812.578°, and dividing that by the MSM yields about 389.107°, indicating that on average Moon moved slightly more than 29.1° in excess of a full 360° orbit per lunation. This explains the long known rule-of-thumb which says that the zodiac sign that is the background for the first visible new lunar crescent at sunset is the same as sign that is the background at the end of the month when the old lunar cresent is last seen in the morning before sunrise, as each zodiac sign spans about 360° / 12 signs of the zodiac = 30° and in the interim Sun advances one zodiac sign eastward. Dividing 389.107° by the MSM we find that the mean sidereal angular motion was about 13.1764° or 13° 10' 35" per day.

Calculating the length of the mean sidereal month is simply a matter of dividing the daily sidereal angular motion into 360°:

MeanSiderealMonth = 360° / MeanSiderealAngularMotionPerDay

Continuing our example, we have 360° / 13.1764° = about 27.3216 days or 27 days 7 hours 43 minutes and nearly 5 seconds.

If one wishes to calculate the mean sidereal month as directly as possible without going through all of the above steps, then the simplified expression is:

MeanSiderealMonth = ( MeanSynodicMonth × SolarYearLength ) / ( MeanSynodicMonth + SolarYearLength )

The Rate of Change of the Mean Synodic Month

Taking the derivative of each LunationToMSM polynomial gives quadratic (2nd order) polynomials for the rate of change of the mean synodic month (MSM), in terms of a miniscule fraction of a second per lunation. Multiplying by 1000000 gives the rate of change in microseconds per lunation:

LunationToMSMrate_Atomic( L ) = 1000000 ( 3.728585E-16 L2 – 4.043092E-11 L + 1.736907E-05 )

LunationToMSMrate_Solar( L ) = 1000000 ( 3.730276E-16 L2 – 4.043358E-11 L – 2.5120395E-05 )

where L is the lunation number relative to J2000.

Although the curves are parabolic over the evaluated range of 19000 years, the trends actually appear gently sinusoidal, with a period more than double as many years. The obvious correlation to check is the Earth axial tilt (obliquity) cycle, which has a period of about 41000 years, shown as a secondary y-axis in the chart below:

Is There a Long-Term Effect of the Mean Earth Orbital Eccentricity?

It has long been held by authors such as Jean Chapront and D.G. Izotov that Earth orbital eccentricity has a long-term secular effect on the length of the lunar cycle, for which they included terms in their analytical expressions for lunar longitude.

This concept originated with Pierre-Simon Laplace, who "proved" it in his early 19th century Mécanique Céleste treatise. His work, however, was based on comparatively crude contemporary observational data, mean lunar position and mean Earth orbital eccentricity expressions of limited accuracy, he and his contemporaries dismissed ancient Babylonia records of higher Earth orbital eccentricity as biased by excess observations near syzygy, he considered only a few ancient eclipses, and he didn't consider variations in the Earth rotation rate (unrecognized at the time). Furthermore his proof concluded with the "explanation" that Moon is gradually approaching Earth, rather than the opposite actual reality.

The following chart shows the direct correlation between the rate of change of the mean synodic month (atomic time) and the mean Earth orbital eccentricity as calculated according to the 18-term astronomical algorithm of Pierre Bretagnon (1984), as given by Jean Meeus in chapter 33, "Long-term variations of the orbit of the Earth" in More Mathematical Astronomy Morsels, pages 201-205, published in 2002 by Willmann-Bell, Richmond VA:

The long-term variations of mean Earth orbital eccentricity has maxima at intervals of around 100000 years, much longer than the 41000-year period that was suggested for the MSM rate of change chart in the section above. For most of the lunations plotted prior to January 2000 AD there appears to be a correlation, but from the present era until 12000 AD two widely separated eccentricity values are associated with the same MSM rate of change. If declining mean orbital eccentricity is truly the explanation for the progressive reduction in MSM rate of change in the past then that trend ought to continue into the future as the mean eccentricity will continue to decline until Earth's orbit becomes nearly circular (29000-30000 AD, according to the Bretagnon algorithm).

Therefore, although mean Earth orbital eccentricity has a well-known effect on the short-term periodic variation of the length of the lunar cycle (find the word "eccentricity" on this page), it doesn't plausibly account for its long-term secular variation.

The Role of Earth Axial Tilt (Obliquity)

Finally, as evidence supporting a relationship to Earth axial tilt, the following chart shows the direct correlation between the rate of change of the mean synodic month (atomic time) and the obliquity angles that are actually used within SOLEX 9.1β:

Indeed there was an exponentially greater rate of change in the lunar revolution time in the remote past when Earth axial tilt was near its periodic maximum. With substantially more angular momentum transferred from Earth to Moon during that era, tidal slowing of the Earth rotation rate must have been much greater than it is today, even more so when one considers the probably much greater polar ice cap masses and lower global sea levels that must also have existed at that time.

The "hook" at the bottom left of the chart does not have the same significance as the broader hook that is seen on the Earth orbital eccentricity chart, because this one begins much further into the future and could be due to inaccuracy of the polynomial used to project obliquity into the future, or it could be due to inherent hysteresis in the Earth-Moon system, possibly due to a reactive change in the lunar orbital inclination. The important point is that the turning back of the curve does coincide with the reversal of the Earth axial tilt cycle.

According to SOLEX, near the present era the lunar orbital inclination varies annually from 5° to 5.3° relative to the ecliptic plane of date. Also according to SOLEX, the mean lunar orbital inclination relative to the Earth equatorial plane of the date is always equal to the mean Earth axial tilt, periodically varying by an amount equal to plus or minus the lunar orbital inclination relative to the ecliptic plane of the date, with one period equal to the time it takes for the lunar orbital nodes to regress westward 360° (approximately 18.6 years in the present era). In other words, during each lunar orbital node westward regression cycle Earth's axis passes through all possible orientations relative to the equatorial lunar orbital inclination, therefore the mean equatorial lunar orbital inclination per lunar orbital node westward regression cycle must nearly equal the Earth axial tilt. More precisely, because of the small annual variation of the ecliptic lunar orbital inclination the mean equatorial lunar orbital inclination more exactly matches the mean Earth axial tilt when averaged over an integral number of years, such as 10 lunar orbital node westward regression cycles (≈186 years).

Lunar secular accleration expressions that employ terms for Earth orbital eccentricity instead of Earth axial tilt most likely appear to work because for several millennia around the present era there happens to be a non-causal but almost perfectly linear correlation between them, as shown below:

In conclusion, it is an oversimplification to describe the tidal acceleration of Moon as a constant rate "secular" change, because it more likely varies as a periodic function of the Earth axial tilt cycle. Earth's oblateness is also a very important factor, both in the short-term and in the long-term, but the most accurate observational data in that regard is based on satellite laser ranging, only available since the 1980s, which is too limited to interpret with regard to long-term trends.

Future Research

SOLEX can output orbital osculating elements for the planets and our Moon, for any range of dates. It would be nice to experiment with this capability to evaluate long-range trends in:

It would be nice to extend the range of lunations further into the remote past and distant future to show at least one full Earth axial tilt cycle (41000 years), but that would require improved expressions for:


Updated Jun 22, 2011 (Symmetry454) = Jun 24, 2011 (Symmetry010) = Jun 27, 2011 (Gregorian) = J2000 Mean Lunation = 141.87