[PATCH] libc/time/gmtime_r.c (gmtime_r): Fixed bug in calculations for dates after year 2069 or before year 1901

[Freddie Chopin <freddie_chopin at op dot pl>]

Hello!

It turns out that the implementation of gmtime_r() with my modifications 
from September has a bug... The problem with this bug is that it shows 
only for dates after 2069 or before 1901, so to actually experience this 
bug time_t must be 64-bits long (or 32-bit unsigned, but I don't think 
such combination exists for any platform). The root cause of the problem 
is that the assumptions that there are 97 leap years in 400-year periods 
and 24 leap years in 100-year periods are valid only if you consider 
"full cycles". So there are 97 leap years between 2000 and 2400, or 24 
leap years between 2000 and 2100, but _NOT_ between 1970 and 2070... The 
fix is to move the epoch from 1st January 1970 to 1st March of 2000 
(beginning of 400-year cycle, right after additional day of leap year - 
29th February 2000), perform the calculations and adjust back at the end.

Attached patch implements this fix. The solution was tested by comparing 
the produced value with gmtime_r() from GLIBC for every hour in 
1602-year period (years 1399 - 2801, inclusive).

The idea for solution comes from implementation of the same 
functionality in musl library, which used almost the same idea that I 
used initially (that was a coincidence).

Unfortunatelly on ARM Cortex-M3 platform these modifications make the 
function run twice as long as previously (2M iterations take 29s vs 15s 
on a 16MHz clock)... I'm not really sure how is that possible (the 
changes don't look that heavy), but that's the reality... As a possible 
optimization, quite a lot of the code can be removed for platforms where 
time_t is 32-bit long, but that's another thing.

BTW - are there any platforms which have time_t as a 64-bit variable 
that use newlib, that would be affected by this problem?

Regards,
FCh

>From c95890ba8d8074d1cb3130669a35f5278e2a088a Sun, 7 Dec 2014 19:11:00 +0100
From: Freddie Chopin <freddie.chopin@gmail.com>
Date: Sun, 7 Dec 2014 19:10:40 +0100
Subject: [PATCH] libc/time/gmtime_r.c (gmtime_r): Fixed bug in calculations for dates after year 2069 or before year 1901. Ideas for solution taken from musl's __secs_to_tm()

diff --git a/newlib/libc/time/gmtime_r.c b/newlib/libc/time/gmtime_r.c
index dddd576..f50c199 100644
--- a/newlib/libc/time/gmtime_r.c
+++ b/newlib/libc/time/gmtime_r.c
@@ -6,6 +6,9 @@
  * - Fixed bug in mday computations - 08/12/04, Alex Mogilnikov <alx@intellectronika.ru>
  * - Fixed bug in __tzcalc_limits - 08/12/04, Alex Mogilnikov <alx@intellectronika.ru>
  * - Move code from _mktm_r() to gmtime_r() - 05/09/14, Freddie Chopin <freddie_chopin@op.pl>
+ * - Fixed bug in calculations for dates after year 2069 or before year 1901. Ideas for
+ *   solution taken from musl's __secs_to_tm() - 07/12/2014, Freddie Chopin
+ *   <freddie_chopin@op.pl>
  *
  * Converts the calendar time pointed to by tim_p into a broken-down time
  * expressed as local time. Returns a pointer to a structure containing the
@@ -14,19 +17,28 @@
 
 #include "local.h"
 
+/* move epoch from 01.01.1970 to 01.03.2000 - this is the first day of new
+ * 400-year long cycle, right after additional day of leap year. This adjustment
+ * is required only for date calculation, so instead of modifying time_t value
+ * (which would require 64-bit operations to work correctly) it's enough to
+ * adjust the calculated number of days since epoch. */
+#define EPOCH_ADJUSTMENT_DAYS	11017
+/* year to which the adjustment was made */
+#define ADJUSTED_EPOCH_YEAR	2000
+/* 1st March of 2000 is Wednesday */
+#define ADJUSTED_EPOCH_WDAY	3
 /* there are 97 leap years in 400-year periods. ((400 - 97) * 365 + 97 * 366) */
 #define DAYS_PER_400_YEARS	146097L
 /* there are 24 leap years in 100-year periods. ((100 - 24) * 365 + 24 * 366) */
 #define DAYS_PER_100_YEARS	36524L
 /* there is one leap year every 4 years */
 #define DAYS_PER_4_YEARS	(3 * 365 + 366)
