原创 Help solve this mind-boggling conundrum

I've just been struggling with a mind-bending conundrum. This is all part of my ongoing Inamorata Prognostication Engine project. I'm currently working on the code to calculate how many days it is until the next full moon, a time that—experience has shown—can be somewhat troublesome with regards to one's Inamorata.

The problem was that my code was producing unexpected results. I was starting to pull my hair out. By yesterday, I was contemplating posting a blog asking for your help. But I was mulling things over in the wee hours of this morning when I thought, "Hang on just a minute. Could it be that...?" It could indeed! The bottom line is that I've tracked down the root of my problem and everything is once again running smoothly in "The House of Max," but for a moment there I thought I was losing my mind.

What? You think you could have sorted this out faster than I? Well, let's see, shall we? The remainder of this blog reflects my musings from yesterday, February 20, 2014, before I'd worked out what I was doing wrong...

Remember that I'm working with an Arduino, so my float variables have only around 6 to 7 decimal digits of precision. The output from my real time clock shows me that today's date (at the time of this writing) is February 20, 2014. Since today is indeed the February 20, I'm happy so far.

I use this date to calculate today's Julian day number and obtain a value of 2456709. I bounce over to the USNO website to confirm that they are in agreement with me and that today does indeed correspond to a Julian day number of 2456709. I'm still smiling.

Next, I go to the MoonPhases.Info website and look up the full moon dates for 2014 as illustrated below.

Purely for the sake of a quick test, I use the date of the most recent full moon—February 14, 2014—as my reference point (any full moon in the past 1,000 years would do, so later I intend to replace this reference with something more appropriate, like a full moon that fell on April 1, or one that occurred on the day that something interesting or relevant happened—do you have any suggestions?).

But we digress... I enter February 14, 2014 as the date of the reference full moon into my code, which calculates that this corresponds to a Julian day number of 2456703. Since this was only six days ago, and since 2456709 – 2456703 = 6, I think we can safely say that we're still on track.

But this is where my smile turns upside down into a frown, because my Arduino now tells me that—after performing some heroic calculations—it has determined that there are 21 days until the next full moon and that the date of that full moon will be March 13, 2014. I would say "Close, but no cigar," except that this is not even particularly close. The actual date, as shown in the image above, should be March 16, 2014, which is 24 days in our future. "Oh dear," I said to myself (or words to that effect).

