From: Juergen Perlinger Date: Tue, 28 May 2019 06:11:59 +0000 (+0200) Subject: [some cleanup of calendar calculations] X-Git-Url: http://git.ipfire.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=255d2e7bf75b40bb8deacccb3f3a2554284db429;p=thirdparty%2Fntp.git [some cleanup of calendar calculations] bk: 5cecd12f6fmazsLK8n1OFSass93pzw --- diff --git a/include/ntp_calendar.h b/include/ntp_calendar.h index 85fffbedc..a13b80306 100644 --- a/include/ntp_calendar.h +++ b/include/ntp_calendar.h @@ -486,6 +486,48 @@ basedate_expand_gpsweek(unsigned short weekno); */ #define GREGORIAN_CYCLE_WEEKS (GREGORIAN_CYCLE_DAYS / 7) -#define is_leapyear(y) (!((y) % 4) && !(!((y) % 100) && (y) % 400)) +/* + * Is a Greogorian calendar year a leap year? The obvious solution is to + * test the expression + * + * (y % 4 == 0) && ((y % 100 != 0) || (y % 400 == 0)) + * + * This needs (in theory) 2 true divisions -- most compilers check the + * (mod 4) condition by doing a bit test. Some compilers have been + * even observed to partially fuse the (mod 100) and (mod 400) test, + * but there is an alternative formula that gives the compiler even + * better chances: + * + * (y % 4 == 0) && ((y % 16 == 0) || (y % 25 != 0)) + * + * The order of checks is chosen so that the shorcut evaluation can fix + * the result as soon as possible. And the compiler has to do only one + * true division here -- the (mod 4) and (mod 16) can be done with + * direct bit tests. *If* the compiler chooses to do so. + * + * The deduction is as follows: rewrite the standard formula as + * (y % 4 == 0) && ((y % 4*25 != 0) || (y % 16*25 == 0)) + * + * then split the congruences: + * (y % 4 == 0) && ((y % 4 != 0 || y % 25 != 0) || (y % 16 == 0 && y % 25 == 0)) + * + * eliminate the 1st inner term, as it is provably false: + * (y % 4 == 0) && (y % 25 != 0 || (y % 16 == 0 && y % 25 == 0)) + * + * Use the distributive laws on the second major group: + * (y % 4 == 0) && ((y % 25 != 0 || y % 16 == 0) && (y % 25 != 0 || y % 25 == 0)) + * + * Eliminate the constant term, reorder, and voila: + */ +static inline int +is_leapyear(int32_t y) { + return !(y % 4) && (!(y % 16) || (y % 25)); +} +/* The (mod 4) test eliminates 3/4 (or 12/16) of all values. + * The (mod 16) test eliminates another 1/16 of all values. + * 3/16 of all values reach the final division. + * Assuming that the true division is the most costly operation, this + * sequence should give most bang for the buck. + */ #endif diff --git a/libntp/ntp_calendar.c b/libntp/ntp_calendar.c index addd50ab5..19b23c826 100644 --- a/libntp/ntp_calendar.c +++ b/libntp/ntp_calendar.c @@ -1136,6 +1136,13 @@ ntpcal_time_to_date( */ #if !defined(HAVE_INT64) +/* multiplication helper. Seconds in days and weeks are multiples of 128, + * and without that factor fit well into 16 bit. So a multiplication + * of 32bit by 16bit and some shifting can be used on pure 32bit machines + * with compilers that do not support 64bit integers. + * + * Calculate ( hi * mul * 128 ) + lo + */ static vint64 _dwjoin( uint16_t mul, @@ -1144,14 +1151,13 @@ _dwjoin( ) { vint64 res; - uint32_t p1, p2; - int sf; + uint32_t p1, p2, sf; - p1 = (uint32_t)hi; - sf = (hi < 0); - p1 = (p1 + sf) ^ sf; /* absolute value if 'hi' */ + /* get sign flag and absolute value of 'hi' in p1 */ + sf = (uint32_t)-(hi < 0); + p1 = ((uint32_t)hi + sf) ^ sf; - /* assemble major units */ + /* assemble major units: res <- |hi| * mul */ res.D_s.lo = (p1 & 0xFFFF) * mul; res.D_s.hi = 0; p1 = (p1 >> 16) * mul; @@ -1159,17 +1165,18 @@ _dwjoin( p1 = p1 << 16; M_ADD(res.D_s.hi, res.D_s.lo, p2, p1); - /* mul by 128, using shift */ + /* mul by 128, using shift: res <-- res << 7 */ res.D_s.hi = (res.D_s.hi << 7) | (res.D_s.lo >> 25); res.D_s.lo = (res.D_s.lo << 7); - /* fix sign */ - if (sf) - M_NEG(res.D_s.hi, res.D_s.lo); + /* fix up sign: res <-- (res + [sf|sf]) ^ [sf|sf] */ + M_ADD(res.D_s.hi, res.D_s.lo, sf, sf); + res.D_s.lo ^= sf; + res.D_s.hi ^= sf; - /* properly add seconds */ + /* properly add seconds: res <-- res + [sx(lo)|lo] */ + p2 = (uint32_t)-(lo < 0); p1 = (uint32_t)lo; - p2 = UINT32_C(0) - (lo < 0); M_ADD(res.D_s.hi, res.D_s.lo, p2, p1); return res; } diff --git a/tests/libntp/calendar.c b/tests/libntp/calendar.c index 97778959f..a54d9d99d 100644 --- a/tests/libntp/calendar.c +++ b/tests/libntp/calendar.c @@ -33,6 +33,7 @@ void test_SplitYearDays2(void); void test_RataDie1(void); void test_LeapYears1(void); void test_LeapYears2(void); +void test_LeapYears3(void); void test_RoundTripDate(void); void test_RoundTripYearStart(void); void test_RoundTripMonthStart(void); @@ -442,6 +443,21 @@ test_LeapYears2(void) return; } +/* check the 'is_leapyear()' implementation for 4400 years */ +void +test_LeapYears3(void) +{ + int32_t year; + int l1, l2; + + for (year = -399; year < 4000; ++year) { + l1 = (year % 4 == 0) && ((year % 100 != 0) || (year % 400 == 0)); + l2 = is_leapyear(year); + snprintf(mbuf, sizeof(mbuf), "y=%d", year); + TEST_ASSERT_EQUAL_MESSAGE(l1, l2, mbuf); + } +} + /* Full roundtrip from 1601-01-01 to 2400-12-31 * checks sequence of rata die numbers and validates date output * (since the input is all nominal days of the calendar in that range diff --git a/tests/libntp/run-calendar.c b/tests/libntp/run-calendar.c index 26089200e..31f196952 100644 --- a/tests/libntp/run-calendar.c +++ b/tests/libntp/run-calendar.c @@ -41,6 +41,7 @@ extern void test_SplitYearDays2(void); extern void test_RataDie1(void); extern void test_LeapYears1(void); extern void test_LeapYears2(void); +extern void test_LeapYears3(void); extern void test_RoundTripDate(void); extern void test_RoundTripYearStart(void); extern void test_RoundTripMonthStart(void); @@ -89,20 +90,21 @@ int main(int argc, char *argv[]) RUN_TEST(test_RataDie1, 33); RUN_TEST(test_LeapYears1, 34); RUN_TEST(test_LeapYears2, 35); - RUN_TEST(test_RoundTripDate, 36); - RUN_TEST(test_RoundTripYearStart, 37); - RUN_TEST(test_RoundTripMonthStart, 38); - RUN_TEST(test_RoundTripWeekStart, 39); - RUN_TEST(test_RoundTripDayStart, 40); - RUN_TEST(test_IsoCalYearsToWeeks, 41); - RUN_TEST(test_IsoCalWeeksToYearStart, 42); - RUN_TEST(test_IsoCalWeeksToYearEnd, 43); - RUN_TEST(test_DaySecToDate, 44); - RUN_TEST(test_GpsRollOver, 45); - RUN_TEST(test_GpsRemapFunny, 46); - RUN_TEST(test_GpsNtpFixpoints, 48); - RUN_TEST(test_NtpToNtp, 49); - RUN_TEST(test_NtpToTime, 50); + RUN_TEST(test_LeapYears3, 36); + RUN_TEST(test_RoundTripDate, 37); + RUN_TEST(test_RoundTripYearStart, 38); + RUN_TEST(test_RoundTripMonthStart, 39); + RUN_TEST(test_RoundTripWeekStart, 40); + RUN_TEST(test_RoundTripDayStart, 41); + RUN_TEST(test_IsoCalYearsToWeeks, 42); + RUN_TEST(test_IsoCalWeeksToYearStart, 43); + RUN_TEST(test_IsoCalWeeksToYearEnd, 44); + RUN_TEST(test_DaySecToDate, 45); + RUN_TEST(test_GpsRollOver, 46); + RUN_TEST(test_GpsRemapFunny, 47); + RUN_TEST(test_GpsNtpFixpoints, 49); + RUN_TEST(test_NtpToNtp, 50); + RUN_TEST(test_NtpToTime, 51); return (UnityEnd()); }