-
-static _CONST int days_per_year[4] = {
-  365,  /* 1970 or 1966 */
-  365,  /* 1971 or 1967 */
-  366,  /* 1972 or 1968 */
-  365   /* 1973 or 1969 */
-} ;
+/* number of days in a non-leap year */
+#define DAYS_PER_YEAR		365
+/* number of days in January */
+#define DAYS_IN_JANUARY		31
+/* number of days in non-leap February */
+#define DAYS_IN_FEBRUARY	28
 
 struct tm *
 _DEFUN (gmtime_r, (tim_p, res),
@@ -35,14 +47,14 @@
 {
   long days, rem;
   _CONST time_t lcltime = *tim_p;
-  int year;
-  int years400, years100, years4;
+  int year, month, yearday, weekday;
+  int years400, years100, years4, remainingyears;
   int yearleap;
   _CONST int *ip;
 
-  days = ((long)lcltime) / SECSPERDAY;
+  days = ((long)lcltime) / SECSPERDAY - EPOCH_ADJUSTMENT_DAYS;
   rem = ((long)lcltime) % SECSPERDAY;
-  while (rem < 0)
+  if (rem < 0)
     {
       rem += SECSPERDAY;
       --days;
@@ -55,45 +67,68 @@
   res->tm_sec = (int) (rem % SECSPERMIN);
 
   /* compute day of week */
-  if ((res->tm_wday = ((EPOCH_WDAY + days) % DAYSPERWEEK)) < 0)
-    res->tm_wday += DAYSPERWEEK;
+  if ((weekday = ((ADJUSTED_EPOCH_WDAY + days) % DAYSPERWEEK)) < 0)
+    weekday += DAYSPERWEEK;
+  res->tm_wday = weekday;
 
   /* compute year & day of year */
   years400 = days / DAYS_PER_400_YEARS;
   days -= years400 * DAYS_PER_400_YEARS;
+  /* simplify by making the values positive */
+  if (days < 0)
+    {
+      days += DAYS_PER_400_YEARS;
+      --years400;
+    }
+
   years100 = days / DAYS_PER_100_YEARS;
+  if (years100 == 4) /* required for proper day of year calculation */
+    --years100;
   days -= years100 * DAYS_PER_100_YEARS;
   years4 = days / DAYS_PER_4_YEARS;
   days -= years4 * DAYS_PER_4_YEARS;
+  remainingyears = days / DAYS_PER_YEAR;
+  if (remainingyears == 4) /* required for proper day of year calculation */
+    --remainingyears;
+  days -= remainingyears * DAYS_PER_YEAR;
 
-  year = EPOCH_YEAR + years400 * 400 + years100 * 100 + years4 * 4;
+  year = ADJUSTED_EPOCH_YEAR + years400 * 400 + years100 * 100 + years4 * 4 +
+      remainingyears;
 
-  /* the value left in days is based in 1970 */
-  if (days < 0)
+  /* If remainingyears is zero, it means that the years were completely
+   * "consumed" by modulo calculations by 400, 100 and 4, so the year is:
+   * 1. a multiple of 4, but not a multiple of 100 or 400 - it's a leap year,
+   * 2. a multiple of 4 and 100, but not a multiple of 400 - it's not a leap
+   * year,
+   * 3. a multiple of 4, 100 and 400 - it's a leap year.
+   * If years4 is non-zero, it means that the year is not a multiple of 100 or
+   * 400 (case 1), so it's a leap year. If years100 is zero (and years4 is zero
+   * - due to short-circuiting), it means that the year is a multiple of 400
+   * (case 3), so it's also a leap year. */
+  yearleap = remainingyears == 0 && (years4 != 0 || years100 == 0);
+
+  /* adjust back to 1st January */
+  yearday = days + DAYS_IN_JANUARY + DAYS_IN_FEBRUARY + yearleap;
+  if (yearday >= DAYS_PER_YEAR + yearleap)
     {
-      ip = &days_per_year[3];
-      while (days < 0)
-        {
-          days += *ip--;
-          --year;
-        }
+      yearday -= DAYS_PER_YEAR + yearleap;
+      ++year;
     }
-  else
-    {
-      ip = &days_per_year[0];
-      while (days >= *ip)
-        {
-          days -= *ip++;
-          ++year;
-        }
-    }
-
+  res->tm_yday = yearday;
   res->tm_year = year - YEAR_BASE;
-  res->tm_yday = days;
-  yearleap = isleap(year);
-  ip = __month_lengths[yearleap];
-  for (res->tm_mon = 0; days >= ip[res->tm_mon]; ++res->tm_mon)
-    days -= ip[res->tm_mon];
+
+  /* Because "days" is the number of days since 1st March, the additional leap
+   * day (29th of February) is the last possible day, so it doesn't matter much
+   * whether the year is actually leap or not. */
+  ip = __month_lengths[1];
+  month = 2;
+  while (days >= ip[month])
+    {
+      days -= ip[month];
+      if (++month >= MONSPERYEAR)
+        month = 0;
+    }
+  res->tm_mon = month;
   res->tm_mday = days + 1;
 
   res->tm_isdst = 0;

Attachment: ChangeLog.txt
Description: Text document