Analysis by Dr. Irv Bromberg, University of Toronto, Canada
[Click here to go back to the Symmetry454 / Kalendis home page]
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.
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 free computer program written by Professor Aldo Vitagliano of the Department of Chemical Sciences at the University of Naples Federico II, 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 latest version 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://www.solexorb.it/>.
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.
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://individual.utoronto.ca/kalendis/seasons.htm>.
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.
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:
lesser peaksoccur when the lunar conjunction is nearly midway between perigee and apogee.
The following chart, adapted from a similar chart received from Jonathan Jay of Hawaii, is a histogram of the durations (rounded to the nearest minute) of more than 60000 lunar cycles over a 5000-year period, showing that the lesser peaks
are indeed the most common. The most common durations are at 29d 11h 0m and 10h 52m, with lots of nearby similarly high frequency durations, and the secondarily most common durations are at 29d 15h 11m and 15h 16m, again with lots of nearby similarly high frequency durations. Click here or on the chart to open a high resolution PDF version.
The next chart takes an even closer view at the periodic variations by zooming into ±70 lunations relative to the J2000 epoch:
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.530851063829787234042553191489361702127659574468... | 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 (52 yerms) | 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.5305869 | ≈ 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.
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.
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):
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 Phase | Maximum >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:52 | 06:33:40 | 13:21:12 | least variations, about 13+1/3 hours |
First Quadrature (1st Quarter = 90°) | 22:12:09 | 04:13:11 | 17:58:58 | almost 18 hours of variations |
Opposition (Full Moon = 180°) | 19:57:48 | 06:34:19 | 13:23:29 | a tad more variable than conjunctions |
Last Quadrature (3rd Quarter = 270°) | 22:12:54 | 04:14:04 | 17:58:51 | almost 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):
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.
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 high-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 high-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 of Jean Meeus' "About the phases of the Moon" in "Mathematical Astronomy Morsels IV", published by Willmann-Bell, Inc., 2007). 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 early and the mean waning quadrature (fraction = 0.75) tends to be about 10 minutes late.
The following chart shows the variation in Earth-Moon distance over 3 present-era years, with the distances at syzygy events (lunar phase 0° with Earth-Moon-Sun approximately aligned at lunar conjunctions = New Moon, or lunar phase 180° with Moon-Earth-Sun approximately aligned at lunar oppositions = Full Moon) highlighted with colored symbols and dashed curves. Click here or on the chart to open a high-resolution PDF version (which also includes the next 2 charts). The y-axis shows the lunar distance in millions of metres and the x-axis shows the Gregorian date in YYYY-MM-DD format, with major tick-mark interval of 59 days (close to 2 lunar cycles) and minor tick-mark interval of 14.75 days (close to one half lunar cycle). Keep in mind that currently Earth passes perihelion (closest to Sun) on January 3rd or 4th, and passes aphelion (furthest from Sun) on July 4th or 5th.
The following chart again shows the variation in Earth-Moon distance over the same 3 present-era years, but with the distances at the first and last quarter moon events (lunar phase 90° or 270°) highlighted with colored symbols and dashed curves. Click here or on the chart to open a high-resolution PDF version (which also includes the previous and next chart). Although Moon looks half-illuminated at these events, we call it a quarter Moon because they mark the end of the first quarter or the beginning of the last quarter of the lunar cycle.
The following chart combines the previous two charts, making it easy to see that although most of the lunar distance extremes occur near syzygy or quarter events, a few don't (about twice as many far as close extremes). Click here or on the chart to open a high-resolution PDF version (which also includes the previous two charts).
In the above lunar distance charts instances of closest distances seem to occur in groups of 3 or 4 events, and likewise instances of furthest distances seem to occur in groups of 3 or 4 events. This is a hint that such extreme events occurs more commonly. To check whether that is correct, the following chart shows the distribution of Earth-Moon distances during the 21st century, including separate cumulative centile curves for all lunar conjunctions, first quarters, oppositions, and last quarters. Also plotted, as dashed curves, is the first derivative (point-by-point slope) of the conjunction and first quarter curves, showing the relative frequencies of the various lunar distances. Click here or on the chart to open a high-resolution PDF version.
What about the lunar distance when Moon is actually at perigee or apogee? The following chart summarizes the distributions of those distances. The distances for perigee are intentionally plotted "backwards" to show the true relationship with apogee: the closest perigee distances go together with the furthest apogee distances (maximum lunar orbital eccentricity, major axis parallel to the Earth-Sun line), and the furthest perigee distances go together with the closest apogee distances (minimum eccentricity, major axis perpendicular to the Earth-Sun line). Click here or on the chart to open a high-resolution PDF version (which also includes the next two charts).
The next two charts zoom in on perigee and apogee separately, each plotted as an ascending cumulative centile curve along with its first derivative (slope), showing yet again that the extreme distances occur more frequently than intermediate distances. This implies that the lunar orbit alternates relatively rapidly between minimum and maximum eccentricity, spending relatively little time at intermediate eccentricities. Click here or on the charts to open a high-resolution PDF version (which also includes the previous chart).
The perigee distance varies over a much greater range (almost 14 × 106 metres) than does the apogee distance (less than 3 × 106 metres). The higher frequency of the perigee derivative curve at the closest distances indicates that these occur most frequently.
The apogee derivative curve, which resembles an inverted Gaussian distribution curve, indicates that the closest and furthest apogee distances occur most frequently and with almost equal frequency.
I calculated the moments of these events using the new-moon-at-or-after or lunar-phase-at-or-after functions, and then calculated the distances between the centers of Earth and Moon using the lunar-distance function as published in the book Calendrical Calculations, 3rd edition, by Nachum Dershowitz and Edward M. Reingold, Cambridge University Press, 2008.
For more discussion and charts about lunar distances and instantaneous lunar orbital elements, see pages 11-25 (chapters 1-4) of Jean Meeus' Mathematical Astronomy Morsels, published by Willmann-Bell Inc., 1997. For calculating the moments of perigee and apogee, I started with the mean perigee from Jean Meeus' Astronomical Algorithms, 2nd edition, published by Willmann-Bell Inc., 1999, chapter 45 "Position of the Moon", page 343, based on Chapront ELP2000-82B, 1998 (mean apogee is the same but add 180° then modulo 360°), then used a binary search within ±2 days to find the actual minimum perigee or maximum apogee distance to within 1 metre precision. The mean-to-actual differences for the 21st century were within ±25.1° for perigee and ±5.6° for apogee.
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 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 ) = ThisMeanNewMoon – DeltaT( 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.
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 = ToLunation – FromLunation
RETURN MSM × ElapsedMonths – EndMeanConjunction + 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://individual.utoronto.ca/kalendis/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 653KB
This workbook is compatible with Excel 2007 or later for Windows. It’s also compatible with recent versions of Excel for macOS (tested in version 16.62) except that the Shift to Past
and Shift to Future
buttons aren’t visible on the charts — you can still use the Setup
worksheet to change the target year numbers and then update the chart(s). It is incompatible with Excel for macOS 2011 (tested in 14.7.7) — although it opens and displays properly, any attempt to update chart(s) causes Excel to crash. It’s also incompatible with LibreOffice CALC because CALC employs a different chart object model.
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.
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://individual.utoronto.ca/kalendis/seasons.htm>), 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 ) )
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/24 – 44/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).
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 )
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:
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.
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° (in the present era 6798.38 taking days with reference to the northward equinox of date or about 18+2/3 years, or 6793.48 days with reference to the stars). The instantaneous nodal regression rate varies with a principal period of 173.3 days, and the nodes are temporarily stationary whenever the major axis of the lunar elliptical orbit is nearly aligned with the Sun-Earth line. 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.
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 Feb 20, 2022 (Symmetry454) = Feb 18, 2022 (Symmetry010) = Feb 19, 2022 (Gregorian)