For lack of other means, I had to resort to empirical evidence. I am not a statistician to draw conclusions from the data, but as far as astropy is to be trusted (and together with Greg Miller answer), we can comfortably say that the average elongation between the Sun and the Moon is indeed less than 90 and about 89.9773. Not less than 1 minute 20 seconds of arc below 90.
using astropy, I was probing the elongation every 100 minutes for a period of 20 Metonic cycles from 1971 (about 380 years).(*). Hopefully I have no errors in the procedure.
from datetime import datetime, timedelta
from astropy.coordinates import get_sun, get_moon
from astropy.time import Time
SYNODIC_MONTH = 29.530588
sign, d, m, s = x.signed_dms
return sign*(((s+60*m) / 3600.0) + d)
def get_mean_elongation(start_time, duration_years, resulation_minutes=60):
full_synodic_months_in_period = int((duration_years * 365.25) / SYNODIC_MONTH)
total_seconds_in_period = full_synodic_months_in_period * SYNODIC_MONTH * 24 * 3600
start_timestamp = (start_time - datetime(1970, 1, 1)).total_seconds()
measurement_dates = numpy.arange(start_timestamp,
int(start_timestamp + total_seconds_in_period),
times = Time(measurement_dates, format='unix')
sun_locations = get_sun(times)
sun_longitudes = sun_locations.geocentrictrueecliptic.lon
moon_locations = get_moon(times)
moon_longitudes = moon_locations.geocentrictrueecliptic.lon
raw_elongations = sun_longitudes - moon_longitudes
elongations = [abs((arc_angle_to_decimal_angle(x) + 180) % 360 - 180) for x in raw_elongations]
initial_time = datetime(year=1971, month=1, day=1)
# randomly using Metonic cycle as our time unit. duration for each run.
start_times = [initial_time + i*timedelta(days=SYNODIC_MONTH*235.0) for i in range(20)]
for start_time in start_times:
mean_elongation_in_cycle = get_mean_elongation(start_time, 19.01, 100)
Overall the mean elongation was only 89.977126428844, with std dev of 0.03 degrees between the cycles. However, every cycle is an average itself of 235 synodic months. the averages of the 20 cycles:
It is interesting to see that the average of each cycle never suppresses 90.02 degrees, but can go down as far 89.93 (three times out of 20). So maybe there is indeed some asymmetry here.
Edit: I later ran the same, again with 20 cycles, but this time of 4 Metonic cycles (~76 years) [though with somewhat decreased resolution but this seems to be insignificant]. In every single one, the mean elongation turned up to be below 90 degrees. with an average of 89.977365350276. here are the numbers:
For fun, I ran it again, this time using start-date 4000 years ahead. the results were about the same. (again, all the cycles below 90, avg: 89.97755)
Now, we are left with the question of why this happens. very far from being sure, I still prefer to stick to my original hunch: that this deviation from 90 degrees is due to the precessions, with their irregular speeds, which are ultimately related to the Sun- so it will somehow, with a complex mechanism, cause it to be so. For I can't find any other reason. (1) When one looks at the inequalities of the Moon and their principal terms, which I believe should be accurate within 1 minute of arc, he can see (unless I'm mistaken) that there is symmetry there (trying to have $D = 180 + D$); though there is asymmetry between the sides of the Full and New Moons (Parallactic inequality) this is still a mirror between them. [the reader should not however that all values there are mean values] (2) the main precession cycles are relatively fast. The most important for our sake here, the apsidal precession, is only 8.85 years; hence I expect that if the precessions and irregularities were really independent, their effect would cancel quite fast.
(*) This is the first time I was using astropy, and had some odd warnings when I run the programs. Besides, for some reason, when I was testing the function I noticed the latitude components of the Sun do not return 0, but can be several (up to 6 seconds of arc). Not sure why is that, and what does it mean about the reliability and accuracy of the longitude component that was used in this software.