In a moment I'll show you the code I'm using to calculate how many days there are to the next full moon. From there, we can calculate the actual date of the next full moon. But first, let me walk you through the reasoning process I went through to determine what algorithm to use (remember that, as for most things in life, I'm making all of this up as I go along).

Purely for the sake of a starting-point example, let's assume that our universe has only been in existence for a little over 100 days. Let's also assume that a full moon occurs every 10 days on the dot—that is, there was full moon on days 10, 20, 30, 40, etc. Let's pick one of these full moons as our reference full moon—say the one that occurred on day 100.

Now, let's assume that we are currently on day number 134 in our hypothetical universe. So, the number of days (the "difference" or 'd') between today and our reference full moon is 134—100 = 34. We know that the period ('p') of our moon's orbit is 10 days. What we want to do is to divide the number of days 'd' by the period 'p' and keep the remainder. In computer terms, when performing integer math, we have two types of divide operation available to us:

// Standard divide operation; remainder is lost
d / p = 34 / 10 = 3

// Modulo divide operation; returns the remainder
d % p = 34 % 10 = 4

In our terms, this remainder of 4 equates to 4 days. Thus, since the period of our moon's orbit is 10 days, the number of days to the next full moon is p—4 = 6. Furthermore, since we know that today's number is 134, we therefore know that the next full moon will occur on day 134 + 6 = 140. Tra la!

Of course, our test as described above was based on integer math operations. If we had been dealing with real numbers, then d / p = 34 / 10 = 3.4. If we now throw away the integer part of this result to leave a remainder of 0.4, we have to multiply this by the period to obtain the number of days; that is, 0.4 x p = 0.4 x 10 = 4 days.

Now, let me show you the code I'm using to generate the number of days to the next full moon and, from there, the actual date of the next full moon. Then we'll see how long it takes you to spot where I went wrong.

I pass my Julian day numbers (both are of type long) as parameters into this function—"jNow" and "jRef" are the Julian day numbers for today and for the day of my reference full moon, respectively. You'll just have to take my word that these values are 2456709 and 2456703 as discussed above.

The first thing I do is to declare a variable called 'd' of type long and assign it the result of today's Julian day number minus the Julian day number of my reference full moon. In the case of the example values I'm using, the result is 2456709—2456703 = 6.

When I was a kid, I was always told that the Moon took 28 days to orbit the Earth. Of course, I know that this is only an approximation, so I performed a Google Search "How long does it take the moon to orbit the earth?" and received the following response: "The sidereal month is the time it takes to make one complete orbit of the earth with respect to the fixed stars—it is about 27.32 days."

The "about" part of this answer raised a small "red flag," so I Googled further and obtained a more accurate value for a sidereal month of 27.321661 days. This is why I declare a float 'p' (for "period") and assign it to a label I defined earlier. Once again, you'll just have to trust me that I did define "siderealMonth" as having a value of 27.321661.

Now, I want to get the remainder from dividing the number of days 'd' by the period of the lunar month 'p,' but the modulus ("%") operator doesn't work with floats. The solution is to perform the operation "((d/p)—floor(d/p))." The idea here is that the "floor()" function throws away the remainder from a floating point value. The "(d/p)" operation gives me the full result from the division in the form of both the integer and fractional parts, and from this I subtract the integer part provided by "floor(d/p)," which therefore leaves me with the fractional part; i.e., the remainder.

Since I come from the days when a division was a computationally "expensive" operation—and since "((d/p)—floor(d/p))" requires me to perform the division "(d/p)" two times—I decided to pre-calculate this using the "float dDp = (float) d / p;" statement. Obviously I know that this saving is miniscule in the scheme of things, but "old habits die hard," as they say.

The rest of the code is pretty much self-explanatory; "(((dDp – floor(dDp)) * p)" is the remainder multiplied by the period, which gives me the number of days since the previous full moon, so subtracting this from 'p' gives me the number of days to the next full moon. Of course, this result is currently stored in the "tmp" variable, which is a float, but I want to cast this (and return it) as an int. The problem is that this casting operation will simply truncate any fractional value, but I would prefer to perform a round operation instead. This explains the reason for my adding 0.5 to the result.

As an aside, if you are unsure what I mean, by the previous point, then think of it this way... If I have a float result of 6.x, for example, then casting this to an int will truncate it to 6. This would be OK if my original number was 6.4, but not OK if my original number was 6.6. In this latter case, I would prefer my final value to be 7. The way to achieve this is to add the 0.5. In the case of 6.4 + 0.5 = 6.9, the cast to an int will still truncate this to 6; by comparison, in the case of 6.6 + 0.5 = 7.1, the cast to an int will truncate this to 7, which is what I want.

So, returning to the problem at hand, when I run the above function, supplying it with values of 2456709 for today's Julian day number, 2456703 for the Julian day number of my reference full moon, and 27.321661 as the value of the sidereal month, my function returns a value of 21 days to the next full moon. I then use this value to calculate that February 20 + 21 days = March 13. But, as we said, the actual date of the next full moon should be March 16, 2014, which is 24 days in our future.

First of all I worked the problem through by hand, but I obtained exactly the same result as my function. On the one hand this made me happy (I hadn't made a stupid coding mistake—you can end up with all sorts of unexpected results if you fail to cast things the right way); on the other hand it made me sad (I still didn't have a clue what was going wrong).

By this time I was "clutching at straws," because I started to wonder if there was any chance the creators of the MoonPhases.Info website had made a mistake, so I bounced over to TimeAndDate.com to check out their version of the 2014 calendar. Bummer; they also have full moons occurring on February 14 and March 16 as shown below. My brain hurts.

So, that's where I left things when I retired to my bed yesterday evening. In the wee hours of this morning, I was lying in bed running through things in my mind, dissecting my algorithm operation by operation, trying different test values to see what would happen, and generally getting nowhere fast. And then came the breakthrough; it gradually dawned on me that there was only one place where the error could be originating, and that place was... [your answer goes here].

So, what do you think? Please post your thoughts in the comments below.