diff --git a/cbp/ModMMDB.cbp b/cbp/ModMMDB.cbp
new file mode 100644
index 00000000..f8a57614
--- /dev/null
+++ b/cbp/ModMMDB.cbp
@@ -0,0 +1,407 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/cbp/Module.cbp b/cbp/Module.cbp
index 745f74ec..c3865803 100644
--- a/cbp/Module.cbp
+++ b/cbp/Module.cbp
@@ -362,6 +362,9 @@
+
+
+
@@ -372,6 +375,7 @@
+
@@ -387,7 +391,6 @@
-
@@ -438,14 +441,14 @@
+
+
-
-
diff --git a/cbp/Sandbox.cbp b/cbp/Sandbox.cbp
index 0e76327d..7567ee8d 100644
--- a/cbp/Sandbox.cbp
+++ b/cbp/Sandbox.cbp
@@ -149,6 +149,7 @@
+
@@ -160,6 +161,7 @@
+
diff --git a/cbp/default.workspace b/cbp/default.workspace
index 09de101a..0d207ea8 100644
--- a/cbp/default.workspace
+++ b/cbp/default.workspace
@@ -8,6 +8,7 @@
+
diff --git a/config/mingw32/maxminddb_config.h b/config/mingw32/maxminddb_config.h
new file mode 100644
index 00000000..61bd5f3e
--- /dev/null
+++ b/config/mingw32/maxminddb_config.h
@@ -0,0 +1,193 @@
+/* include/maxminddb_config.h. Generated from maxminddb_config.h.in by configure. */
+#ifndef MAXMINDDB_CONFIG_H
+#define MAXMINDDB_CONFIG_H
+
+/* Define to 1 if you have the header file. */
+#define HAVE_ARPA_INET_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_ASSERT_H 1
+
+/* Define to 1 if the system has the type `boolean'. */
+/* #undef HAVE_BOOLEAN */
+
+/* Define to 1 if you have the header file. */
+#define HAVE_DLFCN_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_FCNTL_H 1
+
+/* Define to 1 if you have the `getpagesize' function. */
+#define HAVE_GETPAGESIZE 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_INTTYPES_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_LIBGEN_H 1
+
+/* Define to 1 if your system has a GNU libc compatible `malloc' function, and
+ to 0 otherwise. */
+#define HAVE_MALLOC 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_MATH_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_MEMORY_H 1
+
+/* Define to 1 if you have a working `mmap' system call. */
+#define HAVE_MMAP 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_NETDB_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_NETINET_IN_H 1
+
+/* Has an open_memstream() function */
+#define HAVE_OPEN_MEMSTREAM 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_STDARG_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_STDBOOL_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_STDINT_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_STDIO_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_STDLIB_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_STRINGS_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_STRING_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_SYS_MMAN_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_SYS_PARAM_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_SYS_SOCKET_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_SYS_STAT_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_SYS_TIME_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_SYS_TYPES_H 1
+
+/* Define to 1 if you have the header file. */
+#define HAVE_UNISTD_H 1
+
+/* Define to the sub-directory in which libtool stores uninstalled libraries.
+ */
+#define LT_OBJDIR ".libs/"
+
+/* Missing the unsigned __int128 type */
+#define MMDB_UINT128_IS_BYTE_ARRAY 1
+
+/* int128 types are available with __attribute__((mode(TI))) */
+/* #undef MMDB_UINT128_USING_MODE */
+
+/* Name of package */
+#define PACKAGE "libmaxminddb"
+
+/* Define to the address where bug reports for this package should be sent. */
+#define PACKAGE_BUGREPORT "support@maxmind.com"
+
+/* Define to the full name of this package. */
+#define PACKAGE_NAME "libmaxminddb"
+
+/* Define to the full name and version of this package. */
+#define PACKAGE_STRING "libmaxminddb 1.1.4"
+
+/* Define to the one symbol short name of this package. */
+#define PACKAGE_TARNAME "libmaxminddb"
+
+/* Define to the home page for this package. */
+#define PACKAGE_URL ""
+
+/* Define to the version of this package. */
+#define PACKAGE_VERSION "1.1.4"
+
+/* Define to 1 if you have the ANSI C header files. */
+#define STDC_HEADERS 1
+
+/* Version number of package */
+#define VERSION "1.1.4"
+
+/* Define for Solaris 2.5.1 so the uint32_t typedef from ,
+ , or is not used. If the typedef were allowed, the
+ #define below would cause a syntax error. */
+/* #undef _UINT32_T */
+
+/* Define for Solaris 2.5.1 so the uint64_t typedef from ,
+ , or is not used. If the typedef were allowed, the
+ #define below would cause a syntax error. */
+/* #undef _UINT64_T */
+
+/* Define for Solaris 2.5.1 so the uint8_t typedef from ,
+ , or is not used. If the typedef were allowed, the
+ #define below would cause a syntax error. */
+/* #undef _UINT8_T */
+
+/* Define to rpl_malloc if the replacement function should be used. */
+/* #undef malloc */
+
+/* Define to `long int' if does not define. */
+/* #undef off_t */
+
+/* Define to the equivalent of the C99 'restrict' keyword, or to
+ nothing if this is not supported. Do not define if restrict is
+ supported directly. */
+#define restrict __restrict
+/* Work around a bug in Sun C++: it does not support _Restrict or
+ __restrict__, even though the corresponding Sun C compiler ends up with
+ "#define restrict _Restrict" or "#define restrict __restrict__" in the
+ previous line. Perhaps some future version of Sun C++ will work with
+ restrict; if so, hopefully it defines __RESTRICT like Sun C does. */
+#if defined __SUNPRO_CC && !defined __RESTRICT
+# define _Restrict
+# define __restrict__
+#endif
+
+/* Define to `unsigned int' if does not define. */
+/* #undef size_t */
+
+/* Define to `int' if does not define. */
+/* #undef ssize_t */
+
+/* Define to the type of an unsigned integer type of width exactly 32 bits if
+ such a type exists and the standard includes do not define it. */
+/* #undef uint32_t */
+
+/* Define to the type of an unsigned integer type of width exactly 64 bits if
+ such a type exists and the standard includes do not define it. */
+/* #undef uint64_t */
+
+/* Define to the type of an unsigned integer type of width exactly 8 bits if
+ such a type exists and the standard includes do not define it. */
+/* #undef uint8_t */
+
+#ifndef MMDB_UINT128_USING_MODE
+/* Define as 1 if we we use unsigned int __atribute__ ((__mode__(TI))) for uint128 values */
+#define MMDB_UINT128_USING_MODE 0
+#endif
+
+#ifndef MMDB_UINT128_IS_BYTE_ARRAY
+/* Define as 1 if we don't have an unsigned __int128 type */
+#define MMDB_UINT128_IS_BYTE_ARRAY 0
+#endif
+
+#endif /* MAXMINDDB_CONFIG_H */
diff --git a/external/Common/aes256.cpp b/external/Common/aes256.cpp
new file mode 100644
index 00000000..e6d6582e
--- /dev/null
+++ b/external/Common/aes256.cpp
@@ -0,0 +1,407 @@
+/*
+* Byte-oriented AES-256 implementation.
+* All lookup tables replaced with 'on the fly' calculations.
+*
+* Copyright (c) 2007-2011 Ilya O. Levin, http://www.literatecode.com
+* Other contributors: Hal Finney
+*
+* Permission to use, copy, modify, and distribute this software for any
+* purpose with or without fee is hereby granted, provided that the above
+* copyright notice and this permission notice appear in all copies.
+*
+* THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
+* WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
+* MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
+* ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
+* WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
+* ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
+* OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
+*/
+#include "aes256.h"
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+#define FD(x) (((x) >> 1) ^ (((x) & 1) ? 0x8d : 0))
+
+#define BACK_TO_TABLES
+
+static uint8_t rj_xtime(uint8_t);
+static void aes_subBytes(uint8_t *);
+static void aes_subBytes_inv(uint8_t *);
+static void aes_addRoundKey(uint8_t *, uint8_t *);
+static void aes_addRoundKey_cpy(uint8_t *, uint8_t *, uint8_t *);
+static void aes_shiftRows(uint8_t *);
+static void aes_shiftRows_inv(uint8_t *);
+static void aes_mixColumns(uint8_t *);
+static void aes_mixColumns_inv(uint8_t *);
+static void aes_expandEncKey(uint8_t *, uint8_t *);
+static void aes_expandDecKey(uint8_t *, uint8_t *);
+#ifndef BACK_TO_TABLES
+static uint8_t gf_alog(uint8_t);
+static uint8_t gf_log(uint8_t);
+static uint8_t gf_mulinv(uint8_t);
+static uint8_t rj_sbox(uint8_t);
+static uint8_t rj_sbox_inv(uint8_t);
+#endif
+
+#ifdef BACK_TO_TABLES
+
+static const uint8_t sbox[256] =
+{
+ 0x63, 0x7c, 0x77, 0x7b, 0xf2, 0x6b, 0x6f, 0xc5,
+ 0x30, 0x01, 0x67, 0x2b, 0xfe, 0xd7, 0xab, 0x76,
+ 0xca, 0x82, 0xc9, 0x7d, 0xfa, 0x59, 0x47, 0xf0,
+ 0xad, 0xd4, 0xa2, 0xaf, 0x9c, 0xa4, 0x72, 0xc0,
+ 0xb7, 0xfd, 0x93, 0x26, 0x36, 0x3f, 0xf7, 0xcc,
+ 0x34, 0xa5, 0xe5, 0xf1, 0x71, 0xd8, 0x31, 0x15,
+ 0x04, 0xc7, 0x23, 0xc3, 0x18, 0x96, 0x05, 0x9a,
+ 0x07, 0x12, 0x80, 0xe2, 0xeb, 0x27, 0xb2, 0x75,
+ 0x09, 0x83, 0x2c, 0x1a, 0x1b, 0x6e, 0x5a, 0xa0,
+ 0x52, 0x3b, 0xd6, 0xb3, 0x29, 0xe3, 0x2f, 0x84,
+ 0x53, 0xd1, 0x00, 0xed, 0x20, 0xfc, 0xb1, 0x5b,
+ 0x6a, 0xcb, 0xbe, 0x39, 0x4a, 0x4c, 0x58, 0xcf,
+ 0xd0, 0xef, 0xaa, 0xfb, 0x43, 0x4d, 0x33, 0x85,
+ 0x45, 0xf9, 0x02, 0x7f, 0x50, 0x3c, 0x9f, 0xa8,
+ 0x51, 0xa3, 0x40, 0x8f, 0x92, 0x9d, 0x38, 0xf5,
+ 0xbc, 0xb6, 0xda, 0x21, 0x10, 0xff, 0xf3, 0xd2,
+ 0xcd, 0x0c, 0x13, 0xec, 0x5f, 0x97, 0x44, 0x17,
+ 0xc4, 0xa7, 0x7e, 0x3d, 0x64, 0x5d, 0x19, 0x73,
+ 0x60, 0x81, 0x4f, 0xdc, 0x22, 0x2a, 0x90, 0x88,
+ 0x46, 0xee, 0xb8, 0x14, 0xde, 0x5e, 0x0b, 0xdb,
+ 0xe0, 0x32, 0x3a, 0x0a, 0x49, 0x06, 0x24, 0x5c,
+ 0xc2, 0xd3, 0xac, 0x62, 0x91, 0x95, 0xe4, 0x79,
+ 0xe7, 0xc8, 0x37, 0x6d, 0x8d, 0xd5, 0x4e, 0xa9,
+ 0x6c, 0x56, 0xf4, 0xea, 0x65, 0x7a, 0xae, 0x08,
+ 0xba, 0x78, 0x25, 0x2e, 0x1c, 0xa6, 0xb4, 0xc6,
+ 0xe8, 0xdd, 0x74, 0x1f, 0x4b, 0xbd, 0x8b, 0x8a,
+ 0x70, 0x3e, 0xb5, 0x66, 0x48, 0x03, 0xf6, 0x0e,
+ 0x61, 0x35, 0x57, 0xb9, 0x86, 0xc1, 0x1d, 0x9e,
+ 0xe1, 0xf8, 0x98, 0x11, 0x69, 0xd9, 0x8e, 0x94,
+ 0x9b, 0x1e, 0x87, 0xe9, 0xce, 0x55, 0x28, 0xdf,
+ 0x8c, 0xa1, 0x89, 0x0d, 0xbf, 0xe6, 0x42, 0x68,
+ 0x41, 0x99, 0x2d, 0x0f, 0xb0, 0x54, 0xbb, 0x16
+};
+static const uint8_t sboxinv[256] =
+{
+ 0x52, 0x09, 0x6a, 0xd5, 0x30, 0x36, 0xa5, 0x38,
+ 0xbf, 0x40, 0xa3, 0x9e, 0x81, 0xf3, 0xd7, 0xfb,
+ 0x7c, 0xe3, 0x39, 0x82, 0x9b, 0x2f, 0xff, 0x87,
+ 0x34, 0x8e, 0x43, 0x44, 0xc4, 0xde, 0xe9, 0xcb,
+ 0x54, 0x7b, 0x94, 0x32, 0xa6, 0xc2, 0x23, 0x3d,
+ 0xee, 0x4c, 0x95, 0x0b, 0x42, 0xfa, 0xc3, 0x4e,
+ 0x08, 0x2e, 0xa1, 0x66, 0x28, 0xd9, 0x24, 0xb2,
+ 0x76, 0x5b, 0xa2, 0x49, 0x6d, 0x8b, 0xd1, 0x25,
+ 0x72, 0xf8, 0xf6, 0x64, 0x86, 0x68, 0x98, 0x16,
+ 0xd4, 0xa4, 0x5c, 0xcc, 0x5d, 0x65, 0xb6, 0x92,
+ 0x6c, 0x70, 0x48, 0x50, 0xfd, 0xed, 0xb9, 0xda,
+ 0x5e, 0x15, 0x46, 0x57, 0xa7, 0x8d, 0x9d, 0x84,
+ 0x90, 0xd8, 0xab, 0x00, 0x8c, 0xbc, 0xd3, 0x0a,
+ 0xf7, 0xe4, 0x58, 0x05, 0xb8, 0xb3, 0x45, 0x06,
+ 0xd0, 0x2c, 0x1e, 0x8f, 0xca, 0x3f, 0x0f, 0x02,
+ 0xc1, 0xaf, 0xbd, 0x03, 0x01, 0x13, 0x8a, 0x6b,
+ 0x3a, 0x91, 0x11, 0x41, 0x4f, 0x67, 0xdc, 0xea,
+ 0x97, 0xf2, 0xcf, 0xce, 0xf0, 0xb4, 0xe6, 0x73,
+ 0x96, 0xac, 0x74, 0x22, 0xe7, 0xad, 0x35, 0x85,
+ 0xe2, 0xf9, 0x37, 0xe8, 0x1c, 0x75, 0xdf, 0x6e,
+ 0x47, 0xf1, 0x1a, 0x71, 0x1d, 0x29, 0xc5, 0x89,
+ 0x6f, 0xb7, 0x62, 0x0e, 0xaa, 0x18, 0xbe, 0x1b,
+ 0xfc, 0x56, 0x3e, 0x4b, 0xc6, 0xd2, 0x79, 0x20,
+ 0x9a, 0xdb, 0xc0, 0xfe, 0x78, 0xcd, 0x5a, 0xf4,
+ 0x1f, 0xdd, 0xa8, 0x33, 0x88, 0x07, 0xc7, 0x31,
+ 0xb1, 0x12, 0x10, 0x59, 0x27, 0x80, 0xec, 0x5f,
+ 0x60, 0x51, 0x7f, 0xa9, 0x19, 0xb5, 0x4a, 0x0d,
+ 0x2d, 0xe5, 0x7a, 0x9f, 0x93, 0xc9, 0x9c, 0xef,
+ 0xa0, 0xe0, 0x3b, 0x4d, 0xae, 0x2a, 0xf5, 0xb0,
+ 0xc8, 0xeb, 0xbb, 0x3c, 0x83, 0x53, 0x99, 0x61,
+ 0x17, 0x2b, 0x04, 0x7e, 0xba, 0x77, 0xd6, 0x26,
+ 0xe1, 0x69, 0x14, 0x63, 0x55, 0x21, 0x0c, 0x7d
+};
+
+#define rj_sbox(x) sbox[(x)]
+#define rj_sbox_inv(x) sboxinv[(x)]
+
+#else /* tableless subroutines */
+
+/* -------------------------------------------------------------------------- */
+static uint8_t gf_alog(uint8_t x) // calculate anti-logarithm gen 3
+{
+ uint8_t y = 1, i;
+
+ for (i = 0; i < x; i++) y ^= rj_xtime(y);
+
+ return y;
+} /* gf_alog */
+
+/* -------------------------------------------------------------------------- */
+static uint8_t gf_log(uint8_t x) // calculate logarithm gen 3
+{
+ uint8_t y, i = 0;
+
+ if (x)
+ for (i = 1, y = 1; i > 0; i++ )
+ {
+ y ^= rj_xtime(y);
+ if (y == x) break;
+ }
+
+ return i;
+} /* gf_log */
+
+
+/* -------------------------------------------------------------------------- */
+static uint8_t gf_mulinv(uint8_t x) // calculate multiplicative inverse
+{
+ return (x) ? gf_alog(255 - gf_log(x)) : 0;
+} /* gf_mulinv */
+
+/* -------------------------------------------------------------------------- */
+static uint8_t rj_sbox(uint8_t x)
+{
+ uint8_t y, sb;
+
+ sb = y = gf_mulinv(x);
+ y = (uint8_t)(y << 1) | (y >> 7), sb ^= y;
+ y = (uint8_t)(y << 1) | (y >> 7), sb ^= y;
+ y = (uint8_t)(y << 1) | (y >> 7), sb ^= y;
+ y = (uint8_t)(y << 1) | (y >> 7), sb ^= y;
+
+ return (sb ^ 0x63);
+} /* rj_sbox */
+
+/* -------------------------------------------------------------------------- */
+static uint8_t rj_sbox_inv(uint8_t x)
+{
+ uint8_t y, sb;
+
+ y = x ^ 0x63;
+ sb = y = (uint8_t)(y << 1) | (y >> 7);
+ y = (uint8_t)(y << 2) | (y >> 6);
+ sb ^= y;
+ y = (uint8_t)(y << 3) | (y >> 5);
+ sb ^= y;
+
+ return gf_mulinv(sb);
+} /* rj_sbox_inv */
+
+#endif
+
+/* -------------------------------------------------------------------------- */
+static uint8_t rj_xtime(uint8_t x)
+{
+ uint8_t y = (uint8_t)(x << 1);
+ return (x & 0x80) ? (y ^ 0x1b) : y;
+} /* rj_xtime */
+
+/* -------------------------------------------------------------------------- */
+static void aes_subBytes(uint8_t *buf)
+{
+ register uint8_t i = 16;
+
+ while (i--) buf[i] = rj_sbox(buf[i]);
+} /* aes_subBytes */
+
+/* -------------------------------------------------------------------------- */
+static void aes_subBytes_inv(uint8_t *buf)
+{
+ register uint8_t i = 16;
+
+ while (i--) buf[i] = rj_sbox_inv(buf[i]);
+} /* aes_subBytes_inv */
+
+/* -------------------------------------------------------------------------- */
+static void aes_addRoundKey(uint8_t *buf, uint8_t *key)
+{
+ register uint8_t i = 16;
+
+ while (i--) buf[i] ^= key[i];
+} /* aes_addRoundKey */
+
+/* -------------------------------------------------------------------------- */
+static void aes_addRoundKey_cpy(uint8_t *buf, uint8_t *key, uint8_t *cpk)
+{
+ register uint8_t i = 16;
+
+ while (i--) buf[i] ^= (cpk[i] = key[i]), cpk[16 + i] = key[16 + i];
+} /* aes_addRoundKey_cpy */
+
+
+/* -------------------------------------------------------------------------- */
+static void aes_shiftRows(uint8_t *buf)
+{
+ register uint8_t i, j; /* to make it potentially parallelable :) */
+
+ i = buf[1], buf[1] = buf[5], buf[5] = buf[9], buf[9] = buf[13], buf[13] = i;
+ i = buf[10], buf[10] = buf[2], buf[2] = i;
+ j = buf[3], buf[3] = buf[15], buf[15] = buf[11], buf[11] = buf[7], buf[7] = j;
+ j = buf[14], buf[14] = buf[6], buf[6] = j;
+
+} /* aes_shiftRows */
+
+/* -------------------------------------------------------------------------- */
+static void aes_shiftRows_inv(uint8_t *buf)
+{
+ register uint8_t i, j; /* same as above :) */
+
+ i = buf[1], buf[1] = buf[13], buf[13] = buf[9], buf[9] = buf[5], buf[5] = i;
+ i = buf[2], buf[2] = buf[10], buf[10] = i;
+ j = buf[3], buf[3] = buf[7], buf[7] = buf[11], buf[11] = buf[15], buf[15] = j;
+ j = buf[6], buf[6] = buf[14], buf[14] = j;
+
+} /* aes_shiftRows_inv */
+
+/* -------------------------------------------------------------------------- */
+static void aes_mixColumns(uint8_t *buf)
+{
+ register uint8_t i, a, b, c, d, e;
+
+ for (i = 0; i < 16; i += 4)
+ {
+ a = buf[i];
+ b = buf[i + 1];
+ c = buf[i + 2];
+ d = buf[i + 3];
+ e = a ^ b ^ c ^ d;
+ buf[i] ^= e ^ rj_xtime(a ^ b);
+ buf[i + 1] ^= e ^ rj_xtime(b ^ c);
+ buf[i + 2] ^= e ^ rj_xtime(c ^ d);
+ buf[i + 3] ^= e ^ rj_xtime(d ^ a);
+ }
+} /* aes_mixColumns */
+
+/* -------------------------------------------------------------------------- */
+void aes_mixColumns_inv(uint8_t *buf)
+{
+ register uint8_t i, a, b, c, d, e, x, y, z;
+
+ for (i = 0; i < 16; i += 4)
+ {
+ a = buf[i];
+ b = buf[i + 1];
+ c = buf[i + 2];
+ d = buf[i + 3];
+ e = a ^ b ^ c ^ d;
+ z = rj_xtime(e);
+ x = e ^ rj_xtime(rj_xtime(z ^ a ^ c));
+ y = e ^ rj_xtime(rj_xtime(z ^ b ^ d));
+ buf[i] ^= x ^ rj_xtime(a ^ b);
+ buf[i + 1] ^= y ^ rj_xtime(b ^ c);
+ buf[i + 2] ^= x ^ rj_xtime(c ^ d);
+ buf[i + 3] ^= y ^ rj_xtime(d ^ a);
+ }
+} /* aes_mixColumns_inv */
+
+/* -------------------------------------------------------------------------- */
+static void aes_expandEncKey(uint8_t *k, uint8_t *rc)
+{
+ register uint8_t i;
+
+ k[0] ^= rj_sbox(k[29]) ^ (*rc);
+ k[1] ^= rj_sbox(k[30]);
+ k[2] ^= rj_sbox(k[31]);
+ k[3] ^= rj_sbox(k[28]);
+ *rc = rj_xtime( *rc);
+
+ for(i = 4; i < 16; i += 4) k[i] ^= k[i - 4], k[i + 1] ^= k[i - 3],
+ k[i + 2] ^= k[i - 2], k[i + 3] ^= k[i - 1];
+ k[16] ^= rj_sbox(k[12]);
+ k[17] ^= rj_sbox(k[13]);
+ k[18] ^= rj_sbox(k[14]);
+ k[19] ^= rj_sbox(k[15]);
+
+ for(i = 20; i < 32; i += 4) k[i] ^= k[i - 4], k[i + 1] ^= k[i - 3],
+ k[i + 2] ^= k[i - 2], k[i + 3] ^= k[i - 1];
+
+} /* aes_expandEncKey */
+
+/* -------------------------------------------------------------------------- */
+void aes_expandDecKey(uint8_t *k, uint8_t *rc)
+{
+ uint8_t i;
+
+ for(i = 28; i > 16; i -= 4) k[i + 0] ^= k[i - 4], k[i + 1] ^= k[i - 3],
+ k[i + 2] ^= k[i - 2], k[i + 3] ^= k[i - 1];
+
+ k[16] ^= rj_sbox(k[12]);
+ k[17] ^= rj_sbox(k[13]);
+ k[18] ^= rj_sbox(k[14]);
+ k[19] ^= rj_sbox(k[15]);
+
+ for(i = 12; i > 0; i -= 4) k[i + 0] ^= k[i - 4], k[i + 1] ^= k[i - 3],
+ k[i + 2] ^= k[i - 2], k[i + 3] ^= k[i - 1];
+
+ *rc = FD(*rc);
+ k[0] ^= rj_sbox(k[29]) ^ (*rc);
+ k[1] ^= rj_sbox(k[30]);
+ k[2] ^= rj_sbox(k[31]);
+ k[3] ^= rj_sbox(k[28]);
+} /* aes_expandDecKey */
+
+
+/* -------------------------------------------------------------------------- */
+void aes256_init(aes256_context *ctx, uint8_t *k)
+{
+ uint8_t rcon = 1;
+ register uint8_t i;
+
+ for (i = 0; i < sizeof(ctx->key); i++) ctx->enckey[i] = ctx->deckey[i] = k[i];
+ for (i = 8; --i;) aes_expandEncKey(ctx->deckey, &rcon);
+} /* aes256_init */
+
+/* -------------------------------------------------------------------------- */
+void aes256_done(aes256_context *ctx)
+{
+ register uint8_t i;
+
+ for (i = 0; i < sizeof(ctx->key); i++)
+ ctx->key[i] = ctx->enckey[i] = ctx->deckey[i] = 0;
+} /* aes256_done */
+
+/* -------------------------------------------------------------------------- */
+void aes256_encrypt_ecb(aes256_context *ctx, uint8_t *buf)
+{
+ uint8_t i, rcon;
+
+ aes_addRoundKey_cpy(buf, ctx->enckey, ctx->key);
+ for(i = 1, rcon = 1; i < 14; ++i)
+ {
+ aes_subBytes(buf);
+ aes_shiftRows(buf);
+ aes_mixColumns(buf);
+ if( i & 1 ) aes_addRoundKey( buf, &ctx->key[16]);
+ else aes_expandEncKey(ctx->key, &rcon), aes_addRoundKey(buf, ctx->key);
+ }
+ aes_subBytes(buf);
+ aes_shiftRows(buf);
+ aes_expandEncKey(ctx->key, &rcon);
+ aes_addRoundKey(buf, ctx->key);
+} /* aes256_encrypt */
+
+/* -------------------------------------------------------------------------- */
+void aes256_decrypt_ecb(aes256_context *ctx, uint8_t *buf)
+{
+ uint8_t i, rcon;
+
+ aes_addRoundKey_cpy(buf, ctx->deckey, ctx->key);
+ aes_shiftRows_inv(buf);
+ aes_subBytes_inv(buf);
+
+ for (i = 14, rcon = 0x80; --i;)
+ {
+ if( ( i & 1 ) )
+ {
+ aes_expandDecKey(ctx->key, &rcon);
+ aes_addRoundKey(buf, &ctx->key[16]);
+ }
+ else aes_addRoundKey(buf, ctx->key);
+ aes_mixColumns_inv(buf);
+ aes_shiftRows_inv(buf);
+ aes_subBytes_inv(buf);
+ }
+ aes_addRoundKey( buf, ctx->key);
+} /* aes256_decrypt */
+
+#ifdef __cplusplus
+}
+#endif
diff --git a/external/MaxmindDB/LICENSE b/external/MaxmindDB/LICENSE
new file mode 100644
index 00000000..d6456956
--- /dev/null
+++ b/external/MaxmindDB/LICENSE
@@ -0,0 +1,202 @@
+
+ Apache License
+ Version 2.0, January 2004
+ http://www.apache.org/licenses/
+
+ TERMS AND CONDITIONS FOR USE, REPRODUCTION, AND DISTRIBUTION
+
+ 1. Definitions.
+
+ "License" shall mean the terms and conditions for use, reproduction,
+ and distribution as defined by Sections 1 through 9 of this document.
+
+ "Licensor" shall mean the copyright owner or entity authorized by
+ the copyright owner that is granting the License.
+
+ "Legal Entity" shall mean the union of the acting entity and all
+ other entities that control, are controlled by, or are under common
+ control with that entity. For the purposes of this definition,
+ "control" means (i) the power, direct or indirect, to cause the
+ direction or management of such entity, whether by contract or
+ otherwise, or (ii) ownership of fifty percent (50%) or more of the
+ outstanding shares, or (iii) beneficial ownership of such entity.
+
+ "You" (or "Your") shall mean an individual or Legal Entity
+ exercising permissions granted by this License.
+
+ "Source" form shall mean the preferred form for making modifications,
+ including but not limited to software source code, documentation
+ source, and configuration files.
+
+ "Object" form shall mean any form resulting from mechanical
+ transformation or translation of a Source form, including but
+ not limited to compiled object code, generated documentation,
+ and conversions to other media types.
+
+ "Work" shall mean the work of authorship, whether in Source or
+ Object form, made available under the License, as indicated by a
+ copyright notice that is included in or attached to the work
+ (an example is provided in the Appendix below).
+
+ "Derivative Works" shall mean any work, whether in Source or Object
+ form, that is based on (or derived from) the Work and for which the
+ editorial revisions, annotations, elaborations, or other modifications
+ represent, as a whole, an original work of authorship. For the purposes
+ of this License, Derivative Works shall not include works that remain
+ separable from, or merely link (or bind by name) to the interfaces of,
+ the Work and Derivative Works thereof.
+
+ "Contribution" shall mean any work of authorship, including
+ the original version of the Work and any modifications or additions
+ to that Work or Derivative Works thereof, that is intentionally
+ submitted to Licensor for inclusion in the Work by the copyright owner
+ or by an individual or Legal Entity authorized to submit on behalf of
+ the copyright owner. For the purposes of this definition, "submitted"
+ means any form of electronic, verbal, or written communication sent
+ to the Licensor or its representatives, including but not limited to
+ communication on electronic mailing lists, source code control systems,
+ and issue tracking systems that are managed by, or on behalf of, the
+ Licensor for the purpose of discussing and improving the Work, but
+ excluding communication that is conspicuously marked or otherwise
+ designated in writing by the copyright owner as "Not a Contribution."
+
+ "Contributor" shall mean Licensor and any individual or Legal Entity
+ on behalf of whom a Contribution has been received by Licensor and
+ subsequently incorporated within the Work.
+
+ 2. Grant of Copyright License. Subject to the terms and conditions of
+ this License, each Contributor hereby grants to You a perpetual,
+ worldwide, non-exclusive, no-charge, royalty-free, irrevocable
+ copyright license to reproduce, prepare Derivative Works of,
+ publicly display, publicly perform, sublicense, and distribute the
+ Work and such Derivative Works in Source or Object form.
+
+ 3. Grant of Patent License. Subject to the terms and conditions of
+ this License, each Contributor hereby grants to You a perpetual,
+ worldwide, non-exclusive, no-charge, royalty-free, irrevocable
+ (except as stated in this section) patent license to make, have made,
+ use, offer to sell, sell, import, and otherwise transfer the Work,
+ where such license applies only to those patent claims licensable
+ by such Contributor that are necessarily infringed by their
+ Contribution(s) alone or by combination of their Contribution(s)
+ with the Work to which such Contribution(s) was submitted. If You
+ institute patent litigation against any entity (including a
+ cross-claim or counterclaim in a lawsuit) alleging that the Work
+ or a Contribution incorporated within the Work constitutes direct
+ or contributory patent infringement, then any patent licenses
+ granted to You under this License for that Work shall terminate
+ as of the date such litigation is filed.
+
+ 4. Redistribution. You may reproduce and distribute copies of the
+ Work or Derivative Works thereof in any medium, with or without
+ modifications, and in Source or Object form, provided that You
+ meet the following conditions:
+
+ (a) You must give any other recipients of the Work or
+ Derivative Works a copy of this License; and
+
+ (b) You must cause any modified files to carry prominent notices
+ stating that You changed the files; and
+
+ (c) You must retain, in the Source form of any Derivative Works
+ that You distribute, all copyright, patent, trademark, and
+ attribution notices from the Source form of the Work,
+ excluding those notices that do not pertain to any part of
+ the Derivative Works; and
+
+ (d) If the Work includes a "NOTICE" text file as part of its
+ distribution, then any Derivative Works that You distribute must
+ include a readable copy of the attribution notices contained
+ within such NOTICE file, excluding those notices that do not
+ pertain to any part of the Derivative Works, in at least one
+ of the following places: within a NOTICE text file distributed
+ as part of the Derivative Works; within the Source form or
+ documentation, if provided along with the Derivative Works; or,
+ within a display generated by the Derivative Works, if and
+ wherever such third-party notices normally appear. The contents
+ of the NOTICE file are for informational purposes only and
+ do not modify the License. You may add Your own attribution
+ notices within Derivative Works that You distribute, alongside
+ or as an addendum to the NOTICE text from the Work, provided
+ that such additional attribution notices cannot be construed
+ as modifying the License.
+
+ You may add Your own copyright statement to Your modifications and
+ may provide additional or different license terms and conditions
+ for use, reproduction, or distribution of Your modifications, or
+ for any such Derivative Works as a whole, provided Your use,
+ reproduction, and distribution of the Work otherwise complies with
+ the conditions stated in this License.
+
+ 5. Submission of Contributions. Unless You explicitly state otherwise,
+ any Contribution intentionally submitted for inclusion in the Work
+ by You to the Licensor shall be under the terms and conditions of
+ this License, without any additional terms or conditions.
+ Notwithstanding the above, nothing herein shall supersede or modify
+ the terms of any separate license agreement you may have executed
+ with Licensor regarding such Contributions.
+
+ 6. Trademarks. This License does not grant permission to use the trade
+ names, trademarks, service marks, or product names of the Licensor,
+ except as required for reasonable and customary use in describing the
+ origin of the Work and reproducing the content of the NOTICE file.
+
+ 7. Disclaimer of Warranty. Unless required by applicable law or
+ agreed to in writing, Licensor provides the Work (and each
+ Contributor provides its Contributions) on an "AS IS" BASIS,
+ WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or
+ implied, including, without limitation, any warranties or conditions
+ of TITLE, NON-INFRINGEMENT, MERCHANTABILITY, or FITNESS FOR A
+ PARTICULAR PURPOSE. You are solely responsible for determining the
+ appropriateness of using or redistributing the Work and assume any
+ risks associated with Your exercise of permissions under this License.
+
+ 8. Limitation of Liability. In no event and under no legal theory,
+ whether in tort (including negligence), contract, or otherwise,
+ unless required by applicable law (such as deliberate and grossly
+ negligent acts) or agreed to in writing, shall any Contributor be
+ liable to You for damages, including any direct, indirect, special,
+ incidental, or consequential damages of any character arising as a
+ result of this License or out of the use or inability to use the
+ Work (including but not limited to damages for loss of goodwill,
+ work stoppage, computer failure or malfunction, or any and all
+ other commercial damages or losses), even if such Contributor
+ has been advised of the possibility of such damages.
+
+ 9. Accepting Warranty or Additional Liability. While redistributing
+ the Work or Derivative Works thereof, You may choose to offer,
+ and charge a fee for, acceptance of support, warranty, indemnity,
+ or other liability obligations and/or rights consistent with this
+ License. However, in accepting such obligations, You may act only
+ on Your own behalf and on Your sole responsibility, not on behalf
+ of any other Contributor, and only if You agree to indemnify,
+ defend, and hold each Contributor harmless for any liability
+ incurred by, or claims asserted against, such Contributor by reason
+ of your accepting any such warranty or additional liability.
+
+ END OF TERMS AND CONDITIONS
+
+ APPENDIX: How to apply the Apache License to your work.
+
+ To apply the Apache License to your work, attach the following
+ boilerplate notice, with the fields enclosed by brackets "[]"
+ replaced with your own identifying information. (Don't include
+ the brackets!) The text should be enclosed in the appropriate
+ comment syntax for the file format. We also recommend that a
+ file or class name and description of purpose be included on the
+ same "printed page" as the copyright notice for easier
+ identification within third-party archives.
+
+ Copyright [yyyy] [name of copyright owner]
+
+ Licensed under the Apache License, Version 2.0 (the "License");
+ you may not use this file except in compliance with the License.
+ You may obtain a copy of the License at
+
+ http://www.apache.org/licenses/LICENSE-2.0
+
+ Unless required by applicable law or agreed to in writing, software
+ distributed under the License is distributed on an "AS IS" BASIS,
+ WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ See the License for the specific language governing permissions and
+ limitations under the License.
diff --git a/external/MaxmindDB/maxminddb-compat-util.h b/external/MaxmindDB/maxminddb-compat-util.h
new file mode 100644
index 00000000..e3f0320f
--- /dev/null
+++ b/external/MaxmindDB/maxminddb-compat-util.h
@@ -0,0 +1,167 @@
+#include
+#include
+
+/* *INDENT-OFF* */
+
+/* The memmem, strdup, and strndup functions were all copied from the
+ * FreeBSD source, along with the relevant copyright notice.
+ *
+ * It'd be nicer to simply use the functions available on the system if they
+ * exist, but there doesn't seem to be a good way to detect them without also
+ * defining things like _GNU_SOURCE, which we want to avoid, because then we
+ * end up _accidentally_ using GNU features without noticing, which then
+ * breaks on systems like OSX.
+ *
+ * C is fun! */
+
+/* Applies to memmem implementation */
+/*-
+ * Copyright (c) 2005 Pascal Gloor
+ *
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions
+ * are met:
+ * 1. Redistributions of source code must retain the above copyright
+ * notice, this list of conditions and the following disclaimer.
+ * 2. Redistributions in binary form must reproduce the above copyright
+ * notice, this list of conditions and the following disclaimer in the
+ * documentation and/or other materials provided with the distribution.
+ * 3. The name of the author may not be used to endorse or promote
+ * products derived from this software without specific prior written
+ * permission.
+ *
+ * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
+ * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+ * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+ * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
+ * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+ * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
+ * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+ * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+static void *
+mmdb_memmem(const void *l, size_t l_len, const void *s, size_t s_len)
+{
+ register char *cur, *last;
+ const char *cl = (const char *)l;
+ const char *cs = (const char *)s;
+
+ /* we need something to compare */
+ if (l_len == 0 || s_len == 0)
+ return NULL;
+
+ /* "s" must be smaller or equal to "l" */
+ if (l_len < s_len)
+ return NULL;
+
+ /* special case where s_len == 1 */
+ if (s_len == 1)
+ return memchr(l, (int)*cs, l_len);
+
+ /* the last position where its possible to find "s" in "l" */
+ last = (char *)cl + l_len - s_len;
+
+ for (cur = (char *)cl; cur <= last; cur++)
+ if (cur[0] == cs[0] && memcmp(cur, cs, s_len) == 0)
+ return cur;
+
+ return NULL;
+}
+
+/* Applies to strnlen implementation */
+/*-
+ * Copyright (c) 2009 David Schultz
+ * All rights reserved.
+ *
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions
+ * are met:
+ * 1. Redistributions of source code must retain the above copyright
+ * notice, this list of conditions and the following disclaimer.
+ * 2. Redistributions in binary form must reproduce the above copyright
+ * notice, this list of conditions and the following disclaimer in the
+ * documentation and/or other materials provided with the distribution.
+ *
+ * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
+ * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+ * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+ * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
+ * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+ * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
+ * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+ * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+static size_t
+mmdb_strnlen(const char *s, size_t maxlen)
+{
+ size_t len;
+
+ for (len = 0; len < maxlen; len++, s++) {
+ if (!*s)
+ break;
+ }
+ return (len);
+}
+
+/* Applies to strdup and strndup implementation */
+/*
+ * Copyright (c) 1988, 1993
+ * The Regents of the University of California. All rights reserved.
+ *
+ * Redistribution and use in source and binary forms, with or without
+ * modification, are permitted provided that the following conditions
+ * are met:
+ * 1. Redistributions of source code must retain the above copyright
+ * notice, this list of conditions and the following disclaimer.
+ * 2. Redistributions in binary form must reproduce the above copyright
+ * notice, this list of conditions and the following disclaimer in the
+ * documentation and/or other materials provided with the distribution.
+ * 3. Neither the name of the University nor the names of its contributors
+ * may be used to endorse or promote products derived from this software
+ * without specific prior written permission.
+ *
+ * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
+ * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+ * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+ * ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
+ * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+ * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
+ * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+ * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+static char *
+mmdb_strdup(const char *str)
+{
+ size_t len;
+ char *copy;
+
+ len = strlen(str) + 1;
+ if ((copy = malloc(len)) == NULL)
+ return (NULL);
+ memcpy(copy, str, len);
+ return (copy);
+}
+
+static char *
+mmdb_strndup(const char *str, size_t n)
+{
+ size_t len;
+ char *copy;
+
+ len = mmdb_strnlen(str, n);
+ if ((copy = malloc(len + 1)) == NULL)
+ return (NULL);
+ memcpy(copy, str, len);
+ copy[len] = '\0';
+ return (copy);
+}
+/* *INDENT-ON* */
diff --git a/external/MaxmindDB/maxminddb.c b/external/MaxmindDB/maxminddb.c
new file mode 100644
index 00000000..4b98c88e
--- /dev/null
+++ b/external/MaxmindDB/maxminddb.c
@@ -0,0 +1,2045 @@
+#if HAVE_CONFIG_H
+#include
+#endif
+#include "maxminddb.h"
+#include "maxminddb-compat-util.h"
+#include
+#include
+#include
+#include
+#include
+#include
+#include
+#include
+
+#ifdef _WIN32
+#include
+#include
+#else
+#include
+#include
+#include
+#endif
+
+#define MMDB_DATA_SECTION_SEPARATOR (16)
+#define MAXIMUM_DATA_STRUCTURE_DEPTH (512)
+
+#ifdef MMDB_DEBUG
+#define LOCAL
+#define NO_PROTO
+#define DEBUG_FUNC
+#define DEBUG_MSG(msg) fprintf(stderr, msg "\n")
+#define DEBUG_MSGF(fmt, ...) fprintf(stderr, fmt "\n", __VA_ARGS__)
+#define DEBUG_BINARY(fmt, byte) \
+ do { \
+ char *binary = byte_to_binary(byte); \
+ if (NULL == binary) { \
+ fprintf(stderr, "Malloc failed in DEBUG_BINARY\n"); \
+ abort(); \
+ } \
+ fprintf(stderr, fmt "\n", binary); \
+ free(binary); \
+ } while (0)
+#define DEBUG_NL fprintf(stderr, "\n")
+#else
+#define LOCAL static
+#define NO_PROTO static
+#define DEBUG_MSG(...)
+#define DEBUG_MSGF(...)
+#define DEBUG_BINARY(...)
+#define DEBUG_NL
+#endif
+
+#ifdef MMDB_DEBUG
+DEBUG_FUNC char *byte_to_binary(uint8_t byte)
+{
+ char *bits = malloc(sizeof(char) * 9);
+ if (NULL == bits) {
+ return bits;
+ }
+
+ for (uint8_t i = 0; i < 8; i++) {
+ bits[i] = byte & (128 >> i) ? '1' : '0';
+ }
+ bits[8] = '\0';
+
+ return bits;
+}
+
+DEBUG_FUNC char *type_num_to_name(uint8_t num)
+{
+ switch (num) {
+ case 0:
+ return "extended";
+ case 1:
+ return "pointer";
+ case 2:
+ return "utf8_string";
+ case 3:
+ return "double";
+ case 4:
+ return "bytes";
+ case 5:
+ return "uint16";
+ case 6:
+ return "uint32";
+ case 7:
+ return "map";
+ case 8:
+ return "int32";
+ case 9:
+ return "uint64";
+ case 10:
+ return "uint128";
+ case 11:
+ return "array";
+ case 12:
+ return "container";
+ case 13:
+ return "end_marker";
+ case 14:
+ return "boolean";
+ case 15:
+ return "float";
+ default:
+ return "unknown type";
+ }
+}
+#endif
+
+/* None of the values we check on the lhs are bigger than uint32_t, so on
+ * platforms where SIZE_MAX is a 64-bit integer, this would be a no-op, and it
+ * makes the compiler complain if we do the check anyway. */
+#if SIZE_MAX == UINT32_MAX
+#define MAYBE_CHECK_SIZE_OVERFLOW(lhs, rhs, error) \
+ if ((lhs) > (rhs)) { \
+ return error; \
+ }
+#else
+#define MAYBE_CHECK_SIZE_OVERFLOW(...)
+#endif
+
+typedef struct record_info_s {
+ uint16_t record_length;
+ uint32_t (*left_record_getter)(const uint8_t *);
+ uint32_t (*right_record_getter)(const uint8_t *);
+ uint8_t right_record_offset;
+} record_info_s;
+
+#define METADATA_MARKER "\xab\xcd\xefMaxMind.com"
+/* This is 128kb */
+#define METADATA_BLOCK_MAX_SIZE 131072
+
+/* *INDENT-OFF* */
+/* --prototypes automatically generated by dev-bin/regen-prototypes.pl - don't remove this comment */
+LOCAL int map_file(MMDB_s *const mmdb);
+LOCAL const uint8_t *find_metadata(const uint8_t *file_content,
+ ssize_t file_size, uint32_t *metadata_size);
+LOCAL int read_metadata(MMDB_s *mmdb);
+LOCAL MMDB_s make_fake_metadata_db(MMDB_s *mmdb);
+LOCAL int value_for_key_as_uint16(MMDB_entry_s *start, char *key,
+ uint16_t *value);
+LOCAL int value_for_key_as_uint32(MMDB_entry_s *start, char *key,
+ uint32_t *value);
+LOCAL int value_for_key_as_uint64(MMDB_entry_s *start, char *key,
+ uint64_t *value);
+LOCAL int value_for_key_as_string(MMDB_entry_s *start, char *key,
+ char const **value);
+LOCAL int populate_languages_metadata(MMDB_s *mmdb, MMDB_s *metadata_db,
+ MMDB_entry_s *metadata_start);
+LOCAL int populate_description_metadata(MMDB_s *mmdb, MMDB_s *metadata_db,
+ MMDB_entry_s *metadata_start);
+LOCAL int resolve_any_address(const char *ipstr, struct addrinfo **addresses);
+LOCAL int find_address_in_search_tree(MMDB_s *mmdb, uint8_t *address,
+ sa_family_t address_family,
+ MMDB_lookup_result_s *result);
+LOCAL record_info_s record_info_for_database(MMDB_s *mmdb);
+LOCAL int find_ipv4_start_node(MMDB_s *mmdb);
+LOCAL int populate_result(MMDB_s *mmdb, uint32_t node_count, uint32_t value,
+ uint16_t netmask, MMDB_lookup_result_s *result);
+LOCAL uint32_t get_left_28_bit_record(const uint8_t *record);
+LOCAL uint32_t get_right_28_bit_record(const uint8_t *record);
+LOCAL int path_length(va_list va_path);
+LOCAL int lookup_path_in_array(const char *path_elem, MMDB_s *mmdb,
+ MMDB_entry_data_s *entry_data);
+LOCAL int lookup_path_in_map(const char *path_elem, MMDB_s *mmdb,
+ MMDB_entry_data_s *entry_data);
+LOCAL int skip_map_or_array(MMDB_s *mmdb, MMDB_entry_data_s *entry_data);
+LOCAL int decode_one_follow(MMDB_s *mmdb, uint32_t offset,
+ MMDB_entry_data_s *entry_data);
+LOCAL int decode_one(MMDB_s *mmdb, uint32_t offset,
+ MMDB_entry_data_s *entry_data);
+LOCAL int get_ext_type(int raw_ext_type);
+LOCAL uint32_t get_ptr_from(uint8_t ctrl, uint8_t const *const ptr,
+ int ptr_size);
+LOCAL int get_entry_data_list(MMDB_s *mmdb, uint32_t offset,
+ MMDB_entry_data_list_s *const entry_data_list,
+ int depth);
+LOCAL float get_ieee754_float(const uint8_t *restrict p);
+LOCAL double get_ieee754_double(const uint8_t *restrict p);
+LOCAL uint32_t get_uint32(const uint8_t *p);
+LOCAL uint32_t get_uint24(const uint8_t *p);
+LOCAL uint32_t get_uint16(const uint8_t *p);
+LOCAL uint64_t get_uintX(const uint8_t *p, int length);
+LOCAL int32_t get_sintX(const uint8_t *p, int length);
+LOCAL MMDB_entry_data_list_s *new_entry_data_list(void);
+LOCAL void free_mmdb_struct(MMDB_s *const mmdb);
+LOCAL void free_languages_metadata(MMDB_s *mmdb);
+LOCAL void free_descriptions_metadata(MMDB_s *mmdb);
+LOCAL MMDB_entry_data_list_s *dump_entry_data_list(
+ FILE *stream, MMDB_entry_data_list_s *entry_data_list, int indent,
+ int *status);
+LOCAL void print_indentation(FILE *stream, int i);
+LOCAL char *bytes_to_hex(uint8_t *bytes, uint32_t size);
+/* --prototypes end - don't remove this comment-- */
+/* *INDENT-ON* */
+
+#define CHECKED_DECODE_ONE(mmdb, offset, entry_data) \
+ do { \
+ int status = decode_one(mmdb, offset, entry_data); \
+ if (MMDB_SUCCESS != status) { \
+ return status; \
+ } \
+ } while (0)
+
+#define CHECKED_DECODE_ONE_FOLLOW(mmdb, offset, entry_data) \
+ do { \
+ int status = decode_one_follow(mmdb, offset, entry_data); \
+ if (MMDB_SUCCESS != status) { \
+ return status; \
+ } \
+ } while (0)
+
+#define FREE_AND_SET_NULL(p) { free((void *)(p)); (p) = NULL; }
+
+int MMDB_open(const char *const filename, uint32_t flags, MMDB_s *const mmdb)
+{
+ int status = MMDB_SUCCESS;
+
+ mmdb->file_content = NULL;
+ mmdb->data_section = NULL;
+ mmdb->metadata.database_type = NULL;
+ mmdb->metadata.languages.count = 0;
+ mmdb->metadata.description.count = 0;
+
+ mmdb->filename = mmdb_strdup(filename);
+ if (NULL == mmdb->filename) {
+ status = MMDB_OUT_OF_MEMORY_ERROR;
+ goto cleanup;
+ }
+
+ if ((flags & MMDB_MODE_MASK) == 0) {
+ flags |= MMDB_MODE_MMAP;
+ }
+ mmdb->flags = flags;
+
+ if (MMDB_SUCCESS != (status = map_file(mmdb)) ) {
+ goto cleanup;
+ }
+
+#ifdef _WIN32
+ WSADATA wsa;
+ WSAStartup(MAKEWORD(2, 2), &wsa);
+#endif
+
+ uint32_t metadata_size = 0;
+ const uint8_t *metadata = find_metadata(mmdb->file_content, mmdb->file_size,
+ &metadata_size);
+ if (NULL == metadata) {
+ status = MMDB_INVALID_METADATA_ERROR;
+ goto cleanup;
+ }
+
+ mmdb->metadata_section = metadata;
+ mmdb->metadata_section_size = metadata_size;
+
+ status = read_metadata(mmdb);
+ if (MMDB_SUCCESS != status) {
+ goto cleanup;
+ }
+
+ if (mmdb->metadata.binary_format_major_version != 2) {
+ status = MMDB_UNKNOWN_DATABASE_FORMAT_ERROR;
+ goto cleanup;
+ }
+
+ uint32_t search_tree_size = mmdb->metadata.node_count *
+ mmdb->full_record_byte_size;
+
+ mmdb->data_section = mmdb->file_content + search_tree_size;
+ if (search_tree_size > (uint32_t)mmdb->file_size) {
+ status = MMDB_INVALID_METADATA_ERROR;
+ goto cleanup;
+ }
+ mmdb->data_section_size = mmdb->file_size - search_tree_size;
+ mmdb->metadata_section = metadata;
+ mmdb->ipv4_start_node.node_value = 0;
+ mmdb->ipv4_start_node.netmask = 0;
+
+ cleanup:
+ if (MMDB_SUCCESS != status) {
+ int saved_errno = errno;
+ free_mmdb_struct(mmdb);
+ errno = saved_errno;
+ }
+ return status;
+}
+
+#ifdef _WIN32
+
+LOCAL int map_file(MMDB_s *const mmdb)
+{
+ ssize_t size;
+ int status = MMDB_SUCCESS;
+ HANDLE mmh = NULL;
+ HANDLE fd = CreateFileA(mmdb->filename, GENERIC_READ, FILE_SHARE_READ, NULL,
+ OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL, NULL);
+ if (fd == INVALID_HANDLE_VALUE) {
+ status = MMDB_FILE_OPEN_ERROR;
+ goto cleanup;
+ }
+ size = GetFileSize(fd, NULL);
+ if (size == INVALID_FILE_SIZE) {
+ status = MMDB_FILE_OPEN_ERROR;
+ goto cleanup;
+ }
+ mmh = CreateFileMappingA(fd, NULL, PAGE_READONLY, 0, size, NULL);
+ /* Microsoft documentation for CreateFileMapping indicates this returns
+ NULL not INVALID_HANDLE_VALUE on error */
+ if (NULL == mmh) {
+ status = MMDB_IO_ERROR;
+ goto cleanup;
+ }
+ uint8_t *file_content =
+ (uint8_t *)MapViewOfFile(mmh, FILE_MAP_READ, 0, 0, 0);
+ if (file_content == NULL) {
+ status = MMDB_IO_ERROR;
+ goto cleanup;
+ }
+
+ mmdb->file_size = size;
+ mmdb->file_content = file_content;
+
+ cleanup:;
+ int saved_errno = errno;
+ if (INVALID_HANDLE_VALUE != fd) {
+ CloseHandle(fd);
+ }
+ if (NULL != mmh) {
+ CloseHandle(mmh);
+ }
+ errno = saved_errno;
+
+ return status;
+}
+
+#else
+
+LOCAL int map_file(MMDB_s *const mmdb)
+{
+ ssize_t size;
+ int status = MMDB_SUCCESS;
+
+ int fd = open(mmdb->filename, O_RDONLY);
+ struct stat s;
+ if (fd < 0 || fstat(fd, &s)) {
+ status = MMDB_FILE_OPEN_ERROR;
+ goto cleanup;
+ }
+
+ size = s.st_size;
+ if (size < 0 || size != s.st_size) {
+ status = MMDB_OUT_OF_MEMORY_ERROR;
+ goto cleanup;
+ }
+
+ uint8_t *file_content =
+ (uint8_t *)mmap(NULL, size, PROT_READ, MAP_SHARED, fd, 0);
+ if (MAP_FAILED == file_content) {
+ if (ENOMEM == errno) {
+ status = MMDB_OUT_OF_MEMORY_ERROR;
+ } else {
+ status = MMDB_IO_ERROR;
+ }
+ goto cleanup;
+ }
+
+ mmdb->file_size = size;
+ mmdb->file_content = file_content;
+
+ cleanup:;
+ int saved_errno = errno;
+ if (fd >= 0) {
+ close(fd);
+ }
+ errno = saved_errno;
+
+ return status;
+}
+
+#endif
+
+LOCAL const uint8_t *find_metadata(const uint8_t *file_content,
+ ssize_t file_size, uint32_t *metadata_size)
+{
+ const ssize_t marker_len = sizeof(METADATA_MARKER) - 1;
+ ssize_t max_size = file_size >
+ METADATA_BLOCK_MAX_SIZE ? METADATA_BLOCK_MAX_SIZE :
+ file_size;
+
+ uint8_t *search_area = (uint8_t *)(file_content + (file_size - max_size));
+ uint8_t *start = search_area;
+ uint8_t *tmp;
+ do {
+ tmp = mmdb_memmem(search_area, max_size,
+ METADATA_MARKER, marker_len);
+
+ if (NULL != tmp) {
+ max_size -= tmp - search_area;
+ search_area = tmp;
+
+ /* Continue searching just after the marker we just read, in case
+ * there are multiple markers in the same file. This would be odd
+ * but is certainly not impossible. */
+ max_size -= marker_len;
+ search_area += marker_len;
+ }
+ } while (NULL != tmp);
+
+ if (search_area == start) {
+ return NULL;
+ }
+
+ *metadata_size = max_size;
+
+ return search_area;
+}
+
+LOCAL int read_metadata(MMDB_s *mmdb)
+{
+ /* We need to create a fake MMDB_s struct in order to decode values from
+ the metadata. The metadata is basically just like the data section, so we
+ want to use the same functions we use for the data section to get metadata
+ values. */
+ MMDB_s metadata_db = make_fake_metadata_db(mmdb);
+
+ MMDB_entry_s metadata_start = {
+ .mmdb = &metadata_db,
+ .offset = 0
+ };
+
+ int status =
+ value_for_key_as_uint32(&metadata_start, "node_count",
+ &mmdb->metadata.node_count);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ if (!mmdb->metadata.node_count) {
+ DEBUG_MSG("could not find node_count value in metadata");
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+
+ status = value_for_key_as_uint16(&metadata_start, "record_size",
+ &mmdb->metadata.record_size);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ if (!mmdb->metadata.record_size) {
+ DEBUG_MSG("could not find record_size value in metadata");
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+
+ if (mmdb->metadata.record_size != 24 && mmdb->metadata.record_size != 28
+ && mmdb->metadata.record_size != 32) {
+ DEBUG_MSGF("bad record size in metadata: %i",
+ mmdb->metadata.record_size);
+ return MMDB_UNKNOWN_DATABASE_FORMAT_ERROR;
+ }
+
+ status = value_for_key_as_uint16(&metadata_start, "ip_version",
+ &mmdb->metadata.ip_version);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ if (!mmdb->metadata.ip_version) {
+ DEBUG_MSG("could not find ip_version value in metadata");
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+ if (!(mmdb->metadata.ip_version == 4 || mmdb->metadata.ip_version == 6)) {
+ DEBUG_MSGF("ip_version value in metadata is not 4 or 6 - it was %i",
+ mmdb->metadata.ip_version);
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+
+ status = value_for_key_as_string(&metadata_start, "database_type",
+ &mmdb->metadata.database_type);
+ if (MMDB_SUCCESS != status) {
+ DEBUG_MSG("error finding database_type value in metadata");
+ return status;
+ }
+
+ status =
+ populate_languages_metadata(mmdb, &metadata_db, &metadata_start);
+ if (MMDB_SUCCESS != status) {
+ DEBUG_MSG("could not populate languages from metadata");
+ return status;
+ }
+
+ status = value_for_key_as_uint16(
+ &metadata_start, "binary_format_major_version",
+ &mmdb->metadata.binary_format_major_version);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ if (!mmdb->metadata.binary_format_major_version) {
+ DEBUG_MSG(
+ "could not find binary_format_major_version value in metadata");
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+
+ status = value_for_key_as_uint16(
+ &metadata_start, "binary_format_minor_version",
+ &mmdb->metadata.binary_format_minor_version);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+
+ status = value_for_key_as_uint64(&metadata_start, "build_epoch",
+ &mmdb->metadata.build_epoch);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ if (!mmdb->metadata.build_epoch) {
+ DEBUG_MSG("could not find build_epoch value in metadata");
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+
+ status = populate_description_metadata(mmdb, &metadata_db, &metadata_start);
+ if (MMDB_SUCCESS != status) {
+ DEBUG_MSG("could not populate description from metadata");
+ return status;
+ }
+
+ mmdb->full_record_byte_size = mmdb->metadata.record_size * 2 / 8U;
+
+ mmdb->depth = mmdb->metadata.ip_version == 4 ? 32 : 128;
+
+ return MMDB_SUCCESS;
+}
+
+LOCAL MMDB_s make_fake_metadata_db(MMDB_s *mmdb)
+{
+ MMDB_s fake_metadata_db = {
+ .data_section = mmdb->metadata_section,
+ .data_section_size = mmdb->metadata_section_size
+ };
+
+ return fake_metadata_db;
+}
+
+LOCAL int value_for_key_as_uint16(MMDB_entry_s *start, char *key,
+ uint16_t *value)
+{
+ MMDB_entry_data_s entry_data;
+ const char *path[] = { key, NULL };
+ int status = MMDB_aget_value(start, &entry_data, path);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ if (MMDB_DATA_TYPE_UINT16 != entry_data.type) {
+ DEBUG_MSGF("expect uint16 for %s but received %s", key,
+ type_num_to_name(
+ entry_data.type));
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+ *value = entry_data.uint16;
+ return MMDB_SUCCESS;
+}
+
+LOCAL int value_for_key_as_uint32(MMDB_entry_s *start, char *key,
+ uint32_t *value)
+{
+ MMDB_entry_data_s entry_data;
+ const char *path[] = { key, NULL };
+ int status = MMDB_aget_value(start, &entry_data, path);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ if (MMDB_DATA_TYPE_UINT32 != entry_data.type) {
+ DEBUG_MSGF("expect uint32 for %s but received %s", key,
+ type_num_to_name(
+ entry_data.type));
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+ *value = entry_data.uint32;
+ return MMDB_SUCCESS;
+}
+
+LOCAL int value_for_key_as_uint64(MMDB_entry_s *start, char *key,
+ uint64_t *value)
+{
+ MMDB_entry_data_s entry_data;
+ const char *path[] = { key, NULL };
+ int status = MMDB_aget_value(start, &entry_data, path);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ if (MMDB_DATA_TYPE_UINT64 != entry_data.type) {
+ DEBUG_MSGF("expect uint64 for %s but received %s", key,
+ type_num_to_name(
+ entry_data.type));
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+ *value = entry_data.uint64;
+ return MMDB_SUCCESS;
+}
+
+LOCAL int value_for_key_as_string(MMDB_entry_s *start, char *key,
+ char const **value)
+{
+ MMDB_entry_data_s entry_data;
+ const char *path[] = { key, NULL };
+ int status = MMDB_aget_value(start, &entry_data, path);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ if (MMDB_DATA_TYPE_UTF8_STRING != entry_data.type) {
+ DEBUG_MSGF("expect string for %s but received %s", key,
+ type_num_to_name(
+ entry_data.type));
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+ *value = mmdb_strndup((char *)entry_data.utf8_string, entry_data.data_size);
+ if (NULL == *value) {
+ return MMDB_OUT_OF_MEMORY_ERROR;
+ }
+ return MMDB_SUCCESS;
+}
+
+LOCAL int populate_languages_metadata(MMDB_s *mmdb, MMDB_s *metadata_db,
+ MMDB_entry_s *metadata_start)
+{
+ MMDB_entry_data_s entry_data;
+
+ const char *path[] = { "languages", NULL };
+ int status = MMDB_aget_value(metadata_start, &entry_data, path);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ if (MMDB_DATA_TYPE_ARRAY != entry_data.type) {
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+
+ MMDB_entry_s array_start = {
+ .mmdb = metadata_db,
+ .offset = entry_data.offset
+ };
+
+ MMDB_entry_data_list_s *member;
+ status = MMDB_get_entry_data_list(&array_start, &member);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+
+ MMDB_entry_data_list_s *first_member = member;
+
+ uint32_t array_size = member->entry_data.data_size;
+ MAYBE_CHECK_SIZE_OVERFLOW(array_size, SIZE_MAX / sizeof(char *),
+ MMDB_INVALID_METADATA_ERROR);
+
+ mmdb->metadata.languages.count = 0;
+ mmdb->metadata.languages.names = malloc(array_size * sizeof(char *));
+ if (NULL == mmdb->metadata.languages.names) {
+ return MMDB_OUT_OF_MEMORY_ERROR;
+ }
+
+ for (uint32_t i = 0; i < array_size; i++) {
+ member = member->next;
+ if (MMDB_DATA_TYPE_UTF8_STRING != member->entry_data.type) {
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+
+ mmdb->metadata.languages.names[i] =
+ mmdb_strndup((char *)member->entry_data.utf8_string,
+ member->entry_data.data_size);
+
+ if (NULL == mmdb->metadata.languages.names[i]) {
+ return MMDB_OUT_OF_MEMORY_ERROR;
+ }
+ // We assign this as we go so that if we fail a malloc and need to
+ // free it, the count is right.
+ mmdb->metadata.languages.count = i + 1;
+ }
+
+ MMDB_free_entry_data_list(first_member);
+
+ return MMDB_SUCCESS;
+}
+
+LOCAL int populate_description_metadata(MMDB_s *mmdb, MMDB_s *metadata_db,
+ MMDB_entry_s *metadata_start)
+{
+ MMDB_entry_data_s entry_data;
+
+ const char *path[] = { "description", NULL };
+ int status = MMDB_aget_value(metadata_start, &entry_data, path);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+
+ if (MMDB_DATA_TYPE_MAP != entry_data.type) {
+ return MMDB_INVALID_METADATA_ERROR;
+ }
+
+ MMDB_entry_s map_start = {
+ .mmdb = metadata_db,
+ .offset = entry_data.offset
+ };
+
+ MMDB_entry_data_list_s *member;
+ status = MMDB_get_entry_data_list(&map_start, &member);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+
+ MMDB_entry_data_list_s *first_member = member;
+
+ uint32_t map_size = member->entry_data.data_size;
+ mmdb->metadata.description.count = 0;
+ if (0 == map_size) {
+ mmdb->metadata.description.descriptions = NULL;
+ goto cleanup;
+ }
+ MAYBE_CHECK_SIZE_OVERFLOW(map_size, SIZE_MAX / sizeof(MMDB_description_s *),
+ MMDB_INVALID_METADATA_ERROR);
+
+ mmdb->metadata.description.descriptions =
+ malloc(map_size * sizeof(MMDB_description_s *));
+ if (NULL == mmdb->metadata.description.descriptions) {
+ status = MMDB_OUT_OF_MEMORY_ERROR;
+ goto cleanup;
+ }
+
+ for (uint32_t i = 0; i < map_size; i++) {
+ mmdb->metadata.description.descriptions[i] =
+ malloc(sizeof(MMDB_description_s));
+ if (NULL == mmdb->metadata.description.descriptions[i]) {
+ status = MMDB_OUT_OF_MEMORY_ERROR;
+ goto cleanup;
+ }
+
+ mmdb->metadata.description.count = i + 1;
+ mmdb->metadata.description.descriptions[i]->language = NULL;
+ mmdb->metadata.description.descriptions[i]->description = NULL;
+
+ member = member->next;
+
+ if (MMDB_DATA_TYPE_UTF8_STRING != member->entry_data.type) {
+ status = MMDB_INVALID_METADATA_ERROR;
+ goto cleanup;
+ }
+
+ mmdb->metadata.description.descriptions[i]->language =
+ mmdb_strndup((char *)member->entry_data.utf8_string,
+ member->entry_data.data_size);
+
+ if (NULL == mmdb->metadata.description.descriptions[i]->language) {
+ status = MMDB_OUT_OF_MEMORY_ERROR;
+ goto cleanup;
+ }
+
+ member = member->next;
+
+ if (MMDB_DATA_TYPE_UTF8_STRING != member->entry_data.type) {
+ status = MMDB_INVALID_METADATA_ERROR;
+ goto cleanup;
+ }
+
+ mmdb->metadata.description.descriptions[i]->description =
+ mmdb_strndup((char *)member->entry_data.utf8_string,
+ member->entry_data.data_size);
+
+ if (NULL == mmdb->metadata.description.descriptions[i]->description) {
+ status = MMDB_OUT_OF_MEMORY_ERROR;
+ goto cleanup;
+ }
+ }
+
+ cleanup:
+ MMDB_free_entry_data_list(first_member);
+
+ return status;
+}
+
+MMDB_lookup_result_s MMDB_lookup_string(MMDB_s *const mmdb,
+ const char *const ipstr,
+ int *const gai_error,
+ int *const mmdb_error)
+{
+ MMDB_lookup_result_s result = {
+ .found_entry = false,
+ .netmask = 0,
+ .entry = {
+ .mmdb = mmdb,
+ .offset = 0
+ }
+ };
+
+ struct addrinfo *addresses = NULL;
+ *gai_error = resolve_any_address(ipstr, &addresses);
+
+ if (!*gai_error) {
+ result = MMDB_lookup_sockaddr(mmdb, addresses->ai_addr, mmdb_error);
+ }
+
+ if (NULL != addresses) {
+ freeaddrinfo(addresses);
+ }
+
+ return result;
+}
+
+LOCAL int resolve_any_address(const char *ipstr, struct addrinfo **addresses)
+{
+ struct addrinfo hints = {
+ .ai_family = AF_UNSPEC,
+ .ai_flags = AI_NUMERICHOST,
+ // We set ai_socktype so that we only get one result back
+ .ai_socktype = SOCK_STREAM
+ };
+
+ int gai_status = getaddrinfo(ipstr, NULL, &hints, addresses);
+ if (gai_status) {
+ return gai_status;
+ }
+
+ return 0;
+}
+
+MMDB_lookup_result_s MMDB_lookup_sockaddr(
+ MMDB_s *const mmdb,
+ const struct sockaddr *const sockaddr,
+ int *const mmdb_error)
+{
+ MMDB_lookup_result_s result = {
+ .found_entry = false,
+ .netmask = 0,
+ .entry = {
+ .mmdb = mmdb,
+ .offset = 0
+ }
+ };
+
+ uint8_t mapped_address[16], *address;
+ if (mmdb->metadata.ip_version == 4) {
+ if (sockaddr->sa_family == AF_INET6) {
+ *mmdb_error = MMDB_IPV6_LOOKUP_IN_IPV4_DATABASE_ERROR;
+ return result;
+ }
+ address = (uint8_t *)&((struct sockaddr_in *)sockaddr)->sin_addr.s_addr;
+ } else {
+ if (sockaddr->sa_family == AF_INET6) {
+ address =
+ (uint8_t *)&((struct sockaddr_in6 *)sockaddr)->sin6_addr.
+ s6_addr;
+ } else {
+ address = mapped_address;
+ memset(address, 0, 12);
+ memcpy(address + 12,
+ &((struct sockaddr_in *)sockaddr)->sin_addr.s_addr, 4);
+ }
+ }
+
+ *mmdb_error =
+ find_address_in_search_tree(mmdb, address, sockaddr->sa_family,
+ &result);
+
+ return result;
+}
+
+LOCAL int find_address_in_search_tree(MMDB_s *mmdb, uint8_t *address,
+ sa_family_t address_family,
+ MMDB_lookup_result_s *result)
+{
+ record_info_s record_info = record_info_for_database(mmdb);
+ if (0 == record_info.right_record_offset) {
+ return MMDB_UNKNOWN_DATABASE_FORMAT_ERROR;
+ }
+
+ DEBUG_NL;
+ DEBUG_MSG("Looking for address in search tree");
+
+ uint32_t node_count = mmdb->metadata.node_count;
+ uint32_t value = 0;
+ uint16_t max_depth0 = mmdb->depth - 1;
+ uint16_t start_bit = max_depth0;
+
+ if (mmdb->metadata.ip_version == 6 && address_family == AF_INET) {
+ int mmdb_error = find_ipv4_start_node(mmdb);
+ if (MMDB_SUCCESS != mmdb_error) {
+ return mmdb_error;
+ }
+ DEBUG_MSGF("IPv4 start node is %u (netmask %u)",
+ mmdb->ipv4_start_node.node_value,
+ mmdb->ipv4_start_node.netmask);
+ /* We have an IPv6 database with no IPv4 data */
+ if (mmdb->ipv4_start_node.node_value >= node_count) {
+ return populate_result(mmdb, node_count,
+ mmdb->ipv4_start_node.node_value,
+ mmdb->ipv4_start_node.netmask,
+ result);
+ }
+
+ value = mmdb->ipv4_start_node.node_value;
+ start_bit -= mmdb->ipv4_start_node.netmask;
+ }
+
+ const uint8_t *search_tree = mmdb->file_content;
+ const uint8_t *record_pointer;
+ for (int current_bit = start_bit; current_bit >= 0; current_bit--) {
+ uint8_t bit_is_true =
+ address[(max_depth0 - current_bit) >> 3]
+ & (1U << (~(max_depth0 - current_bit) & 7)) ? 1 : 0;
+
+ DEBUG_MSGF("Looking at bit %i - bit's value is %i", current_bit,
+ bit_is_true);
+ DEBUG_MSGF(" current node = %u", value);
+
+ record_pointer = &search_tree[value * record_info.record_length];
+ if (record_pointer + record_info.record_length > mmdb->data_section) {
+ return MMDB_CORRUPT_SEARCH_TREE_ERROR;
+ }
+ if (bit_is_true) {
+ record_pointer += record_info.right_record_offset;
+ value = record_info.right_record_getter(record_pointer);
+ } else {
+ value = record_info.left_record_getter(record_pointer);
+ }
+
+ /* Ideally we'd check to make sure that a record never points to a
+ * previously seen value, but that's more complicated. For now, we can
+ * at least check that we don't end up at the top of the tree again. */
+ if (0 == value) {
+ DEBUG_MSGF(" %s record has a value of 0",
+ bit_is_true ? "right" : "left");
+ return MMDB_CORRUPT_SEARCH_TREE_ERROR;
+ }
+
+ if (value >= node_count) {
+ return populate_result(mmdb, node_count, value, current_bit, result);
+ } else {
+ DEBUG_MSGF(" proceeding to search tree node %i", value);
+ }
+ }
+
+ DEBUG_MSG(
+ "Reached the end of the address bits without leaving the search tree");
+
+ // We should not be able to reach this return. If we do, something very bad happened.
+ return MMDB_CORRUPT_SEARCH_TREE_ERROR;
+}
+
+LOCAL record_info_s record_info_for_database(MMDB_s *mmdb)
+{
+ record_info_s record_info = {
+ .record_length = mmdb->full_record_byte_size,
+ .right_record_offset = 0
+ };
+
+ if (record_info.record_length == 6) {
+ record_info.left_record_getter = &get_uint24;
+ record_info.right_record_getter = &get_uint24;
+ record_info.right_record_offset = 3;
+ } else if (record_info.record_length == 7) {
+ record_info.left_record_getter = &get_left_28_bit_record;
+ record_info.right_record_getter = &get_right_28_bit_record;
+ record_info.right_record_offset = 3;
+ } else if (record_info.record_length == 8) {
+ record_info.left_record_getter = &get_uint32;
+ record_info.right_record_getter = &get_uint32;
+ record_info.right_record_offset = 4;
+ } else {
+ assert(false);
+ }
+
+ return record_info;
+}
+
+LOCAL int find_ipv4_start_node(MMDB_s *mmdb)
+{
+ /* In a pathological case of a database with a single node search tree,
+ * this check will be true even after we've found the IPv4 start node, but
+ * that doesn't seem worth trying to fix. */
+ if (mmdb->ipv4_start_node.node_value != 0) {
+ return MMDB_SUCCESS;
+ }
+
+ record_info_s record_info = record_info_for_database(mmdb);
+
+ const uint8_t *search_tree = mmdb->file_content;
+ uint32_t node_value = 0;
+ const uint8_t *record_pointer;
+ uint32_t netmask;
+ for (netmask = 0; netmask < 96; netmask++) {
+ record_pointer = &search_tree[node_value * record_info.record_length];
+ if (record_pointer + record_info.record_length > mmdb->data_section) {
+ return MMDB_CORRUPT_SEARCH_TREE_ERROR;
+ }
+ node_value = record_info.left_record_getter(record_pointer);
+ /* This can happen if there's no IPv4 data _or_ if there is a subnet
+ * with data that contains the entire IPv4 range (like ::/64) */
+ if (node_value >= mmdb->metadata.node_count) {
+ break;
+ }
+ }
+
+ mmdb->ipv4_start_node.node_value = node_value;
+ mmdb->ipv4_start_node.netmask = netmask;
+
+ return MMDB_SUCCESS;
+}
+
+LOCAL int populate_result(MMDB_s *mmdb, uint32_t node_count, uint32_t value,
+ uint16_t netmask, MMDB_lookup_result_s *result)
+{
+ uint32_t offset = value - node_count;
+ DEBUG_MSGF(" data section offset is %i (record value = %i)", offset, value);
+
+ if (offset > mmdb->data_section_size) {
+ return MMDB_CORRUPT_SEARCH_TREE_ERROR;
+ }
+
+ result->netmask = mmdb->depth - netmask;
+ result->entry.offset = offset;
+ result->found_entry = result->entry.offset > 0 ? true : false;
+ return MMDB_SUCCESS;
+}
+
+LOCAL uint32_t get_left_28_bit_record(const uint8_t *record)
+{
+ return record[0] * 65536 + record[1] * 256 + record[2] +
+ ((record[3] & 0xf0) << 20);
+}
+
+LOCAL uint32_t get_right_28_bit_record(const uint8_t *record)
+{
+ uint32_t value = get_uint32(record);
+ return value & 0xfffffff;
+}
+
+int MMDB_read_node(MMDB_s *const mmdb, uint32_t node_number,
+ MMDB_search_node_s *const node)
+{
+ record_info_s record_info = record_info_for_database(mmdb);
+ if (0 == record_info.right_record_offset) {
+ return MMDB_UNKNOWN_DATABASE_FORMAT_ERROR;
+ }
+
+ if (node_number > mmdb->metadata.node_count) {
+ return MMDB_INVALID_NODE_NUMBER_ERROR;
+ }
+
+ const uint8_t *search_tree = mmdb->file_content;
+ const uint8_t *record_pointer =
+ &search_tree[node_number * record_info.record_length];
+ node->left_record = record_info.left_record_getter(record_pointer);
+ record_pointer += record_info.right_record_offset;
+ node->right_record = record_info.right_record_getter(record_pointer);
+
+ return MMDB_SUCCESS;
+}
+
+int MMDB_get_value(MMDB_entry_s *const start,
+ MMDB_entry_data_s *const entry_data,
+ ...)
+{
+ va_list path;
+ va_start(path, entry_data);
+ int status = MMDB_vget_value(start, entry_data, path);
+ va_end(path);
+ return status;
+}
+
+int MMDB_vget_value(MMDB_entry_s *const start,
+ MMDB_entry_data_s *const entry_data,
+ va_list va_path)
+{
+ int length = path_length(va_path);
+ const char *path_elem;
+ int i = 0;
+
+ MAYBE_CHECK_SIZE_OVERFLOW(length, SIZE_MAX / sizeof(const char *) - 1,
+ MMDB_INVALID_METADATA_ERROR);
+
+ const char **path = malloc((length + 1) * sizeof(const char *));
+ if (NULL == path) {
+ return MMDB_OUT_OF_MEMORY_ERROR;
+ }
+
+ while (NULL != (path_elem = va_arg(va_path, char *))) {
+ path[i] = path_elem;
+ i++;
+ }
+ path[i] = NULL;
+
+ int status = MMDB_aget_value(start, entry_data, path);
+
+ free((char **)path);
+
+ return status;
+}
+
+LOCAL int path_length(va_list va_path)
+{
+ int i = 0;
+ const char *ignore;
+ va_list path_copy;
+ va_copy(path_copy, va_path);
+
+ while (NULL != (ignore = va_arg(path_copy, char *))) {
+ i++;
+ }
+
+ va_end(path_copy);
+
+ return i;
+}
+
+int MMDB_aget_value(MMDB_entry_s *const start,
+ MMDB_entry_data_s *const entry_data,
+ const char *const *const path)
+{
+ MMDB_s *mmdb = start->mmdb;
+ uint32_t offset = start->offset;
+
+ memset(entry_data, 0, sizeof(MMDB_entry_data_s));
+ DEBUG_NL;
+ DEBUG_MSG("looking up value by path");
+
+ CHECKED_DECODE_ONE_FOLLOW(mmdb, offset, entry_data);
+
+ DEBUG_NL;
+ DEBUG_MSGF("top level element is a %s", type_num_to_name(entry_data->type));
+
+ /* Can this happen? It'd probably represent a pathological case under
+ * normal use, but there's nothing preventing someone from passing an
+ * invalid MMDB_entry_s struct to this function */
+ if (!entry_data->has_data) {
+ return MMDB_INVALID_LOOKUP_PATH_ERROR;
+ }
+
+ const char *path_elem;
+ int i = 0;
+ while (NULL != (path_elem = path[i++])) {
+ DEBUG_NL;
+ DEBUG_MSGF("path elem = %s", path_elem);
+
+ /* XXX - it'd be good to find a quicker way to skip through these
+ entries that doesn't involve decoding them
+ completely. Basically we need to just use the size from the
+ control byte to advance our pointer rather than calling
+ decode_one(). */
+ if (entry_data->type == MMDB_DATA_TYPE_ARRAY) {
+ int status = lookup_path_in_array(path_elem, mmdb, entry_data);
+ if (MMDB_SUCCESS != status) {
+ memset(entry_data, 0, sizeof(MMDB_entry_data_s));
+ return status;
+ }
+ } else if (entry_data->type == MMDB_DATA_TYPE_MAP) {
+ int status = lookup_path_in_map(path_elem, mmdb, entry_data);
+ if (MMDB_SUCCESS != status) {
+ memset(entry_data, 0, sizeof(MMDB_entry_data_s));
+ return status;
+ }
+ } else {
+ /* Once we make the code traverse maps & arrays without calling
+ * decode_one() we can get rid of this. */
+ memset(entry_data, 0, sizeof(MMDB_entry_data_s));
+ return MMDB_LOOKUP_PATH_DOES_NOT_MATCH_DATA_ERROR;
+ }
+ }
+
+ return MMDB_SUCCESS;
+}
+
+LOCAL int lookup_path_in_array(const char *path_elem, MMDB_s *mmdb,
+ MMDB_entry_data_s *entry_data)
+{
+ uint32_t size = entry_data->data_size;
+ char *first_invalid;
+
+ int saved_errno = errno;
+ errno = 0;
+ int array_index = strtol(path_elem, &first_invalid, 10);
+ if (array_index < 0 || ERANGE == errno) {
+ errno = saved_errno;
+ return MMDB_INVALID_LOOKUP_PATH_ERROR;
+ }
+ errno = saved_errno;
+
+ if (*first_invalid || (uint32_t)array_index >= size) {
+ return MMDB_LOOKUP_PATH_DOES_NOT_MATCH_DATA_ERROR;
+ }
+
+ for (int i = 0; i < array_index; i++) {
+ /* We don't want to follow a pointer here. If the next element is a
+ * pointer we simply skip it and keep going */
+ CHECKED_DECODE_ONE(mmdb, entry_data->offset_to_next, entry_data);
+ int status = skip_map_or_array(mmdb, entry_data);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ }
+
+ MMDB_entry_data_s value;
+ CHECKED_DECODE_ONE_FOLLOW(mmdb, entry_data->offset_to_next, &value);
+ memcpy(entry_data, &value, sizeof(MMDB_entry_data_s));
+
+ return MMDB_SUCCESS;
+}
+
+LOCAL int lookup_path_in_map(const char *path_elem, MMDB_s *mmdb,
+ MMDB_entry_data_s *entry_data)
+{
+ uint32_t size = entry_data->data_size;
+ uint32_t offset = entry_data->offset_to_next;
+ size_t path_elem_len = strlen(path_elem);
+
+ while (size-- > 0) {
+ MMDB_entry_data_s key, value;
+ CHECKED_DECODE_ONE_FOLLOW(mmdb, offset, &key);
+
+ uint32_t offset_to_value = key.offset_to_next;
+
+ if (MMDB_DATA_TYPE_UTF8_STRING != key.type) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+
+ if (key.data_size == path_elem_len &&
+ !memcmp(path_elem, key.utf8_string, path_elem_len)) {
+
+ DEBUG_MSG("found key matching path elem");
+
+ CHECKED_DECODE_ONE_FOLLOW(mmdb, offset_to_value, &value);
+ memcpy(entry_data, &value, sizeof(MMDB_entry_data_s));
+ return MMDB_SUCCESS;
+ } else {
+ /* We don't want to follow a pointer here. If the next element is
+ * a pointer we simply skip it and keep going */
+ CHECKED_DECODE_ONE(mmdb, offset_to_value, &value);
+ int status = skip_map_or_array(mmdb, &value);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ offset = value.offset_to_next;
+ }
+ }
+
+ memset(entry_data, 0, sizeof(MMDB_entry_data_s));
+ return MMDB_LOOKUP_PATH_DOES_NOT_MATCH_DATA_ERROR;
+}
+
+LOCAL int skip_map_or_array(MMDB_s *mmdb, MMDB_entry_data_s *entry_data)
+{
+ if (entry_data->type == MMDB_DATA_TYPE_MAP) {
+ uint32_t size = entry_data->data_size;
+ while (size-- > 0) {
+ CHECKED_DECODE_ONE(mmdb, entry_data->offset_to_next, entry_data); // key
+ CHECKED_DECODE_ONE(mmdb, entry_data->offset_to_next, entry_data); // value
+ int status = skip_map_or_array(mmdb, entry_data);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ }
+ } else if (entry_data->type == MMDB_DATA_TYPE_ARRAY) {
+ uint32_t size = entry_data->data_size;
+ while (size-- > 0) {
+ CHECKED_DECODE_ONE(mmdb, entry_data->offset_to_next, entry_data); // value
+ int status = skip_map_or_array(mmdb, entry_data);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ }
+ }
+
+ return MMDB_SUCCESS;
+}
+
+LOCAL int decode_one_follow(MMDB_s *mmdb, uint32_t offset,
+ MMDB_entry_data_s *entry_data)
+{
+ CHECKED_DECODE_ONE(mmdb, offset, entry_data);
+ if (entry_data->type == MMDB_DATA_TYPE_POINTER) {
+ uint32_t next = entry_data->offset_to_next;
+ CHECKED_DECODE_ONE(mmdb, entry_data->pointer, entry_data);
+ /* Pointers to pointers are illegal under the spec */
+ if (entry_data->type == MMDB_DATA_TYPE_POINTER) {
+ DEBUG_MSG("pointer points to another pointer");
+ return MMDB_INVALID_DATA_ERROR;
+ }
+
+ /* The pointer could point to any part of the data section but the
+ * next entry for this particular offset may be the one after the
+ * pointer, not the one after whatever the pointer points to. This
+ * depends on whether the pointer points to something that is a simple
+ * value or a compound value. For a compound value, the next one is
+ * the one after the pointer result, not the one after the pointer. */
+ if (entry_data->type != MMDB_DATA_TYPE_MAP
+ && entry_data->type != MMDB_DATA_TYPE_ARRAY) {
+
+ entry_data->offset_to_next = next;
+ }
+ }
+
+ return MMDB_SUCCESS;
+}
+
+#if !MMDB_UINT128_IS_BYTE_ARRAY
+NO_PROTO mmdb_uint128_t get_uint128(const uint8_t *p, int length)
+{
+ mmdb_uint128_t value = 0;
+ while (length-- > 0) {
+ value <<= 8;
+ value += *p++;
+ }
+ return value;
+}
+#endif
+
+LOCAL int decode_one(MMDB_s *mmdb, uint32_t offset,
+ MMDB_entry_data_s *entry_data)
+{
+ const uint8_t *mem = mmdb->data_section;
+
+ if (offset + 1 > mmdb->data_section_size) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+
+ entry_data->offset = offset;
+ entry_data->has_data = true;
+
+ DEBUG_NL;
+ DEBUG_MSGF("Offset: %i", offset);
+
+ uint8_t ctrl = mem[offset++];
+ DEBUG_BINARY("Control byte: %s", ctrl);
+
+ int type = (ctrl >> 5) & 7;
+ DEBUG_MSGF("Type: %i (%s)", type, type_num_to_name(type));
+
+ if (type == MMDB_DATA_TYPE_EXTENDED) {
+ if (offset + 1 > mmdb->data_section_size) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ type = get_ext_type(mem[offset++]);
+ DEBUG_MSGF("Extended type: %i (%s)", type, type_num_to_name(type));
+ }
+
+ entry_data->type = type;
+
+ if (type == MMDB_DATA_TYPE_POINTER) {
+ int psize = (ctrl >> 3) & 3;
+ DEBUG_MSGF("Pointer size: %i", psize);
+
+ if (offset + psize + 1 > mmdb->data_section_size) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ entry_data->pointer = get_ptr_from(ctrl, &mem[offset], psize);
+ DEBUG_MSGF("Pointer to: %i", entry_data->pointer);
+
+ entry_data->data_size = psize + 1;
+ entry_data->offset_to_next = offset + psize + 1;
+ return MMDB_SUCCESS;
+ }
+
+ uint32_t size = ctrl & 31;
+ switch (size) {
+ case 29:
+ if (offset + 1 > mmdb->data_section_size) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ size = 29 + mem[offset++];
+ break;
+ case 30:
+ if (offset + 2 > mmdb->data_section_size) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ size = 285 + get_uint16(&mem[offset]);
+ offset += 2;
+ break;
+ case 31:
+ if (offset + 3 > mmdb->data_section_size) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ size = 65821 + get_uint24(&mem[offset]);
+ offset += 3;
+ default:
+ break;
+ }
+
+ DEBUG_MSGF("Size: %i", size);
+
+ if (type == MMDB_DATA_TYPE_MAP || type == MMDB_DATA_TYPE_ARRAY) {
+ entry_data->data_size = size;
+ entry_data->offset_to_next = offset;
+ return MMDB_SUCCESS;
+ }
+
+ if (type == MMDB_DATA_TYPE_BOOLEAN) {
+ entry_data->boolean = size ? true : false;
+ entry_data->data_size = 0;
+ entry_data->offset_to_next = offset;
+ DEBUG_MSGF("boolean value: %s", entry_data->boolean ? "true" : "false");
+ return MMDB_SUCCESS;
+ }
+
+ // check that the data doesn't extend past the end of the memory
+ // buffer
+ if (offset + size > mmdb->data_section_size) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+
+ if (type == MMDB_DATA_TYPE_UINT16) {
+ if (size > 2) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ entry_data->uint16 = (uint16_t)get_uintX(&mem[offset], size);
+ DEBUG_MSGF("uint16 value: %u", entry_data->uint16);
+ } else if (type == MMDB_DATA_TYPE_UINT32) {
+ if (size > 4) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ entry_data->uint32 = (uint32_t)get_uintX(&mem[offset], size);
+ DEBUG_MSGF("uint32 value: %u", entry_data->uint32);
+ } else if (type == MMDB_DATA_TYPE_INT32) {
+ if (size > 4) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ entry_data->int32 = get_sintX(&mem[offset], size);
+ DEBUG_MSGF("int32 value: %i", entry_data->int32);
+ } else if (type == MMDB_DATA_TYPE_UINT64) {
+ if (size > 8) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ entry_data->uint64 = get_uintX(&mem[offset], size);
+ DEBUG_MSGF("uint64 value: %" PRIu64, entry_data->uint64);
+ } else if (type == MMDB_DATA_TYPE_UINT128) {
+ if (size > 16) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+#if MMDB_UINT128_IS_BYTE_ARRAY
+ memset(entry_data->uint128, 0, 16);
+ if (size > 0) {
+ memcpy(entry_data->uint128 + 16 - size, &mem[offset], size);
+ }
+#else
+ entry_data->uint128 = get_uint128(&mem[offset], size);
+#endif
+ } else if (type == MMDB_DATA_TYPE_FLOAT) {
+ if (size != 4) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ size = 4;
+ entry_data->float_value = get_ieee754_float(&mem[offset]);
+ DEBUG_MSGF("float value: %f", entry_data->float_value);
+ } else if (type == MMDB_DATA_TYPE_DOUBLE) {
+ if (size != 8) {
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ size = 8;
+ entry_data->double_value = get_ieee754_double(&mem[offset]);
+ DEBUG_MSGF("double value: %f", entry_data->double_value);
+ } else if (type == MMDB_DATA_TYPE_UTF8_STRING) {
+ entry_data->utf8_string = size == 0 ? "" : (char *)&mem[offset];
+ entry_data->data_size = size;
+#ifdef MMDB_DEBUG
+ char *string = mmdb_strndup(entry_data->utf8_string,
+ size > 50 ? 50 : size);
+ if (NULL == string) {
+ abort();
+ }
+ DEBUG_MSGF("string value: %s", string);
+ free(string);
+#endif
+ } else if (type == MMDB_DATA_TYPE_BYTES) {
+ entry_data->bytes = &mem[offset];
+ entry_data->data_size = size;
+ }
+
+ entry_data->offset_to_next = offset + size;
+
+ return MMDB_SUCCESS;
+}
+
+LOCAL int get_ext_type(int raw_ext_type)
+{
+ return 7 + raw_ext_type;
+}
+
+LOCAL uint32_t get_ptr_from(uint8_t ctrl, uint8_t const *const ptr,
+ int ptr_size)
+{
+ uint32_t new_offset;
+ switch (ptr_size) {
+ case 0:
+ new_offset = (ctrl & 7) * 256 + ptr[0];
+ break;
+ case 1:
+ new_offset = 2048 + (ctrl & 7) * 65536 + ptr[0] * 256 + ptr[1];
+ break;
+ case 2:
+ new_offset = 2048 + 524288 + (ctrl & 7) * 16777216 + get_uint24(ptr);
+ break;
+ case 3:
+ default:
+ new_offset = get_uint32(ptr);
+ break;
+ }
+ return MMDB_DATA_SECTION_SEPARATOR + new_offset;
+}
+
+int MMDB_get_metadata_as_entry_data_list(
+ MMDB_s *const mmdb, MMDB_entry_data_list_s **const entry_data_list)
+{
+ MMDB_s metadata_db = make_fake_metadata_db(mmdb);
+
+ MMDB_entry_s metadata_start = {
+ .mmdb = &metadata_db,
+ .offset = 0
+ };
+
+ return MMDB_get_entry_data_list(&metadata_start, entry_data_list);
+}
+
+int MMDB_get_entry_data_list(
+ MMDB_entry_s *start, MMDB_entry_data_list_s **const entry_data_list)
+{
+ *entry_data_list = new_entry_data_list();
+ if (NULL == *entry_data_list) {
+ return MMDB_OUT_OF_MEMORY_ERROR;
+ }
+ return get_entry_data_list(start->mmdb, start->offset, *entry_data_list, 0);
+}
+
+LOCAL int get_entry_data_list(MMDB_s *mmdb, uint32_t offset,
+ MMDB_entry_data_list_s *const entry_data_list,
+ int depth)
+{
+ if (depth >= MAXIMUM_DATA_STRUCTURE_DEPTH) {
+ DEBUG_MSG("reached the maximum data structure depth");
+ return MMDB_INVALID_DATA_ERROR;
+ }
+ depth++;
+ CHECKED_DECODE_ONE(mmdb, offset, &entry_data_list->entry_data);
+
+ switch (entry_data_list->entry_data.type) {
+ case MMDB_DATA_TYPE_POINTER:
+ {
+ uint32_t next_offset = entry_data_list->entry_data.offset_to_next;
+ uint32_t last_offset;
+ CHECKED_DECODE_ONE(mmdb, last_offset =
+ entry_data_list->entry_data.pointer,
+ &entry_data_list->entry_data);
+
+ /* Pointers to pointers are illegal under the spec */
+ if (entry_data_list->entry_data.type == MMDB_DATA_TYPE_POINTER) {
+ DEBUG_MSG("pointer points to another pointer");
+ return MMDB_INVALID_DATA_ERROR;
+ }
+
+ if (entry_data_list->entry_data.type == MMDB_DATA_TYPE_ARRAY
+ || entry_data_list->entry_data.type == MMDB_DATA_TYPE_MAP) {
+
+ int status =
+ get_entry_data_list(mmdb, last_offset, entry_data_list,
+ depth);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+ }
+ entry_data_list->entry_data.offset_to_next = next_offset;
+ }
+ break;
+ case MMDB_DATA_TYPE_ARRAY:
+ {
+ uint32_t array_size = entry_data_list->entry_data.data_size;
+ uint32_t array_offset = entry_data_list->entry_data.offset_to_next;
+ MMDB_entry_data_list_s *previous = entry_data_list;
+ while (array_size-- > 0) {
+ MMDB_entry_data_list_s *entry_data_list_to = previous->next =
+ new_entry_data_list();
+ if (NULL == entry_data_list_to) {
+ return MMDB_OUT_OF_MEMORY_ERROR;
+ }
+
+ int status =
+ get_entry_data_list(mmdb, array_offset, entry_data_list_to,
+ depth);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+
+ array_offset = entry_data_list_to->entry_data.offset_to_next;
+ while (previous->next) {
+ previous = previous->next;
+ }
+ }
+ entry_data_list->entry_data.offset_to_next = array_offset;
+
+ }
+ break;
+ case MMDB_DATA_TYPE_MAP:
+ {
+ uint32_t size = entry_data_list->entry_data.data_size;
+
+ offset = entry_data_list->entry_data.offset_to_next;
+ MMDB_entry_data_list_s *previous = entry_data_list;
+ while (size-- > 0) {
+ MMDB_entry_data_list_s *entry_data_list_to = previous->next =
+ new_entry_data_list();
+ if (NULL == entry_data_list_to) {
+ return MMDB_OUT_OF_MEMORY_ERROR;
+ }
+
+ int status =
+ get_entry_data_list(mmdb, offset, entry_data_list_to,
+ depth);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+
+ while (previous->next) {
+ previous = previous->next;
+ }
+
+ offset = entry_data_list_to->entry_data.offset_to_next;
+ entry_data_list_to = previous->next =
+ new_entry_data_list();
+
+ if (NULL == entry_data_list_to) {
+ return MMDB_OUT_OF_MEMORY_ERROR;
+ }
+
+ status = get_entry_data_list(mmdb, offset, entry_data_list_to,
+ depth);
+ if (MMDB_SUCCESS != status) {
+ return status;
+ }
+
+ while (previous->next) {
+ previous = previous->next;
+ }
+ offset = entry_data_list_to->entry_data.offset_to_next;
+ }
+ entry_data_list->entry_data.offset_to_next = offset;
+ }
+ break;
+ default:
+ break;
+ }
+
+ return MMDB_SUCCESS;
+}
+
+LOCAL float get_ieee754_float(const uint8_t *restrict p)
+{
+ volatile float f;
+ uint8_t *q = (void *)&f;
+#if __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
+ q[3] = p[0];
+ q[2] = p[1];
+ q[1] = p[2];
+ q[0] = p[3];
+#else
+ memcpy(q, p, 4);
+#endif
+ return f;
+}
+
+LOCAL double get_ieee754_double(const uint8_t *restrict p)
+{
+ volatile double d;
+ uint8_t *q = (void *)&d;
+#if __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
+ q[7] = p[0];
+ q[6] = p[1];
+ q[5] = p[2];
+ q[4] = p[3];
+ q[3] = p[4];
+ q[2] = p[5];
+ q[1] = p[6];
+ q[0] = p[7];
+#else
+ memcpy(q, p, 8);
+#endif
+
+ return d;
+}
+
+LOCAL uint32_t get_uint32(const uint8_t *p)
+{
+ return p[0] * 16777216U + p[1] * 65536 + p[2] * 256 + p[3];
+}
+
+LOCAL uint32_t get_uint24(const uint8_t *p)
+{
+ return p[0] * 65536U + p[1] * 256 + p[2];
+}
+
+LOCAL uint32_t get_uint16(const uint8_t *p)
+{
+ return p[0] * 256U + p[1];
+}
+
+LOCAL uint64_t get_uintX(const uint8_t *p, int length)
+{
+ uint64_t value = 0;
+ while (length-- > 0) {
+ value <<= 8;
+ value += *p++;
+ }
+ return value;
+}
+
+LOCAL int32_t get_sintX(const uint8_t *p, int length)
+{
+ return (int32_t)get_uintX(p, length);
+}
+
+LOCAL MMDB_entry_data_list_s *new_entry_data_list(void)
+{
+ /* We need calloc here in order to ensure that the ->next pointer in the
+ * struct doesn't point to some random address. */
+ return calloc(1, sizeof(MMDB_entry_data_list_s));
+}
+
+void MMDB_free_entry_data_list(MMDB_entry_data_list_s *const entry_data_list)
+{
+ if (entry_data_list == NULL) {
+ return;
+ }
+ if (entry_data_list->next) {
+ MMDB_free_entry_data_list(entry_data_list->next);
+ }
+ free(entry_data_list);
+}
+
+void MMDB_close(MMDB_s *const mmdb)
+{
+ free_mmdb_struct(mmdb);
+}
+
+LOCAL void free_mmdb_struct(MMDB_s *const mmdb)
+{
+ if (!mmdb) {
+ return;
+ }
+
+ if (NULL != mmdb->filename) {
+ FREE_AND_SET_NULL(mmdb->filename);
+ }
+ if (NULL != mmdb->file_content) {
+#ifdef _WIN32
+ UnmapViewOfFile(mmdb->file_content);
+ /* Winsock is only initialized if open was successful so we only have
+ * to cleanup then. */
+ WSACleanup();
+#else
+ munmap((void *)mmdb->file_content, mmdb->file_size);
+#endif
+ }
+
+ if (NULL != mmdb->metadata.database_type) {
+ FREE_AND_SET_NULL(mmdb->metadata.database_type);
+ }
+
+ free_languages_metadata(mmdb);
+ free_descriptions_metadata(mmdb);
+}
+
+LOCAL void free_languages_metadata(MMDB_s *mmdb)
+{
+ if (!mmdb->metadata.languages.count) {
+ return;
+ }
+
+ for (size_t i = 0; i < mmdb->metadata.languages.count; i++) {
+ FREE_AND_SET_NULL(mmdb->metadata.languages.names[i]);
+ }
+ FREE_AND_SET_NULL(mmdb->metadata.languages.names);
+}
+
+LOCAL void free_descriptions_metadata(MMDB_s *mmdb)
+{
+ if (!mmdb->metadata.description.count) {
+ return;
+ }
+
+ for (size_t i = 0; i < mmdb->metadata.description.count; i++) {
+ if (NULL != mmdb->metadata.description.descriptions[i]) {
+ if (NULL !=
+ mmdb->metadata.description.descriptions[i]->language) {
+ FREE_AND_SET_NULL(
+ mmdb->metadata.description.descriptions[i]->language);
+ }
+
+ if (NULL !=
+ mmdb->metadata.description.descriptions[i]->description) {
+ FREE_AND_SET_NULL(
+ mmdb->metadata.description.descriptions[i]->description);
+ }
+ FREE_AND_SET_NULL(mmdb->metadata.description.descriptions[i]);
+ }
+ }
+
+ FREE_AND_SET_NULL(mmdb->metadata.description.descriptions);
+}
+
+const char *MMDB_lib_version(void)
+{
+ return PACKAGE_VERSION;
+}
+
+int MMDB_dump_entry_data_list(FILE *const stream,
+ MMDB_entry_data_list_s *const entry_data_list,
+ int indent)
+{
+ int status;
+ dump_entry_data_list(stream, entry_data_list, indent, &status);
+ return status;
+}
+
+LOCAL MMDB_entry_data_list_s *dump_entry_data_list(
+ FILE *stream, MMDB_entry_data_list_s *entry_data_list, int indent,
+ int *status)
+{
+ switch (entry_data_list->entry_data.type) {
+ case MMDB_DATA_TYPE_MAP:
+ {
+ uint32_t size = entry_data_list->entry_data.data_size;
+
+ print_indentation(stream, indent);
+ fprintf(stream, "{\n");
+ indent += 2;
+
+ for (entry_data_list = entry_data_list->next;
+ size && entry_data_list; size--) {
+
+ if (MMDB_DATA_TYPE_UTF8_STRING !=
+ entry_data_list->entry_data.type) {
+ *status = MMDB_INVALID_DATA_ERROR;
+ return NULL;
+ }
+ char *key =
+ mmdb_strndup(
+ (char *)entry_data_list->entry_data.utf8_string,
+ entry_data_list->entry_data.data_size);
+ if (NULL == key) {
+ *status = MMDB_OUT_OF_MEMORY_ERROR;
+ return NULL;
+ }
+
+ print_indentation(stream, indent);
+ fprintf(stream, "\"%s\": \n", key);
+ free(key);
+
+ entry_data_list = entry_data_list->next;
+ entry_data_list =
+ dump_entry_data_list(stream, entry_data_list, indent + 2,
+ status);
+
+ if (MMDB_SUCCESS != *status) {
+ return NULL;
+ }
+ }
+
+ indent -= 2;
+ print_indentation(stream, indent);
+ fprintf(stream, "}\n");
+ }
+ break;
+ case MMDB_DATA_TYPE_ARRAY:
+ {
+ uint32_t size = entry_data_list->entry_data.data_size;
+
+ print_indentation(stream, indent);
+ fprintf(stream, "[\n");
+ indent += 2;
+
+ for (entry_data_list = entry_data_list->next;
+ size && entry_data_list; size--) {
+ entry_data_list =
+ dump_entry_data_list(stream, entry_data_list, indent,
+ status);
+ if (MMDB_SUCCESS != *status) {
+ return NULL;
+ }
+ }
+
+ indent -= 2;
+ print_indentation(stream, indent);
+ fprintf(stream, "]\n");
+ }
+ break;
+ case MMDB_DATA_TYPE_UTF8_STRING:
+ {
+ char *string =
+ mmdb_strndup((char *)entry_data_list->entry_data.utf8_string,
+ entry_data_list->entry_data.data_size);
+ if (NULL == string) {
+ *status = MMDB_OUT_OF_MEMORY_ERROR;
+ return NULL;
+ }
+ print_indentation(stream, indent);
+ fprintf(stream, "\"%s\" \n", string);
+ free(string);
+ entry_data_list = entry_data_list->next;
+ }
+ break;
+ case MMDB_DATA_TYPE_BYTES:
+ {
+ char *hex_string =
+ bytes_to_hex((uint8_t *)entry_data_list->entry_data.bytes,
+ entry_data_list->entry_data.data_size);
+ if (NULL == hex_string) {
+ *status = MMDB_OUT_OF_MEMORY_ERROR;
+ return NULL;
+ }
+
+ print_indentation(stream, indent);
+ fprintf(stream, "%s \n", hex_string);
+ free(hex_string);
+
+ entry_data_list = entry_data_list->next;
+ }
+ break;
+ case MMDB_DATA_TYPE_DOUBLE:
+ print_indentation(stream, indent);
+ fprintf(stream, "%f \n",
+ entry_data_list->entry_data.double_value);
+ entry_data_list = entry_data_list->next;
+ break;
+ case MMDB_DATA_TYPE_FLOAT:
+ print_indentation(stream, indent);
+ fprintf(stream, "%f \n",
+ entry_data_list->entry_data.float_value);
+ entry_data_list = entry_data_list->next;
+ break;
+ case MMDB_DATA_TYPE_UINT16:
+ print_indentation(stream, indent);
+ fprintf(stream, "%u \n", entry_data_list->entry_data.uint16);
+ entry_data_list = entry_data_list->next;
+ break;
+ case MMDB_DATA_TYPE_UINT32:
+ print_indentation(stream, indent);
+ fprintf(stream, "%u \n", entry_data_list->entry_data.uint32);
+ entry_data_list = entry_data_list->next;
+ break;
+ case MMDB_DATA_TYPE_BOOLEAN:
+ print_indentation(stream, indent);
+ fprintf(stream, "%s \n",
+ entry_data_list->entry_data.boolean ? "true" : "false");
+ entry_data_list = entry_data_list->next;
+ break;
+ case MMDB_DATA_TYPE_UINT64:
+ print_indentation(stream, indent);
+ fprintf(stream, "%" PRIu64 " \n",
+ entry_data_list->entry_data.uint64);
+ entry_data_list = entry_data_list->next;
+ break;
+ case MMDB_DATA_TYPE_UINT128:
+ print_indentation(stream, indent);
+#if MMDB_UINT128_IS_BYTE_ARRAY
+ char *hex_string =
+ bytes_to_hex((uint8_t *)entry_data_list->entry_data.uint128, 16);
+ if (NULL == hex_string) {
+ *status = MMDB_OUT_OF_MEMORY_ERROR;
+ return NULL;
+ }
+ fprintf(stream, "0x%s \n", hex_string);
+ free(hex_string);
+#else
+ uint64_t high = entry_data_list->entry_data.uint128 >> 64;
+ uint64_t low = (uint64_t)entry_data_list->entry_data.uint128;
+ fprintf(stream, "0x%016" PRIX64 "%016" PRIX64 " \n", high,
+ low);
+#endif
+ entry_data_list = entry_data_list->next;
+ break;
+ case MMDB_DATA_TYPE_INT32:
+ print_indentation(stream, indent);
+ fprintf(stream, "%d \n", entry_data_list->entry_data.int32);
+ entry_data_list = entry_data_list->next;
+ break;
+ default:
+ *status = MMDB_INVALID_DATA_ERROR;
+ return NULL;
+ }
+
+ *status = MMDB_SUCCESS;
+ return entry_data_list;
+}
+
+LOCAL void print_indentation(FILE *stream, int i)
+{
+ char buffer[1024];
+ int size = i >= 1024 ? 1023 : i;
+ memset(buffer, 32, size);
+ buffer[size] = '\0';
+ fputs(buffer, stream);
+}
+
+LOCAL char *bytes_to_hex(uint8_t *bytes, uint32_t size)
+{
+ char *hex_string;
+ MAYBE_CHECK_SIZE_OVERFLOW(size, SIZE_MAX / 2 - 1, NULL);
+
+ hex_string = malloc((size * 2) + 1);
+ if (NULL == hex_string) {
+ return NULL;
+ }
+
+ for (uint32_t i = 0; i < size; i++) {
+ sprintf(hex_string + (2 * i), "%02X", bytes[i]);
+ }
+
+ return hex_string;
+}
+
+const char *MMDB_strerror(int error_code)
+{
+ switch (error_code) {
+ case MMDB_SUCCESS:
+ return "Success (not an error)";
+ case MMDB_FILE_OPEN_ERROR:
+ return "Error opening the specified MaxMind DB file";
+ case MMDB_CORRUPT_SEARCH_TREE_ERROR:
+ return "The MaxMind DB file's search tree is corrupt";
+ case MMDB_INVALID_METADATA_ERROR:
+ return "The MaxMind DB file contains invalid metadata";
+ case MMDB_IO_ERROR:
+ return "An attempt to read data from the MaxMind DB file failed";
+ case MMDB_OUT_OF_MEMORY_ERROR:
+ return "A memory allocation call failed";
+ case MMDB_UNKNOWN_DATABASE_FORMAT_ERROR:
+ return
+ "The MaxMind DB file is in a format this library can't handle (unknown record size or binary format version)";
+ case MMDB_INVALID_DATA_ERROR:
+ return
+ "The MaxMind DB file's data section contains bad data (unknown data type or corrupt data)";
+ case MMDB_INVALID_LOOKUP_PATH_ERROR:
+ return
+ "The lookup path contained an invalid value (like a negative integer for an array index)";
+ case MMDB_LOOKUP_PATH_DOES_NOT_MATCH_DATA_ERROR:
+ return
+ "The lookup path does not match the data (key that doesn't exist, array index bigger than the array, expected array or map where none exists)";
+ case MMDB_INVALID_NODE_NUMBER_ERROR:
+ return
+ "The MMDB_read_node function was called with a node number that does not exist in the search tree";
+ case MMDB_IPV6_LOOKUP_IN_IPV4_DATABASE_ERROR:
+ return
+ "You attempted to look up an IPv6 address in an IPv4-only database";
+ default:
+ return "Unknown error code";
+ }
+}
diff --git a/external/RandomLib/LICENSE.txt b/external/RandomLib/LICENSE.txt
deleted file mode 100644
index 78f3414c..00000000
--- a/external/RandomLib/LICENSE.txt
+++ /dev/null
@@ -1,23 +0,0 @@
-This license applies to RandomLib, versions 1.2 and later.
-
-Copyright (c) 2006-2013, Charles Karney
-
-Permission is hereby granted, free of charge, to any person
-obtaining a copy of this software and associated documentation
-files (the "Software"), to deal in the Software without
-restriction, including without limitation the rights to use, copy,
-modify, merge, publish, distribute, sublicense, and/or sell copies
-of the Software, and to permit persons to whom the Software is
-furnished to do so, subject to the following conditions:
-
-The above copyright notice and this permission notice shall be
-included in all copies or substantial portions of the Software.
-
-THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
-EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
-MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
-NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
-HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
-WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
-OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
-DEALINGS IN THE SOFTWARE.
diff --git a/external/RandomLib/Random.cpp b/external/RandomLib/Random.cpp
deleted file mode 100644
index 1b752095..00000000
--- a/external/RandomLib/Random.cpp
+++ /dev/null
@@ -1,1405 +0,0 @@
-/**
- * \file Random.cpp
- * \brief Implementation code for %RandomLib
- *
- * Copyright (c) Charles Karney (2006-2011) and licensed
- * under the MIT/X11 License. For more information, see
- * http://randomlib.sourceforge.net/
- *
- * \brief Code for MixerMT0, MixerMT1, MixerSFMT.
- *
- * MixerMT0 is adapted from MT19937 (init_by_array) and MT19937_64
- * (init_by_array64) by Makoto Matsumoto and Takuji Nishimura. See
- * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c and
- * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/VERSIONS/C-LANG/mt19937-64.c
- *
- * MixerMT1 contains modifications to MixerMT0 by Charles Karney to
- * correct defects in MixerMT0. This is described in W. E. Brown,
- * M. Fischler, J. Kowalkowski, M. Paterno, Random Number Generation in C++0X:
- * A Comprehensive Proposal, version 3, Sept 2006, Sec. 26.4.7.1,
- * http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2006/n2079.pdf
- * This has been replaced in the C++11 standard by MixerSFMT.
- *
- * MixerSFMT is adapted from SFMT19937's init_by_array Mutsuo Saito given in
- * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/SFMT-src-1.2.tar.gz and
- * is part of the C++11 standard; see P. Becker, Working Draft, Standard for
- * Programming Language C++, Oct. 2007, Sec. 26.4.7.1,
- * http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2007/n2461.pdf
- *
- * The adaption to the C++ is copyright (c) Charles Karney (2006-2011)
- * and licensed under the MIT/X11 License. For more
- * information, see http://randomlib.sourceforge.net/
- *
- * \brief Code for MT19937 and SFMT19937.
- *
- * MT19937 is adapted from MT19937 and MT19937_64 by Makoto Matsumoto and
- * Takuji Nishimura. See
- * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c and
- * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/VERSIONS/C-LANG/mt19937-64.c
- *
- * The code for stepping MT19937 backwards is adapted (and simplified) from
- * revrand() by Katsumi Hagita. See
- * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/VERSIONS/FORTRAN/REVmt19937b.f
- *
- * SFMT19937 is adapted from SFMT19937 Mutsuo Saito given in
- * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/M062821.pdf and
- * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/SFMT-src-1.2.tar.gz
- *
- * The code for stepping SFMT19937 backwards is by Charles Karney.
- *
- * The adaption to the C++ is copyright (c) Charles Karney (2006-2011)
- * and licensed under the MIT/X11 License. For more
- * information, see http://randomlib.sourceforge.net/
- **********************************************************************/
-
-#define RANDOMLIB_RANDOM_CPP 1
-
-/**
- * Let the header file know that the library is being built.
- **********************************************************************/
-#define RANDOMLIB_BUILDING_LIBRARY 1
-
-#include
-
-#if defined(_MSC_VER) || defined(__MINGW32__)
-#define RANDOMLIB_WINDOWS 1
-#else
-#define RANDOMLIB_WINDOWS 0
-#endif
-
-#include // For SeedWord reading /dev/urandom
-#include // For SeedWord calling time()
-#include // For formatting in Write32/Read32
-#include // For formatting in Write32/Read32
-#if !RANDOMLIB_WINDOWS
-#include // For SeedWord calling gettimeofday
-#include // For SeedWord calling getpid(), gethostid()
-#else
-#include // For SeedWord calling high prec timer
-#include
-#include // For SeedWord calling getpid()
-#define getpid _getpid
-#define gmtime_r(t,g) gmtime_s(g,t)
-#endif
-
-#if RANDOMLIB_WINDOWS || defined(__CYGWIN__)
-#define strtoull strtoul
-#endif
-
-#if defined(_MSC_VER)
-// Squelch warnings about constant conditional expressions
-# pragma warning (disable: 4127)
-#endif
-
-namespace RandomLib {
-
- // RandomType implementation
-
- template<>
- void Random_u32::Write32(std::ostream& os, bool bin, int& cnt,
- Random_u32::type x) {
- if (bin) {
- unsigned char buf[4];
- // Use network order -- most significant byte first
- buf[3] = (unsigned char)(x);
- buf[2] = (unsigned char)(x >>= 8);
- buf[1] = (unsigned char)(x >>= 8);
- buf[0] = (unsigned char)(x >>= 8);
- os.write(reinterpret_cast(buf), 4);
- } else {
- const int longsperline = 72/9;
- // Use hexadecimal to minimize storage together with stringstream to
- // isolate the effect of changing the base.
- std::ostringstream str;
- // No spacing before or after
- if (cnt > 0)
- // Newline every longsperline longs
- str << (cnt % longsperline ? ' ' : '\n');
- str << std::hex << x;
- os << str.str();
- ++cnt;
- }
- }
-
- template<>
- void Random_u32::Read32(std::istream& is, bool bin, Random_u32::type& x) {
- if (bin) {
- unsigned char buf[4];
- is.read(reinterpret_cast(buf), 4);
- // Use network order -- most significant byte first
- x = Random_u32::type(buf[0]) << 24 | Random_u32::type(buf[1]) << 16 |
- Random_u32::type(buf[2]) << 8 | Random_u32::type(buf[3]);
- } else {
- std::string s;
- is >> std::ws >> s;
- // Use hexadecimal to minimize storage together with stringstream to
- // isolate the effect of changing the base.
- std::istringstream str(s);
- str >> std::hex >> x;
- }
- x &= Random_u32::mask;
- }
-
- template<>
- void Random_u64::Write32(std::ostream& os, bool bin, int& cnt,
- Random_u64::type x) {
- Random_u32::Write32(os, bin, cnt, Random_u32::cast(x >> 32));
- Random_u32::Write32(os, bin, cnt, Random_u32::cast(x ));
- }
-
- template<>
- void Random_u64::Read32(std::istream& is, bool bin, Random_u64::type& x) {
- Random_u32::type t;
- Random_u32::Read32(is, bin, t);
- x = Random_u64::type(t) << 32;
- Random_u32::Read32(is, bin, t);
- x |= Random_u64::type(t);
- }
-
- // RandomSeed implementation
-
- RandomSeed::seed_type RandomSeed::SeedWord() {
- // Check that the assumptions made about the capabilities of the number
- // system are valid.
- STATIC_ASSERT(std::numeric_limits::radix == 2 &&
- !std::numeric_limits::is_signed &&
- std::numeric_limits::digits >= 32,
- "seed_type is a bad type");
- u32::type t = 0;
- // Linux has /dev/urandom to initialize the seed randomly. (Use
- // /dev/urandom instead of /dev/random because it does not block.)
- {
- std::ifstream f("/dev/urandom", std::ios::binary | std::ios::in);
- if (f.good()) {
- // Read 32 bits from /dev/urandom
- f.read(reinterpret_cast(&t), sizeof(t));
- }
- }
- std::vector v = SeedVector();
- for (size_t i = v.size(); i--;)
- u32::CheckSum(u32::type(v[i]), t);
- return seed_t::cast(t);
- }
-
- std::vector RandomSeed::SeedVector() {
- std::vector v;
- {
- // fine-grained timer
-#if !RANDOMLIB_WINDOWS
- timeval tv;
- if (gettimeofday(&tv, 0) == 0)
- v.push_back(seed_t::cast(tv.tv_usec));
-#else
- LARGE_INTEGER taux;
- if (QueryPerformanceCounter((LARGE_INTEGER *)&taux)) {
- v.push_back(seed_t::cast(taux.LowPart));
- v.push_back(seed_t::cast(taux.HighPart));
- }
-#endif
- }
- // seconds
- const time_t tim = std::time(0);
- v.push_back(seed_t::cast(seed_type(tim)));
- // PID
- v.push_back(seed_t::cast(getpid()));
-#if !RANDOMLIB_WINDOWS
- // host ID
- v.push_back(seed_t::cast(gethostid()));
-#endif
- {
- // year
-#if !defined(__MINGW32__)
- tm gt;
- gmtime_r(&tim, >);
- v.push_back((seed_type(1900) + seed_t::cast(gt.tm_year)));
-#else
- tm* gt = gmtime(&tim);
- v.push_back((seed_type(1900) + seed_t::cast(gt->tm_year)));
-#endif
- }
- // Candidates for additional elements:
- // ip address(es) of computer, thread index.
- std::transform(v.begin(), v.end(), v.begin(), seed_t::cast);
- return v;
- }
-
- std::vector
- RandomSeed::StringToVector(const std::string& s) {
- std::vector v(0);
- const char* c = s.c_str();
- char* q;
- std::string::size_type p = 0;
- while (true) {
- p = s.find_first_of("0123456789", p);
- if (p == std::string::npos)
- break;
- v.push_back(seed_t::cast(std::strtoull(c + p, &q, 0)));
- p = q - c;
- }
- return v;
- }
-
- // RandomEngine implementation
-
- template
- void RandomEngine::Init() throw() {
- // On exit we have _ptr == N.
-
- STATIC_ASSERT(std::numeric_limits::radix == 2 &&
- !std::numeric_limits::is_signed &&
- std::numeric_limits::digits >=
- int(mixer_t::width),
- "mixer_type is a bad type");
-
- STATIC_ASSERT(std::numeric_limits::radix == 2 &&
- !std::numeric_limits::is_signed &&
- std::numeric_limits::digits >= width,
- "engine_type is a bad type");
-
- STATIC_ASSERT(mixer_t::width == 32 || mixer_t::width == 64,
- "Mixer width must be 32 or 64");
-
- STATIC_ASSERT(width == 32 || width == 64,
- "Algorithm width must be 32 or 64");
-
- // If the bit-widths are the same then the data sizes must be the same.
- STATIC_ASSERT(!(mixer_t::width == width) ||
- sizeof(_stateu) == sizeof(_state),
- "Same bit-widths but different storage");
-
- // Repacking assumes that narrower data type is at least as wasteful than
- // the broader one.
- STATIC_ASSERT(!(mixer_t::width < width) ||
- sizeof(_stateu) >= sizeof(_state),
- "Narrow data type uses less storage");
-
- STATIC_ASSERT(!(mixer_t::width > width) ||
- sizeof(_stateu) <= sizeof(_state),
- "Narrow data type uses less storage");
-
- // Require that _statev and _state are aligned since no repacking is done
- // when calling Transition
- STATIC_ASSERT(sizeof(_statev) == sizeof(_state),
- "Storage mismatch with internal engine data type");
-
- // Convert the seed into state
- Mixer::SeedToState(_seed, _stateu, NU);
-
- // Pack into _state
- if (mixer_t::width < width) {
- for (size_t i = 0; i < N; ++i)
- // Assume 2:1 LSB packing
- _state[i] = result_type(_stateu[2*i]) |
- result_type(_stateu[2*i + 1]) <<
- (mixer_t::width < width ? mixer_t::width : 0);
- } else if (mixer_t::width > width) {
- for (size_t i = N; i--;)
- // Assume 1:2 LSB packing
- _state[i] = result_t::cast(_stateu[i>>1] >> width * (i&1u));
- } // Otherwise the union takes care of it
-
- Algorithm::NormalizeState(_state);
-
- _rounds = -1;
- _ptr = N;
- }
-
- template Random_u32::type
- RandomEngine::Check(u64::type v, u32::type e,
- u32::type m) const {
- if (v != version)
- throw RandomErr(Name() + ": Unknown version");
- if (e != Algorithm::version)
- throw RandomErr(Name() + ": Algorithm mismatch");
- if (m != Mixer::version)
- throw RandomErr(Name() + ": Mixer mismatch");
- u32::type check = 0;
- u64::CheckSum(v, check);
- u32::CheckSum(e, check);
- u32::CheckSum(m, check);
- u32::CheckSum(u32::type(_seed.size()), check);
- for (std::vector::const_iterator n = _seed.begin();
- n != _seed.end(); ++n) {
- if (*n != seed_t::cast(*n))
- throw RandomErr(Name() + ": Illegal seed value");
- u32::CheckSum(u32::type(*n), check);
- }
- u32::CheckSum(_ptr, check);
- if (_stride == 0 || _stride > UNINIT/2)
- throw RandomErr(Name() + ": Invalid stride");
- u32::CheckSum(_stride, check);
- if (_ptr != UNINIT) {
- if (_ptr >= N + _stride)
- throw RandomErr(Name() + ": Invalid pointer");
- u64::CheckSum(_rounds, check);
- Algorithm::CheckState(_state, check);
- }
- return check;
- }
-
- template
- RandomEngine::RandomEngine(std::istream& is, bool bin) {
- u64::type versionr;
- u32::type versione, versionm, t;
- u64::Read32(is, bin, versionr);
- u32::Read32(is, bin, versione);
- u32::Read32(is, bin, versionm);
- u32::Read32(is, bin, t);
- _seed.resize(size_t(t));
- for (std::vector::iterator n = _seed.begin();
- n != _seed.end(); ++n) {
- u32::Read32(is, bin, t);
- *n = seed_type(t);
- }
- u32::Read32(is, bin, t);
- // Don't need to worry about sign extension because _ptr is unsigned.
- _ptr = unsigned(t);
- u32::Read32(is, bin, t);
- _stride = unsigned(t);
- if (_ptr != UNINIT) {
- u64::type p;
- u64::Read32(is, bin, p);
- _rounds = (long long)(p);
- // Sign extension in case long long is bigger than 64 bits.
- _rounds <<= 63 - std::numeric_limits::digits;
- _rounds >>= 63 - std::numeric_limits::digits;
- for (unsigned i = 0; i < N; ++i)
- result_t::Read32(is, bin, _state[i]);
- }
- u32::Read32(is, bin, t);
- if (t != Check(versionr, versione, versionm))
- throw RandomErr(Name() + ": Checksum failure");
- }
-
- template
- void RandomEngine::Save(std::ostream& os,
- bool bin) const {
- u32::type check = Check(version, Algorithm::version, Mixer::version);
- int c = 0;
- u64::Write32(os, bin, c, version);
- u32::Write32(os, bin, c, Algorithm::version);
- u32::Write32(os, bin, c, Mixer::version);
- u32::Write32(os, bin, c, u32::type(_seed.size()));
- for (std::vector::const_iterator n = _seed.begin();
- n != _seed.end(); ++n)
- u32::Write32(os, bin, c, u32::type(*n));
- u32::Write32(os, bin, c, _ptr);
- u32::Write32(os, bin, c, _stride);
- if (_ptr != UNINIT) {
- u64::Write32(os, bin, c, u64::type(_rounds));
- for (unsigned i = 0; i < N; ++i)
- result_t::Write32(os, bin, c, _state[i]);
- }
- u32::Write32(os, bin, c, check);
- }
-
- template
- void RandomEngine::StepCount(long long n) throw() {
- // On exit we have 0 <= _ptr <= N.
- if (_ptr == UNINIT)
- Init();
- const long long ncount = n + Count(); // new Count()
- long long nrounds = ncount / N;
- int nptr = int(ncount - nrounds * N);
- // We pick _ptr = N or _ptr = 0 depending on which choice involves the
- // least work. We thus avoid doing one (potentially unneeded) call to
- // Transition.
- if (nptr < 0) {
- --nrounds;
- nptr += N;
- } else if (nptr == 0 && nrounds > _rounds) {
- nptr = N;
- --nrounds;
- }
- if (nrounds != _rounds)
- Algorithm::Transition(nrounds - _rounds, _statev);
- _rounds = nrounds;
- _ptr = nptr;
- }
-
- template
- void RandomEngine::SelfTest() {
- RandomEngine g(std::vector(0));
- g.SetCount(10000-1);
- result_type x = g();
- if (SelfTestResult(0) && x != SelfTestResult(1))
- throw RandomErr(Name() + ": Incorrect result with seed " +
- g.SeedString());
- seed_type s[] = {0x1234U, 0x5678U, 0x9abcU, 0xdef0U};
- // seed_type s[] = {1, 2, 3, 4};
- g.Reseed(s, s+4);
- g.StepCount(-20000);
- std::string save;
- {
- std::ostringstream stream;
- stream << g << "\n";
- save = stream.str();
- }
- g.Reset();
- {
- std::istringstream stream(save);
- stream >> g;
- }
- g.SetCount(10000);
- {
- std::ostringstream stream;
- g.Save(stream, true);
- save = stream.str();
- }
- {
- std::istringstream stream(save);
- RandomEngine h(std::vector(0));
- h.Load(stream, true);
- h.SetCount(1000000-1);
- x = h();
- if (SelfTestResult(0) && x != SelfTestResult(2))
- throw RandomErr(Name() + ": Incorrect result with seed " +
- h.SeedString());
- g.SetCount(1000000);
- if (h != g)
- throw RandomErr(Name() + ": Comparison failure");
- }
- }
-
- template<> Random_u32::type
- RandomEngine, MixerMT0 >::
- SelfTestResult(unsigned i) throw() {
- return i == 0 ? 1 :
- i == 1 ? 4123659995UL : 3016432305UL;
- }
-
- template<> Random_u64::type
- RandomEngine, MixerMT0 >::
- SelfTestResult(unsigned i) throw() {
- return i == 0 ? 1 :
- i == 1 ? 9981545732273789042ULL : 1384037754719008581ULL;
- }
-
- template<> Random_u32::type
- RandomEngine, MixerMT1 >::
- SelfTestResult(unsigned i) throw() {
- return i == 0 ? 1 :
- i == 1 ? 4123659995UL : 2924523180UL;
- }
-
- template<> Random_u64::type
- RandomEngine, MixerMT1 >::
- SelfTestResult(unsigned i) throw() {
- return i == 0 ? 1 :
- i == 1 ? 9981545732273789042ULL : 5481486777409645478ULL;
- }
-
- template<> Random_u32::type
- RandomEngine, MixerSFMT>::
- SelfTestResult(unsigned i) throw() {
- return i == 0 ? 1 :
- i == 1 ? 666528879UL : 2183745132UL;
- }
-
- template<> Random_u64::type
- RandomEngine, MixerSFMT>::
- SelfTestResult(unsigned i) throw() {
- return i == 0 ? 1 :
- i == 1 ? 12176471137395770412ULL : 66914054428611861ULL;
- }
-
- template<> Random_u32::type
- RandomEngine, MixerSFMT>::
- SelfTestResult(unsigned i) throw() {
- return i == 0 ? 1 :
- i == 1 ? 2695024307UL : 782200760UL;
- }
-
- template<> Random_u64::type
- RandomEngine, MixerSFMT>::
- SelfTestResult(unsigned i) throw() {
- return i == 0 ? 1 :
- i == 1 ? 1464461649847485149ULL : 5050640804923595109ULL;
- }
-
- // RandomMixer implementation
-
- template void MixerMT0::
- SeedToState(const std::vector& seed,
- mixer_type state[], unsigned n) throw() {
- // Adapted from
- // http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c
- // http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/VERSIONS/C-LANG/mt19937-64.c
- const unsigned s = unsigned(seed.size());
- const unsigned w = mixer_t::width;
-
- mixer_type r = s ? a1 + mixer_type(0) : a0 + mixer_type(0);
- state[0] = r;
- for (unsigned k = 1; k < n; ++k) {
- r = b * (r ^ r >> (w - 2)) + k;
- r &= mask;
- state[k] = r;
- }
- if (s > 0) {
- const unsigned m = mixer_t::width / 32,
- s2 = (s + m - 1)/m;
- unsigned i1 = 1;
- r = state[0];
- for (unsigned k = (n > s2 ? n : s2), j = 0;
- k; --k, i1 = i1 == n - 1 ? 1 : i1 + 1, // i1 = i1 + 1 mod n - 1
- j = j == s2 - 1 ? 0 : j + 1 ) { // j = j+1 mod s2
- r = state[i1] ^ c * (r ^ r >> (w - 2));
- r += j + mixer_type(seed[m * j]) +
- (m == 1 || 2 * j + 1 == s ? mixer_type(0) :
- mixer_type(seed[m * j + 1]) << (w - 32));
- r &= mask;
- state[i1] = r;
- }
- for (unsigned k = n - 1; k; --k,
- i1 = i1 == n - 1 ? 1 : i1 + 1) { // i1 = i1 + 1 mod n - 1
- r = state[i1] ^ d * (r ^ r >> (w - 2));
- r -= i1;
- r &= mask;
- state[i1] = r;
- }
- state[0] = typename mixer_t::type(1) << (w - 1);
- }
- }
-
- template void MixerMT1::
- SeedToState(const std::vector& seed,
- mixer_type state[], unsigned n) throw() {
- // This is the algorithm given in the seed_seq class described in
- // http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2006/n2079.pdf It is
- // a modification of
- // http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/CODES/mt19937ar.c
- // http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/VERSIONS/C-LANG/mt19937-64.c
- const unsigned s = unsigned(seed.size());
- const unsigned w = mixer_t::width;
-
- mixer_type r = (a + s) & mask;
- state[0] = r;
- for (unsigned k = 1; k < n; ++k) {
- r = b * (r ^ r >> (w - 2)) + k;
- r &= mask;
- state[k] = r;
- }
- if (s > 0) {
- const unsigned m = mixer_t::width / 32,
- s2 = (s + m - 1)/m;
- unsigned i1 = 0;
- for (unsigned k = (n > s2 ? n : s2), j = 0;
- k; --k, i1 = i1 == n - 1 ? 0 : i1 + 1, // i1 = i1 + 1 mod n
- j = j == s2 - 1 ? 0 : j + 1 ) { // j = j+1 mod s2
- r = state[i1] ^ c * (r ^ r >> (w - 2));
- r += j + mixer_type(seed[m * j]) +
- (m == 1 || 2 * j + 1 == s ? mixer_type(0) :
- mixer_type(seed[m * j + 1]) << (w - 32));
- r &= mask;
- state[i1] = r;
- }
- for (unsigned k = n; k; --k,
- i1 = i1 == n - 1 ? 0 : i1 + 1) { // i1 = i1 + 1 mod n
- r = state[i1] ^ d * (r ^ r >> (w - 2));
- r -= i1;
- r &= mask;
- state[i1] = r;
- }
- }
- }
-
- void MixerSFMT::SeedToState(const std::vector& seed,
- mixer_type state[], unsigned n) throw() {
- // This is adapted from the routine init_by_array by Mutsuo Saito given in
- // http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/SFMT-src-1.2.tar.gz
-
- if (n == 0)
- return; // Nothing to do
-
- const unsigned s = unsigned(seed.size()),
- // Add treatment of small n with lag = (n - 1)/2 for n <= 7. In
- // particular, the first operation (xor or plus) in each for loop
- // involves three distinct indices for n > 2.
- lag = n >= 623 ? 11 : (n >= 68 ? 7 : (n >= 39 ? 5 :
- (n >= 7 ? 3 : (n - 1)/2))),
- // count = max( s + 1, n )
- count = s + 1 > n ? s + 1 : n;
-
- std::fill(state, state + n, mixer_type(a));
- const unsigned w = mixer_t::width;
-
- unsigned i = 0, k = (n - lag) / 2, l = k + lag;
- mixer_type r = state[n - 1];
- for (unsigned j = 0; j < count; ++j,
- i = i == n - 1 ? 0 : i + 1,
- k = k == n - 1 ? 0 : k + 1,
- l = l == n - 1 ? 0 : l + 1) {
- // Here r = state[(j - 1) mod n]
- // i = j mod n
- // k = (j + (n - lag)/2) mod n
- // l = (j + (n - lag)/2 + lag) mod n
- r ^= state[i] ^ state[k];
- r &= mask;
- r = b * (r ^ r >> (w - 5));
- state[k] += r;
- r += i + (j > s ? 0 : (j ? mixer_type(seed[j - 1]) : s));
- state[l] += r;
- state[i] = r;
- }
-
- for (unsigned j = n; j; --j,
- i = i == n - 1 ? 0 : i + 1,
- k = k == n - 1 ? 0 : k + 1,
- l = l == n - 1 ? 0 : l + 1) {
- // Here r = state[(i - 1) mod n]
- // k = (i + (n - lag)/2) mod n
- // l = (i + (n - lag)/2 + lag) mod n
- r += state[i] + state[k];
- r &= mask;
- r = c * (r ^ r >> (w - 5));
- r &= mask;
- state[k] ^= r;
- r -= i;
- r &= mask;
- state[l] ^= r;
- state[i] = r;
- }
- }
-
- // RandomAlgorithm implementation
-
- // Here, input is I, J = I + 1, K = I + M; output is I = I + N (mod N)
-
-#define MT19937_STEP(I, J, K) statev[I] = statev[K] ^ \
- (statev[J] & engine_type(1) ? magic : engine_type(0)) ^ \
- ((statev[I] & upper) | (statev[J] & lower)) >> 1
-
- // The code is cleaned up a little from Hagita's Fortran version by getting
- // rid of the unnecessary masking by YMASK and by using simpler logic to
- // restore the correct value of _state[0].
- //
- // Here input is J = I + N - 1, K = I + M - 1, and p = y[I] (only the high
- // bits are used); output _state[I] and p = y[I - 1].
-
-#define MT19937_REVSTEP(I, J, K) { \
- engine_type q = statev[J] ^ statev[K], s = q >> (width - 1); \
- q = (q ^ (s ? magic : engine_type(0))) << 1 | s; \
- statev[I] = (p & upper) | (q & lower); \
- p = q; \
- }
-
- template
- void MT19937::Transition(long long count, internal_type statev[])
- throw() {
- if (count > 0)
- for (; count; --count) {
- // This ONLY uses high bit of statev[0]
- unsigned i = 0;
- for (; i < N - M; ++i) MT19937_STEP(i, i + 1, i + M );
- for (; i < N - 1; ++i) MT19937_STEP(i, i + 1, i + M - N);
- MT19937_STEP(N - 1, 0, M - 1); // i = N - 1
- }
- else if (count < 0)
- for (; count; ++count) {
- // This ONLY uses high bit of statev[0]
- engine_type p = statev[0];
- // Fix low bits of statev[0] and compute y[-1]
- MT19937_REVSTEP(0, N - 1, M - 1); // i = N
- unsigned i = N - 1;
- for (; i > N - M; --i) MT19937_REVSTEP(i, i - 1, i + M - 1 - N);
- for (; i ; --i) MT19937_REVSTEP(i, i - 1, i + M - 1 );
- MT19937_REVSTEP(0, N - 1, M - 1); // i = 0
- }
- }
-
-#undef MT19937_STEP
-#undef MT19937_REVSTEP
-
- template
- void MT19937::NormalizeState(engine_type state[]) throw() {
-
- // Perform the MT-specific sanity check on the resulting state ensuring
- // that the significant 19937 bits are not all zero.
- state[0] &= upper; // Mask out unused bits
- unsigned i = 0;
- while (i < N && state[i] == 0)
- ++i;
- if (i >= N)
- state[0] = engine_type(1) << (width - 1); // with prob 2^-19937
-
- // This sets the low R bits of _state[0] consistent with the rest of the
- // state. Needed to handle SetCount(-N); Ran32(); immediately following
- // reseeding. This wasn't required in the original code because a
- // Transition was always done first.
- engine_type q = state[N - 1] ^ state[M - 1], s = q >> (width - 1);
- q = (q ^ (s ? magic : engine_type(0))) << 1 | s;
- state[0] = (state[0] & upper) | (q & lower);
- }
-
- template
- void MT19937::CheckState(const engine_type state[],
- Random_u32::type& check) {
- engine_type x = 0;
- Random_u32::type c = check;
- for (unsigned i = 0; i < N; ++i) {
- engine_t::CheckSum(state[i], c);
- x |= state[i];
- }
- if (x == 0)
- throw RandomErr("MT19937: All-zero state");
-
- // There are only width*(N-1) + 1 = 19937 independent bits of state. Thus
- // the low width-1 bits of _state[0] are derivable from the other bits in
- // state. Verify that the redundant bits bits are consistent.
- engine_type q = state[N - 1] ^ state[M - 1], s = q >> (width - 1);
- q = (q ^ (s ? magic : engine_type(0))) << 1 | s;
- if ((q ^ state[0]) & lower)
- throw RandomErr("MT19937: Invalid state");
-
- check = c;
- }
-
-#if defined(HAVE_SSE2) && HAVE_SSE2
-
- // Transition is from Saito's Master's Thesis
- // http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/M062821.pdf
- //
- // This implements
- //
- // w_{i+N} = w_i A xor w_M B xor w_{i+N-2} C xor w_{i+N-1} D
- //
- // where w_i is a 128-bit word and
- //
- // w A = (w << 8)_128 xor w
- // w B = (w >> 11)_32 & MSK
- // w C = (w >> 8)_128
- // w D = (w << 18)_32
- //
- // Here the _128 means shift is of whole 128-bit word. _32 means the shifts
- // are independently done on each 32-bit word.
- //
- // In SFMT19937_STEP32 and SFMT19937_STEP64 input is I, J = I + M and output
- // is I = I + N (mod N). On input, s and r give state for I + N - 2 and I +
- // N - 1; on output s and r give state for I + N - 1 and I + N. The
- // implementation of 128-bit operations is open-coded in a portable fashion
- // (with LSB ordering).
- //
- // N.B. Here N and M are the lags in units of BitWidth words and so are 4
- // (for u32 implementation) or 2 (for u64 implementation) times bigger than
- // defined in Saito's thesis.
-
- // This is adapted from SFMT-sse.c in the SFMT 1.2 distribution.
- // The order of instructions has been rearranged to increase the
- // speed slightly
-
-#define SFMT19937_STEP128(I, J) { \
- internal_type x = _mm_load_si128(statev + I), \
- y = _mm_srli_epi32(statev[J], 11), \
- z = _mm_srli_si128(s, 1); \
- s = _mm_slli_epi32(r, 18); \
- z = _mm_xor_si128(z, x); \
- x = _mm_slli_si128(x, 1); \
- z = _mm_xor_si128(z, s); \
- y = _mm_and_si128(y, m); \
- z = _mm_xor_si128(z, x); \
- s = r; \
- r = _mm_xor_si128(z, y); \
- _mm_store_si128(statev + I, r); \
- }
-
- // This undoes SFMT19937_STEP. Trivially, we have
- //
- // w_i A = w_{i+N} xor w_{i+M} B xor w_{i+N-2} C xor w_{i+N-1} D
- //
- // Given w_i A we can determine w_i from the observation that A^16 =
- // identity, thus
- //
- // w_i = (w_i A) A^15
- //
- // Because x A^(2^n) = x << (8*2^n) xor x, the operation y = x A^15 can be
- // implemented as
- //
- // y' = (x << 64)_128 xor x = x A^8
- // y'' = (y' << 32)_128 xor y' = y' A^4 = x A^12
- // y''' = (y'' << 16)_128 xor y'' = y'' A^2 = x A^14
- // y = (y''' << 8)_128 xor y''' = y''' A = x A^15
- //
- // Here input is I = I + N, J = I + M, K = I + N - 2, L = I + N -1, and
- // output is I = I.
- //
- // This is about 15-35% times slower than SFMT19937_STEPNN, because (1) there
- // doesn't appear to be a straightforward way of saving intermediate results
- // across calls as with SFMT19937_STEPNN and (2) w A^15 is slower to compute
- // than w A.
-
-#define SFMT19937_REVSTEP128(I, J, K, L) { \
- internal_type x = _mm_load_si128(statev + I), \
- y = _mm_srli_epi32(statev[J], 11), \
- z = _mm_slli_epi32(statev[L], 18); \
- y = _mm_and_si128(y, m); \
- x = _mm_xor_si128(x, _mm_srli_si128(statev[K], 1)); \
- x = _mm_xor_si128(x, z); \
- x = _mm_xor_si128(x, y); \
- x = _mm_xor_si128(_mm_slli_si128(x, 8), x); \
- x = _mm_xor_si128(_mm_slli_si128(x, 4), x); \
- x = _mm_xor_si128(_mm_slli_si128(x, 2), x); \
- x = _mm_xor_si128(_mm_slli_si128(x, 1), x); \
- _mm_store_si128(statev + I, x); \
- }
-
- template
- void SFMT19937::Transition(long long count,
- internal_type statev[])
- throw() {
- const internal_type m = _mm_set_epi32(magic3, magic2, magic1, magic0);
- if (count > 0) {
- internal_type s = _mm_load_si128(statev + N128 - 2),
- r = _mm_load_si128(statev + N128 - 1);
- for (; count; --count) {
- unsigned i = 0;
- for (; i + M128 < N128; ++i) SFMT19937_STEP128(i, i + M128 );
- for (; i < N128 ; ++i) SFMT19937_STEP128(i, i + M128 - N128);
- }
- } else if (count < 0)
- for (; count; ++count) {
- unsigned i = N128;
- for (; i + M128 > N128;) {
- --i; SFMT19937_REVSTEP128(i, i + M128 - N128, i - 2, i - 1);
- }
- for (; i > 2;) {
- --i; SFMT19937_REVSTEP128(i, i + M128, i - 2, i - 1);
- }
- SFMT19937_REVSTEP128(1, M128 + 1, N128 - 1, 0 ); // i = 1
- SFMT19937_REVSTEP128(0, M128 , N128 - 2, N128 - 1); // i = 0
- }
- }
-
-#undef SFMT19937_STEP128
-#undef SFMT19937_REVSTEP128
-
-#elif defined(HAVE_ALTIVEC) && HAVE_ALTIVEC
-
- // The Altivec versions of SFMT19937_{,REV}STEP128 are simply translated from
- // the SSE2 versions. The only significant differences arise because of the
- // MSB ordering of the PowerPC. This means that the 32-bit and 64-bit
- // versions are no different because 32-bit and 64-bit words don't pack
- // together in the same way as on an SSE2 machine (see the two definitions of
- // magic). This also means that the 128-bit byte shifts on an LSB machine
- // change into more complicated byte permutations.
-
-#define ALTIVEC_PERM(X, P) vec_perm(X, P, P)
-
-#define SFMT19937_STEP128(I, J) { \
- internal_type x = statev[I], \
- z = vec_xor(vec_xor(ALTIVEC_PERM(s, right1), x), \
- vec_sl(r, bitleft)); \
- s = r; \
- r = vec_xor(z, \
- vec_xor(ALTIVEC_PERM(x, left1), \
- vec_and(vec_sr(statev[J], bitright), \
- magic))); \
- statev[I] = r; \
- }
-
-#define SFMT19937_REVSTEP128(I, J, K, L) { \
- internal_type x = statev[I], \
- y = vec_sr(statev[J], bitright), \
- z = vec_sl(statev[L], bitleft); \
- y = vec_and(y, magic); \
- x = vec_xor(x, ALTIVEC_PERM(statev[K], right1)); \
- x = vec_xor(x, z); \
- x = vec_xor(x, y); \
- x = vec_xor(ALTIVEC_PERM(x, left8), x); \
- x = vec_xor(ALTIVEC_PERM(x, left4), x); \
- x = vec_xor(ALTIVEC_PERM(x, left2), x); \
- statev[I] = vec_xor(ALTIVEC_PERM(x, left1), x); \
- }
-
- template
- void SFMT19937::Transition(long long count,
- internal_type statev[])
- throw() {
- const internal_type magic = width == 32 ?
- (vector unsigned)(magic0, magic1, magic2, magic3) :
- (vector unsigned)(magic1, magic0, magic3, magic2),
- bitleft = (vector unsigned)(18, 18, 18, 18),
- bitright = (vector unsigned)(11, 11, 11, 11);
- // Shift left and right by 1 byte. Note that vec_perm(X, Y, P) glues X and
- // Y together into a 32-byte quantity and then the 16-byte permutation
- // vector P specifies which bytes to put into the 16-byte output. We
- // follow here the convention of using Y = P and using the zero entries in
- // P to allow zero bytes to be introduces into the shifted output. The
- // following describes how the left1 table (32-bit version) is produced:
- //
- // Byte layout of original with LSB ordering
- // 33 32 31 30 23 22 21 20 13 12 11 10 03 02 01 00
- // shift left by 1 byte (z means zeros enter)
- // 32 31 30 23 22 21 20 13 12 11 10 03 02 01 00 zz
- //
- // Rearrange original to LSB order in 4-byte units
- // 03 02 01 00 13 12 11 10 23 22 21 20 33 32 31 30
- // with sequential MSB byte indices
- // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
- //
- // Rearrange shift left verion to LSB order in 4-byte units
- // 02 01 00 zz 12 11 10 03 22 21 20 13 32 31 30 23
- // with corresponding MSB byte indices
- // 1 2 3 z 5 6 7 0 9 10 11 4 13 14 15 8
- //
- // Replace byte index at x by 16 + index of 0 = 16 + 7 = 23 to give
- // 1 2 3 23 5 6 7 0 9 10 11 4 13 14 15 8
- const vector unsigned char left1 = width == 32 ?
- (vector unsigned char)(1,2,3,23, 5,6,7,0, 9,10,11,4, 13,14,15,8) :
- (vector unsigned char)(1,2,3,4,5,6,7,31, 9,10,11,12,13,14,15,0),
- right1 = width == 32 ?
- (vector unsigned char)(7,0,1,2, 11,4,5,6, 15,8,9,10, 17,12,13,14) :
- (vector unsigned char)(15,0,1,2,3,4,5,6, 17,8,9,10,11,12,13,14);
- if (count > 0) {
- internal_type s = statev[N128 - 2],
- r = statev[N128 - 1];
- for (; count; --count) {
- unsigned i = 0;
- for (; i + M128 < N128; ++i) SFMT19937_STEP128(i, i + M128 );
- for (; i < N128 ; ++i) SFMT19937_STEP128(i, i + M128 - N128);
- }
- } else if (count < 0) {
- // leftN shifts left by N bytes.
- const vector unsigned char left2 = width == 32 ?
- (vector unsigned char)(2,3,22,22, 6,7,0,1, 10,11,4,5, 14,15,8,9) :
- (vector unsigned char)(2,3,4,5,6,7,30,30, 10,11,12,13,14,15,0,1),
- left4 = width == 32 ?
- (vector unsigned char)(20,20,20,20, 0,1,2,3, 4,5,6,7, 8,9,10,11) :
- (vector unsigned char)(4,5,6,7,28,28,28,28, 12,13,14,15,0,1,2,3),
- left8 = (vector unsigned char)(24,24,24,24,24,24,24,24,0,1,2,3,4,5,6,7);
- for (; count; ++count) {
- unsigned i = N128;
- for (; i + M128 > N128;) {
- --i; SFMT19937_REVSTEP128(i, i + M128 - N128, i - 2, i - 1);
- }
- for (; i > 2;) {
- --i; SFMT19937_REVSTEP128(i, i + M128, i - 2, i - 1);
- }
- SFMT19937_REVSTEP128(1, M128 + 1, N128 - 1, 0 ); // i = 1
- SFMT19937_REVSTEP128(0, M128 , N128 - 2, N128 - 1); // i = 0
- }
- }
- }
-
-#undef SFMT19937_STEP128
-#undef SFMT19937_REVSTEP128
-#undef ALTIVEC_PERM
-
-#else // neither HAVE_SSE2 or HAVE_ALTIVEC
-
-#define SFMT19937_STEP32(I, J) { \
- internal_type t = statev[I] ^ statev[I] << 8 ^ \
- (statev[J] >> 11 & magic0) ^ \
- (s0 >> 8 | s1 << 24) ^ r0 << 18; \
- s0 = r0; r0 = t & mask; \
- t = statev[I + 1] ^ \
- (statev[I + 1] << 8 | statev[I] >> 24) ^ \
- (statev[J + 1] >> 11 & magic1) ^ \
- (s1 >> 8 | s2 << 24) ^ r1 << 18; \
- s1 = r1; r1 = t & mask; \
- t = statev[I + 2] ^ \
- (statev[I + 2] << 8 | statev[I + 1] >> 24) ^ \
- (statev[J + 2] >> 11 & magic2) ^ \
- (s2 >> 8 | s3 << 24) ^ r2 << 18; \
- s2 = r2; r2 = t & mask; \
- t = statev[I + 3] ^ \
- (statev[I + 3] << 8 | statev[I + 2] >> 24) ^ \
- (statev[J + 3] >> 11 & magic3) ^ s3 >> 8 ^ r3 << 18; \
- s3 = r3; r3 = t & mask; \
- statev[I ] = r0; statev[I + 1] = r1; \
- statev[I + 2] = r2; statev[I + 3] = r3; \
- }
-
-#define SFMT19937_REVSTEP32(I, J, K, L) { \
- internal_type \
- t0 = (statev[I] ^ (statev[J] >> 11 & magic0) ^ \
- (statev[K] >> 8 | statev[K + 1] << 24) ^ \
- statev[L] << 18) & mask, \
- t1 = (statev[I + 1] ^ \
- (statev[J + 1] >> 11 & magic1) ^ \
- (statev[K + 1] >> 8 | statev[K + 2] << 24) ^ \
- statev[L + 1] << 18) & mask, \
- t2 = (statev[I + 2] ^ \
- (statev[J + 2] >> 11 & magic2) ^ \
- (statev[K + 2] >> 8 | statev[K + 3] << 24) ^ \
- statev[L + 2] << 18) & mask, \
- t3 = (statev[I + 3] ^ \
- (statev[J + 3] >> 11 & magic3) ^ \
- statev[K + 3] >> 8 ^ \
- statev[L + 3] << 18) & mask; \
- t3 ^= t1; t2 ^= t0; t3 ^= t2; t2 ^= t1; t1 ^= t0; \
- t3 ^= t2 >> 16 | (t3 << 16 & mask); \
- t2 ^= t1 >> 16 | (t2 << 16 & mask); \
- t1 ^= t0 >> 16 | (t1 << 16 & mask); \
- t0 ^= t0 << 16 & mask; \
- statev[I ] = t0 ^ (t0 << 8 & mask); \
- statev[I + 1] = t1 ^ (t0 >> 24 | (t1 << 8 & mask)); \
- statev[I + 2] = t2 ^ (t1 >> 24 | (t2 << 8 & mask)); \
- statev[I + 3] = t3 ^ (t2 >> 24 | (t3 << 8 & mask)); \
- }
-
- template<>
- void SFMT19937::Transition(long long count,
- internal_type statev[])
- throw() {
- if (count > 0) {
- // x[i+N] = g(x[i], x[i+M], x[i+N-2], x[i,N-1])
- internal_type
- s0 = statev[N - 8], s1 = statev[N - 7],
- s2 = statev[N - 6], s3 = statev[N - 5],
- r0 = statev[N - 4], r1 = statev[N - 3],
- r2 = statev[N - 2], r3 = statev[N - 1];
- for (; count; --count) {
- unsigned i = 0;
- for (; i + M < N; i += R) SFMT19937_STEP32(i, i + M );
- for (; i < N ; i += R) SFMT19937_STEP32(i, i + M - N);
- }
- } else if (count < 0)
- for (; count; ++count) {
- unsigned i = N;
- for (; i + M > N;) {
- i -= R; SFMT19937_REVSTEP32(i, i + M - N, i - 2 * R, i - R);
- }
- for (; i > 2 * R;) {
- i -= R; SFMT19937_REVSTEP32(i, i + M , i - 2 * R, i - R);
- }
- SFMT19937_REVSTEP32(R, M + R, N - R, 0 ); // i = R
- SFMT19937_REVSTEP32(0, M , N - 2 * R, N - R); // i = 0
- }
- }
-
-#undef SFMT19937_STEP32
-#undef SFMT19937_REVSTEP32
-
-#define SFMT19937_STEP64(I, J) { \
- internal_type t = statev[I] ^ statev[I] << 8 ^ \
- (statev[J] >> 11 & magic0) ^ \
- (s0 >> 8 | s1 << 56) ^ (r0 << 18 & mask18); \
- s0 = r0; r0 = t & mask; \
- t = statev[I + 1] ^ \
- (statev[I + 1] << 8 | statev[I] >> 56) ^ \
- (statev[J + 1] >> 11 & magic1) ^ \
- s1 >> 8 ^ (r1 << 18 & mask18); \
- s1 = r1; r1 = t & mask; \
- statev[I] = r0; statev[I + 1] = r1; \
- }
-
- // In combining the left and right shifts to simulate a 128-bit shift we
- // usually use or. However we can equivalently use xor (e.g., t1 << 8 ^ t0
- // >> 56 instead of t1 ^ t1 << 8 | t0 >> 56) and this speeds up the code if
- // used in some places.
-
-#define SFMT19937_REVSTEP64(I, J, K, L) { \
- internal_type \
- t0 = statev[I] ^ (statev[J] >> 11 & magic0) ^ \
- (statev[K] >> 8 | (statev[K + 1] << 56 & mask)) ^ \
- (statev[L] << 18 & mask18), \
- t1 = statev[I + 1] ^ (statev[J + 1] >> 11 & magic1) ^ \
- statev[K + 1] >> 8 ^ (statev[L + 1] << 18 & mask18); \
- t1 ^= t0; \
- t1 ^= t0 >> 32 ^ (t1 << 32 & mask); \
- t0 ^= t0 << 32 & mask; \
- t1 ^= t0 >> 48 ^ (t1 << 16 & mask); \
- t0 ^= t0 << 16 & mask; \
- statev[I ] = t0 ^ (t0 << 8 & mask); \
- statev[I + 1] = t1 ^ t0 >> 56 ^ (t1 << 8 & mask); \
- }
-
- template<>
- void SFMT19937::Transition(long long count,
- internal_type statev[])
- throw() {
- // x[i+N] = g(x[i], x[i+M], x[i+N-2], x[i,N-1])
- if (count > 0) {
- internal_type
- s0 = statev[N - 4], s1 = statev[N - 3],
- r0 = statev[N - 2], r1 = statev[N - 1];
- for (; count; --count) {
- unsigned i = 0;
- for (; i + M < N; i += R) SFMT19937_STEP64(i, i + M );
- for (; i < N ; i += R) SFMT19937_STEP64(i, i + M - N);
- }
- } else if (count < 0)
- for (; count; ++count) {
- unsigned i = N;
- for (; i + M > N;) {
- i -= R; SFMT19937_REVSTEP64(i, i + M - N, i - 2 * R, i - R);
- }
- for (; i > 2 * R;) {
- i -= R; SFMT19937_REVSTEP64(i, i + M , i - 2 * R, i - R);
- }
- SFMT19937_REVSTEP64(R, M + R, N - R, 0 ); // i = R
- SFMT19937_REVSTEP64(0, M , N - 2 * R, N - R); // i = 0
- }
- }
-
-#undef SFMT19937_STEP64
-#undef SFMT19937_REVSTEP64
-
-#endif // HAVE_SSE2 and HAVE_ALTIVEC
-
- template<>
- void SFMT19937::NormalizeState(engine_type state[]) throw() {
- // Carry out the Period Certification for SFMT19937
- engine_type inner = (state[0] & PARITY0) ^ (state[1] & PARITY1) ^
- (state[2] & PARITY2) ^ (state[3] & PARITY3);
- for (unsigned s = 16; s; s >>= 1)
- inner ^= inner >> s;
- STATIC_ASSERT(PARITY_LSB < 32 && PARITY0 & 1u << PARITY_LSB,
- "inconsistent PARITY_LSB or PARITY0");
- // Now inner & 1 is the parity of the number of 1 bits in w_0 & p.
- if ((inner & 1u) == 0)
- // Change bit of w_0 corresponding to LSB of PARITY
- state[PARITY_LSB >> 5] ^= engine_type(1u) << (PARITY_LSB & 31u);
- }
-
- template<>
- void SFMT19937::NormalizeState(engine_type state[]) throw() {
- // Carry out the Period Certification for SFMT19937
- engine_type inner = (state[0] & PARITY0) ^ (state[1] & PARITY1);
- for (unsigned s = 32; s; s >>= 1)
- inner ^= inner >> s;
- STATIC_ASSERT(PARITY_LSB < 64 && PARITY0 & 1u << PARITY_LSB,
- "inconsistent PARITY_LSB or PARITY0");
- // Now inner & 1 is the parity of the number of 1 bits in w_0 & p.
- if ((inner & 1u) == 0)
- // Change bit of w_0 corresponding to LSB of PARITY
- state[PARITY_LSB >> 6] ^= engine_type(1u) << (PARITY_LSB & 63u);
- }
-
- template
- void SFMT19937::CheckState(const engine_type state[],
- Random_u32::type& check) {
- engine_type x = 0;
- Random_u32::type c = check;
- for (unsigned i = 0; i < N; ++i) {
- engine_t::CheckSum(state[i], c);
- x |= state[i];
- }
- if (x == 0)
- throw RandomErr("SFMT19937: All-zero state");
- check = c;
- }
-
- // RandomPower2 implementation
-
-#if RANDOMLIB_POWERTABLE
- // Powers of two. Just use floats here. As long as there's no overflow
- // or underflow these are exact. In particular they can be cast to
- // doubles or long doubles with no error.
- const float RandomPower2::power2[maxpow - minpow + 1] = {
-#if RANDOMLIB_LONGDOUBLEPREC > 64
- // It would be nice to be able to use the C99 notation of 0x1.0p-120
- // for 2^-120 here.
- 1/1329227995784915872903807060280344576.f, // 2^-120
- 1/664613997892457936451903530140172288.f, // 2^-119
- 1/332306998946228968225951765070086144.f, // 2^-118
- 1/166153499473114484112975882535043072.f, // 2^-117
- 1/83076749736557242056487941267521536.f, // 2^-116
- 1/41538374868278621028243970633760768.f, // 2^-115
- 1/20769187434139310514121985316880384.f, // 2^-114
- 1/10384593717069655257060992658440192.f, // 2^-113
- 1/5192296858534827628530496329220096.f, // 2^-112
- 1/2596148429267413814265248164610048.f, // 2^-111
- 1/1298074214633706907132624082305024.f, // 2^-110
- 1/649037107316853453566312041152512.f, // 2^-109
- 1/324518553658426726783156020576256.f, // 2^-108
- 1/162259276829213363391578010288128.f, // 2^-107
- 1/81129638414606681695789005144064.f, // 2^-106
- 1/40564819207303340847894502572032.f, // 2^-105
- 1/20282409603651670423947251286016.f, // 2^-104
- 1/10141204801825835211973625643008.f, // 2^-103
- 1/5070602400912917605986812821504.f, // 2^-102
- 1/2535301200456458802993406410752.f, // 2^-101
- 1/1267650600228229401496703205376.f, // 2^-100
- 1/633825300114114700748351602688.f, // 2^-99
- 1/316912650057057350374175801344.f, // 2^-98
- 1/158456325028528675187087900672.f, // 2^-97
- 1/79228162514264337593543950336.f, // 2^-96
- 1/39614081257132168796771975168.f, // 2^-95
- 1/19807040628566084398385987584.f, // 2^-94
- 1/9903520314283042199192993792.f, // 2^-93
- 1/4951760157141521099596496896.f, // 2^-92
- 1/2475880078570760549798248448.f, // 2^-91
- 1/1237940039285380274899124224.f, // 2^-90
- 1/618970019642690137449562112.f, // 2^-89
- 1/309485009821345068724781056.f, // 2^-88
- 1/154742504910672534362390528.f, // 2^-87
- 1/77371252455336267181195264.f, // 2^-86
- 1/38685626227668133590597632.f, // 2^-85
- 1/19342813113834066795298816.f, // 2^-84
- 1/9671406556917033397649408.f, // 2^-83
- 1/4835703278458516698824704.f, // 2^-82
- 1/2417851639229258349412352.f, // 2^-81
- 1/1208925819614629174706176.f, // 2^-80
- 1/604462909807314587353088.f, // 2^-79
- 1/302231454903657293676544.f, // 2^-78
- 1/151115727451828646838272.f, // 2^-77
- 1/75557863725914323419136.f, // 2^-76
- 1/37778931862957161709568.f, // 2^-75
- 1/18889465931478580854784.f, // 2^-74
- 1/9444732965739290427392.f, // 2^-73
- 1/4722366482869645213696.f, // 2^-72
- 1/2361183241434822606848.f, // 2^-71
- 1/1180591620717411303424.f, // 2^-70
- 1/590295810358705651712.f, // 2^-69
- 1/295147905179352825856.f, // 2^-68
- 1/147573952589676412928.f, // 2^-67
- 1/73786976294838206464.f, // 2^-66
- 1/36893488147419103232.f, // 2^-65
-#endif
- 1/18446744073709551616.f, // 2^-64
- 1/9223372036854775808.f, // 2^-63
- 1/4611686018427387904.f, // 2^-62
- 1/2305843009213693952.f, // 2^-61
- 1/1152921504606846976.f, // 2^-60
- 1/576460752303423488.f, // 2^-59
- 1/288230376151711744.f, // 2^-58
- 1/144115188075855872.f, // 2^-57
- 1/72057594037927936.f, // 2^-56
- 1/36028797018963968.f, // 2^-55
- 1/18014398509481984.f, // 2^-54
- 1/9007199254740992.f, // 2^-53
- 1/4503599627370496.f, // 2^-52
- 1/2251799813685248.f, // 2^-51
- 1/1125899906842624.f, // 2^-50
- 1/562949953421312.f, // 2^-49
- 1/281474976710656.f, // 2^-48
- 1/140737488355328.f, // 2^-47
- 1/70368744177664.f, // 2^-46
- 1/35184372088832.f, // 2^-45
- 1/17592186044416.f, // 2^-44
- 1/8796093022208.f, // 2^-43
- 1/4398046511104.f, // 2^-42
- 1/2199023255552.f, // 2^-41
- 1/1099511627776.f, // 2^-40
- 1/549755813888.f, // 2^-39
- 1/274877906944.f, // 2^-38
- 1/137438953472.f, // 2^-37
- 1/68719476736.f, // 2^-36
- 1/34359738368.f, // 2^-35
- 1/17179869184.f, // 2^-34
- 1/8589934592.f, // 2^-33
- 1/4294967296.f, // 2^-32
- 1/2147483648.f, // 2^-31
- 1/1073741824.f, // 2^-30
- 1/536870912.f, // 2^-29
- 1/268435456.f, // 2^-28
- 1/134217728.f, // 2^-27
- 1/67108864.f, // 2^-26
- 1/33554432.f, // 2^-25
- 1/16777216.f, // 2^-24
- 1/8388608.f, // 2^-23
- 1/4194304.f, // 2^-22
- 1/2097152.f, // 2^-21
- 1/1048576.f, // 2^-20
- 1/524288.f, // 2^-19
- 1/262144.f, // 2^-18
- 1/131072.f, // 2^-17
- 1/65536.f, // 2^-16
- 1/32768.f, // 2^-15
- 1/16384.f, // 2^-14
- 1/8192.f, // 2^-13
- 1/4096.f, // 2^-12
- 1/2048.f, // 2^-11
- 1/1024.f, // 2^-10
- 1/512.f, // 2^-9
- 1/256.f, // 2^-8
- 1/128.f, // 2^-7
- 1/64.f, // 2^-6
- 1/32.f, // 2^-5
- 1/16.f, // 2^-4
- 1/8.f, // 2^-3
- 1/4.f, // 2^-2
- 1/2.f, // 2^-1
- 1.f, // 2^0
- 2.f, // 2^1
- 4.f, // 2^2
- 8.f, // 2^3
- 16.f, // 2^4
- 32.f, // 2^5
- 64.f, // 2^6
- 128.f, // 2^7
- 256.f, // 2^8
- 512.f, // 2^9
- 1024.f, // 2^10
- 2048.f, // 2^11
- 4096.f, // 2^12
- 8192.f, // 2^13
- 16384.f, // 2^14
- 32768.f, // 2^15
- 65536.f, // 2^16
- 131072.f, // 2^17
- 262144.f, // 2^18
- 524288.f, // 2^19
- 1048576.f, // 2^20
- 2097152.f, // 2^21
- 4194304.f, // 2^22
- 8388608.f, // 2^23
- 16777216.f, // 2^24
- 33554432.f, // 2^25
- 67108864.f, // 2^26
- 134217728.f, // 2^27
- 268435456.f, // 2^28
- 536870912.f, // 2^29
- 1073741824.f, // 2^30
- 2147483648.f, // 2^31
- 4294967296.f, // 2^32
- 8589934592.f, // 2^33
- 17179869184.f, // 2^34
- 34359738368.f, // 2^35
- 68719476736.f, // 2^36
- 137438953472.f, // 2^37
- 274877906944.f, // 2^38
- 549755813888.f, // 2^39
- 1099511627776.f, // 2^40
- 2199023255552.f, // 2^41
- 4398046511104.f, // 2^42
- 8796093022208.f, // 2^43
- 17592186044416.f, // 2^44
- 35184372088832.f, // 2^45
- 70368744177664.f, // 2^46
- 140737488355328.f, // 2^47
- 281474976710656.f, // 2^48
- 562949953421312.f, // 2^49
- 1125899906842624.f, // 2^50
- 2251799813685248.f, // 2^51
- 4503599627370496.f, // 2^52
- 9007199254740992.f, // 2^53
- 18014398509481984.f, // 2^54
- 36028797018963968.f, // 2^55
- 72057594037927936.f, // 2^56
- 144115188075855872.f, // 2^57
- 288230376151711744.f, // 2^58
- 576460752303423488.f, // 2^59
- 1152921504606846976.f, // 2^60
- 2305843009213693952.f, // 2^61
- 4611686018427387904.f, // 2^62
- 9223372036854775808.f, // 2^63
- 18446744073709551616.f, // 2^64
- };
-#endif
-
- // RandomEngine (and implicitly RandomAlgorithm and RandomMixer)
- // instantiations. The first 4 (using MixerMT[01]) are not recommended.
- template class RandomEngine< MT19937, MixerMT0 >;
- template class RandomEngine< MT19937, MixerMT0 >;
- template class RandomEngine< MT19937, MixerMT1 >;
- template class RandomEngine< MT19937, MixerMT1 >;
-
- template class RandomEngine< MT19937, MixerSFMT>;
- template class RandomEngine< MT19937, MixerSFMT>;
- template class RandomEngine, MixerSFMT>;
- template class RandomEngine, MixerSFMT>;
-
- // RandomCanonial instantiations
-
- template<> RandomCanonical
- RandomCanonical::Global = RandomCanonical();
- template<> RandomCanonical
- RandomCanonical::Global = RandomCanonical();
- template<> RandomCanonical
- RandomCanonical::Global = RandomCanonical();
- template<> RandomCanonical
- RandomCanonical::Global = RandomCanonical();
-
-} // namespace RandomLib
diff --git a/include/RandomLib/Config.h b/include/RandomLib/Config.h
deleted file mode 100644
index 94d35c2f..00000000
--- a/include/RandomLib/Config.h
+++ /dev/null
@@ -1,13 +0,0 @@
-#define RANDOMLIB_VERSION_STRING "Unconfigured"
-#define RANDOMLIB_VERSION_MAJOR -1
-#define RANDOMLIB_VERSION_MINOR -1
-#define RANDOMLIB_VERSION_PATCH -1
-
-// Define HAVE_SSE2 to be 1 if Intel/AMD CPU with SSE2 support
-/* #undef HAVE_SSE2 */
-
-// Define HAVE_ALTIVEC to be 1 if Power PC CPU with AltiVec support
-/* #undef HAVE_ALTIVEC */
-
-// Undefine HAVE_LONG_DOUBLE if this type is unknown to the compiler
-#define HAVE_LONG_DOUBLE 1
diff --git a/include/RandomLib/DiscreteNormal.hpp b/include/RandomLib/DiscreteNormal.hpp
deleted file mode 100644
index 6775bae4..00000000
--- a/include/RandomLib/DiscreteNormal.hpp
+++ /dev/null
@@ -1,432 +0,0 @@
-/**
- * \file DiscreteNormal.hpp
- * \brief Header for DiscreteNormal
- *
- * Sample exactly from the discrete normal distribution.
- *
- * Copyright (c) Charles Karney (2013) and licensed
- * under the MIT/X11 License. For more information, see
- * http://randomlib.sourceforge.net/
- **********************************************************************/
-
-#if !defined(RANDOMLIB_DISCRETENORMAL_HPP)
-#define RANDOMLIB_DISCRETENORMAL_HPP 1
-
-#include
-#include
-
-namespace RandomLib {
- /**
- * \brief The discrete normal distribution.
- *
- * Sample integers \e i with probability proportional to
- * \f[
- * \exp\biggl[-\frac12\biggl(\frac{i-\mu}{\sigma}\biggr)^2\biggr],
- * \f]
- * where σ and μ are given as rationals (the ratio of two integers).
- * The sampling is exact (provided that the random generator is ideal). For
- * example
- * \code
- #include
- #include
- #include
-
- int main() {
- RandomLib::Random r; // Create r
- r.Reseed(); // and give it a unique seed
- int sigma_num = 7, sigma_den = 1, mu_num = 1, mu_den = 3;
- RandomLib::DiscreteNormal d(sigma_num, sigma_den,
- mu_num, mu_den);
- for (int i = 0; i < 100; ++i)
- std::cout << d(r) << "\n";
- }
- \endcode
- * prints out 100 samples with σ = 7 and μ = 1/3.
- *
- * The algorithm is much the same as for ExactNormal; for details see
- * - C. F. F. Karney, Sampling exactly from the normal distribution,
- * http://arxiv.org/abs/1303.6257 (Mar. 2013).
- * .
- * That algorithm samples the integer part of the result \e k, samples \e x
- * in [0,1], and (unless rejected) returns s(\e k + \e x), where \e s
- * = ±1. For the discrete case, we sample \e x in [0,1) such that
- * \f[
- * s(k + x) = (i - \mu)/\sigma,
- * \f]
- * or
- * \f[
- * x = s(i - \mu)/\sigma - k
- * \f]
- * The value of \e i which results in the smallest \e x ≥ 0 is
- * \f[
- * i_0 = s\lceil k \sigma + s \mu\rceil
- * \f]
- * so sample
- * \f[
- * i = i_0 + sj
- * \f]
- * where \e j is uniformly distributed in [0, ⌈σ⌉). The
- * corresponding value of \e x is
- * \f[
- * \begin{aligned}
- * x &= \bigl(si_0 - (k\sigma + s\mu)\bigr)/\sigma + j/\sigma\\
- * &= x_0 + j/\sigma,\\
- * x_0 &= \bigl(\lceil k \sigma + s \mu\rceil -
- * (k \sigma + s \mu)\bigr)/\sigma.
- * \end{aligned}
- * \f]
- * After \e x is sampled in this way, it should be rejected if \e x ≥ 1
- * (this is counted with the next larger value of \e k) or if \e x = 0, \e k
- * = 0, and \e s = −1 (to avoid double counting the origin). If \e x
- * is accepted (in Step 4 of the ExactNormal algorithm), then return \e i.
- *
- * When σ and μ are given as rationals, all the arithmetic outlined
- * above can be carried out exactly. The basic rejection techniques used by
- * ExactNormal are exact. Thus the result of this discrete form of the
- * algorithm is also exact.
- *
- * RandomLib provides two classes to sample from this distribution:
- * - DiscreteNormal which is tuned for speed on a typical general purpose
- * computer. This assumes that random samples can be generated relatively
- * quickly.
- * - DiscreteNormalAlt, which is a prototype for what might be needed on a
- * small device used for cryptography which is using a hardware generator
- * for obtaining truly random bits. This assumption here is that the
- * random bits are relatively expensive to obtain.
- * .
-
- * The basic algorithm is the same in the two cases. The main advantages of
- * this method are:
- * - exact sampling (provided that the source of random numbers is ideal),
- * - no need to cut off the tails of the distribution,
- * - a short program involving simple integer operations only,
- * - no dependence on external libraries (except to generate random bits),
- * - no large tables of constants needed,
- * - minimal time to set up for a new σ and μ (roughly comparable to
- * the time it takes to generate one sample),
- * - only about 5–20 times slower than standard routines to sample from
- * a normal distribution using plain double-precision arithmetic.
- * - DiscreteNormalAlt exhibits ideal scaling for the consumption of random
- * bits, namely a constant + log2σ, for large σ,
- * where the constant is about 31.
- * .
- * The possible drawbacks of this method are:
- * - σ and μ are restricted to rational numbers with sufficiently
- * small numerators and denominators to avoid overflow (this is unlikely to
- * be a severe restriction especially if the template parameter IntType is
- * set to long long),
- * - the running time is unbounded (but not in any practical sense),
- * - the memory consumption is unbounded (but not in any practical sense),
- * - the toll, about 30 bits, is considerably worse than that obtained using
- * the Knuth-Yao algorithm, for which the toll is no more than 2 (but this
- * requires a large table which is expensive to compute and requires a lot
- * of memory to store).
- *
- * This class uses a mutable private vector. So a single DiscreteNormal
- * object cannot safely be used by multiple threads. In a multi-processing
- * environment, each thread should use a thread-specific DiscreteNormal
- * object.
- *
- * Some timing results for IntType = int, μ = 0, and 108
- * samples (time = time per sample, including setup time, rv = mean number of
- * random variables per sample)
- * - σ = 10, time = 219 ns, rv = 17.52
- * - σ = 32, time = 223 ns, rv = 17.82
- * - σ = 1000, time = 225 ns, rv = 17.95
- * - σ = 160000, time = 226 ns, rv = 17.95
- *
- * @tparam IntType the integer type to use (default int).
- **********************************************************************/
- template class DiscreteNormal {
- public:
- /**
- * Constructor.
- *
- * @param[in] sigma_num the numerator of σ.
- * @param[in] sigma_den the denominator of σ (default 1).
- * @param[in] mu_num the numerator of μ (default 0).
- * @param[in] mu_den the denominator of μ (default 1).
- *
- * The constructor creates a DiscreteNormal objects for sampling with
- * specific values of σ and μ. This may throw an exception if the
- * parameters are such that overflow is possible. Internally σ and
- * μ are expressed with a common denominator, so it may be possible to
- * avoid overflow by picking the fractions of these quantities so that \e
- * sigma_den and \e mu_den have many common factors.
- **********************************************************************/
- DiscreteNormal(IntType sigma_num, IntType sigma_den = 1,
- IntType mu_num = 0, IntType mu_den = 1);
- /**
- * Return a sample.
- *
- * @tparam Random the type of the random generator.
- * @param[in,out] r a random generator.
- * @return discrete normal integer.
- **********************************************************************/
- template
- IntType operator()(Random& r) const;
- private:
- /**
- * sigma = _sig / _d, mu = _imu + _mu / _d, _isig = floor(sigma)
- **********************************************************************/
- IntType _sig, _mu, _d, _isig, _imu;
- typedef unsigned short word;
- /**
- * Holds as much of intermediate uniform deviates as needed.
- **********************************************************************/
- mutable std::vector _v;
- mutable unsigned _m, _l;
- /**
- * Increment on size of _v.
- **********************************************************************/
- static const unsigned alloc_incr = 16;
-
- // ceil(n/d) for d > 0
- static IntType iceil(IntType n, IntType d);
- // abs(n) needed because Visual Studio's std::abs has problems
- static IntType iabs(IntType n);
- static IntType gcd(IntType u, IntType v);
-
- // After x = LeadingDigit(p), p/_sig = (x + p'/_sig)/b where p and p' are
- // in [0, _sig) and b = 1 + max(word).
- word LeadingDigit(IntType& p) const;
-
- /**
- * Implement outcomes for choosing with prob (\e x + 2\e k) / (2\e k + 2);
- * return:
- * - 1 (succeed unconditionally) with prob (\e m − 2) / \e m,
- * - 0 (succeed with probability x) with prob 1 / \e m,
- * - −1 (fail unconditionally) with prob 1 / \e m.
- **********************************************************************/
- template static int Choose(Random& r, int m);
-
- // Compute v' < v. If true set v = v'.
- template bool less_than(Random& r) const;
-
- // Compute v < (x + p/_sig)/base (updating v)
- template bool less_than(Random& r, word x, IntType p) const;
-
- // true with prob (x + p/_sig)/base
- template bool bernoulli(Random& r, word x, IntType p) const;
-
- /**
- * Return true with probability exp(−1/2).
- **********************************************************************/
- template bool ExpProbH(Random& r) const;
-
- /**
- * Return true with probability exp(−n/2).
- **********************************************************************/
- template bool ExpProb(Random& r, int n) const;
-
- /**
- * Return \e n with probability exp(−n/2)
- * (1−exp(−1/2)).
- **********************************************************************/
- template int ExpProbN(Random& r) const;
-
- /**
- * Algorithm B: true with prob exp(-x * (2*k + x) / (2*k + 2)) where
- * x = (x0 + xn / _sig)/b.
- **********************************************************************/
- template
- bool B(Random& r, int k, word x0, IntType xn) const;
- };
-
- template DiscreteNormal::DiscreteNormal
- (IntType sigma_num, IntType sigma_den,
- IntType mu_num, IntType mu_den)
- : _v(std::vector(alloc_incr)), _m(0), _l(alloc_incr) {
- STATIC_ASSERT(std::numeric_limits::is_integer,
- "DiscreteNormal: invalid integer type IntType");
- STATIC_ASSERT(std::numeric_limits::is_signed,
- "DiscreteNormal: IntType must be a signed type");
- STATIC_ASSERT(!std::numeric_limits::is_signed,
- "DiscreteNormal: word must be an unsigned type");
- STATIC_ASSERT(std::numeric_limits::digits + 1 >=
- std::numeric_limits::digits,
- "DiscreteNormal: IntType must be at least as wide as word");
- if (!( sigma_num > 0 && sigma_den > 0 && mu_den > 0 ))
- throw RandomErr("DiscreteNormal: need sigma > 0");
- _imu = mu_num / mu_den;
- if (_imu == (std::numeric_limits::min)())
- throw RandomErr("DiscreteNormal: abs(mu) too large");
- mu_num -= _imu * mu_den;
- IntType l;
- l = gcd(sigma_num, sigma_den); sigma_num /= l; sigma_den /= l;
- l = gcd(mu_num, mu_den); mu_num /= l; mu_den /= l;
- _isig = iceil(sigma_num, sigma_den);
- l = gcd(sigma_den, mu_den);
- _sig = sigma_num * (mu_den / l);
- _mu = mu_num * (sigma_den / l);
- _d = sigma_den * (mu_den / l);
- // The rest of the constructor tests for possible overflow
- // Check for overflow in computing member variables
- IntType maxint = (std::numeric_limits::max)();
- if (!( mu_den / l <= maxint / sigma_num &&
- mu_num <= maxint / (sigma_den / l) &&
- mu_den / l <= maxint / sigma_den ))
- throw RandomErr("DiscreteNormal: sigma or mu overflow");
- // The probability that k = kmax is about 10^-543.
- int kmax = 50;
- // Check that max plausible result fits in an IntType, i.e.,
- // _isig * (kmax + 1) + abs(_imu) does not lead to overflow.
- if (!( kmax + 1 <= maxint / _isig &&
- _isig * (kmax + 1) <= maxint - iabs(_imu) ))
- throw RandomErr("DiscreteNormal: possible overflow a");
- // Check xn0 = _sig * k + s * _mu;
- if (!( kmax <= maxint / _sig &&
- _sig * kmax <= maxint - iabs(_mu) ))
- throw RandomErr("DiscreteNormal: possible overflow b");
- // Check for overflow in LeadingDigit
- // p << bits, p = _sig - 1, bits = 8
- if (!( _sig <= (maxint >> 8) ))
- throw RandomErr("DiscreteNormal: overflow in LeadingDigit");
- }
-
- template template
- IntType DiscreteNormal::operator()(Random& r) const {
- for (;;) {
- int k = ExpProbN(r);
- if (!ExpProb(r, k * (k - 1))) continue;
- IntType
- s = r.Boolean() ? -1 : 1,
- xn = _sig * IntType(k) + s * _mu,
- i = iceil(xn, _d) + r.template Integer(_isig);
- xn = i * _d - xn;
- if (xn >= _sig || (k == 0 && s < 0 && xn <= 0)) continue;
- if (xn > 0) {
- word x0 = LeadingDigit(xn); // Find first digit in expansion in words
- int h = k + 1; while (h-- && B(r, k, x0, xn));
- if (!(h < 0)) continue;
- }
- return s * i + _imu;
- }
- }
-
- template
- IntType DiscreteNormal::iceil(IntType n, IntType d)
- { IntType k = n / d; return k + (k * d < n ? 1 : 0); }
-
- template IntType DiscreteNormal::iabs(IntType n)
- { return n < 0 ? -n : n; }
-
- template
- IntType DiscreteNormal::gcd(IntType u, IntType v) {
- // Knuth, TAOCP, vol 2, 4.5.2, Algorithm A
- u = iabs(u); v = iabs(v);
- while (v > 0) { IntType r = u % v; u = v; v = r; }
- return u;
- }
-
- template typename DiscreteNormal::word
- DiscreteNormal::LeadingDigit(IntType& p) const {
- static const unsigned bits = 8;
- static const unsigned num = std::numeric_limits::digits / bits;
- STATIC_ASSERT(bits * num == std::numeric_limits::digits,
- "Number of digits in word must be multiple of 8");
- word s = 0;
- for (unsigned c = num; c--;) {
- p <<= bits; s <<= bits;
- word d = word(p / _sig);
- s += d;
- p -= IntType(d) * _sig;
- }
- return s;
- }
-
- template template
- int DiscreteNormal::Choose(Random& r, int m) {
- int k = r.template Integer(m);
- return k == 0 ? 0 : (k == 1 ? -1 : 1);
- }
-
- template template
- bool DiscreteNormal::less_than(Random& r) const {
- for (unsigned j = 0; ; ++j) {
- if (j == _m) {
- // Need more bits in the old V
- if (_l == _m) _v.resize(_l += alloc_incr);
- _v[_m++] = r.template Integer();
- }
- word w = r.template Integer();
- if (w > _v[j])
- return false; // New V is bigger, so exit
- else if (w < _v[j]) {
- _v[j] = w; // New V is smaller, update _v
- _m = j + 1; // adjusting its size
- return true; // and generate the next V
- }
- // Else w == _v[j] and we need to check the next word
- }
- }
-
- template template
- bool DiscreteNormal::less_than(Random& r, word x, IntType p) const {
- if (_m == 0) _v[_m++] = r.template Integer();
- if (_v[0] != x) return _v[0] < x;
- for (unsigned j = 1; ; ++j) {
- if (p == 0) return false;
- if (j == _m) {
- if (_l == _m) _v.resize(_l += alloc_incr);
- _v[_m++] = r.template Integer();
- }
- x = LeadingDigit(p);
- if (_v[j] != x) return _v[j] < x;
- }
- }
-
- template template
- bool DiscreteNormal::bernoulli(Random& r, word x, IntType p) const {
- word w = r.template Integer();
- if (w != x) return w < x;
- for (;;) {
- if (p == 0) return false;
- x = LeadingDigit(p);
- w = r.template Integer();
- if (w != x) return w < x;
- }
- }
-
- template template
- bool DiscreteNormal::ExpProbH(Random& r) const {
- static const word half = word(1) << (std::numeric_limits::digits - 1);
- _m = 0;
- if ((_v[_m++] = r.template Integer()) & half) return true;
- // Here _v < 1/2. Now loop finding decreasing V. Exit when first
- // increasing one is found.
- for (unsigned s = 0; ; s ^= 1) { // Parity of loop count
- if (!less_than(r)) return s != 0u;
- }
- }
-
- template template
- bool DiscreteNormal::ExpProb(Random& r, int n) const {
- while (n-- > 0) { if (!ExpProbH(r)) return false; }
- return true;
- }
-
- template template
- int DiscreteNormal::ExpProbN(Random& r) const {
- int n = 0;
- while (ExpProbH(r)) ++n;
- return n;
- }
-
- template template
- bool DiscreteNormal::B(Random& r, int k, word x0, IntType xn)
- const {
- int n = 0, h = 2 * k + 2, f;
- _m = 0;
- for (;; ++n) {
- if ( ((f = k ? 0 : Choose(r, h)) < 0) ||
- !(n ? less_than(r) : less_than(r, x0, xn)) ||
- ((f = k ? Choose(r, h) : f) < 0) ||
- (f == 0 && !bernoulli(r, x0, xn)) ) break;
- }
- return (n % 2) == 0;
- }
-
-} // namespace RandomLib
-
-#endif // RANDOMLIB_DISCRETENORMAL_HPP
diff --git a/include/RandomLib/DiscreteNormalAlt.hpp b/include/RandomLib/DiscreteNormalAlt.hpp
deleted file mode 100644
index b92e4974..00000000
--- a/include/RandomLib/DiscreteNormalAlt.hpp
+++ /dev/null
@@ -1,328 +0,0 @@
-/**
- * \file DiscreteNormalAlt.hpp
- * \brief Header for DiscreteNormalAlt
- *
- * Sample exactly from the discrete normal distribution (alternate version).
- *
- * Copyright (c) Charles Karney (2013)