Use PCG as random number generator

liquid_container
Lukas Werling 2016-04-20 20:46:55 +02:00
parent 250bbf9aa2
commit 7005eae55d
13 changed files with 3379 additions and 43 deletions

201
licenses/apache.txt 100644
View File

@ -0,0 +1,201 @@
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.

View File

@ -34,7 +34,7 @@ DWORD GenerateRandomPlayerColor(int32_t iTry) // generate a random player color
{
// generate a random one biased towards max channel luminance
// (for greater color difference and less gray-ish colors)
return C4RGB(std::min(SafeRandom(302), 256), std::min(SafeRandom(302), 256), std::min(SafeRandom(302), 256));
return C4RGB(std::min<int>(SafeRandom(302), 256), std::min<int>(SafeRandom(302), 256), std::min<int>(SafeRandom(302), 256));
}
bool IsColorConflict(DWORD dwClr1, DWORD dwClr2) // return whether dwClr1 and dwClr2 are closely together

View File

@ -140,8 +140,8 @@ struct C4RCMassMover
struct C4RCRandom
{
int Cnt; // index in seed
int Range; // random range query
int Val; // random value
uint32_t Range; // random range query
uint32_t Val; // random value
};
struct C4RCCreateObj

View File

@ -3095,7 +3095,7 @@ C4Player *C4Game::JoinPlayer(const char *szFilename, int32_t iAtClient, const ch
return pPlr;
}
void C4Game::FixRandom(int32_t iSeed)
void C4Game::FixRandom(uint64_t iSeed)
{
FixedRandom(iSeed);
}

View File

@ -148,7 +148,7 @@ public:
bool DoKeyboardInput(C4KeyCode vk_code, C4KeyEventType eEventType, bool fAlt, bool fCtrl, bool fShift, bool fRepeated, class C4GUI::Dialog *pForDialog=NULL, bool fPlrCtrlOnly=false, int32_t iStrength=-1);
bool DoKeyboardInput(C4KeyCodeEx Key, C4KeyEventType eEventType, class C4GUI::Dialog *pForDialog=NULL, bool fPlrCtrlOnly=false, int32_t iStrength=-1);
void DrawCrewOverheadText(C4TargetFacet &cgo, int32_t iPlayer);
void FixRandom(int32_t iSeed);
void FixRandom(uint64_t iSeed);
bool Init();
bool PreInit();
void SetScenarioFilename(const char*);

View File

@ -441,8 +441,8 @@ static std::vector<int32_t> GetRoundPolygon(int32_t x, int32_t y, int32_t size,
vertices.reserve(count * 2);
// varying phase of the sin/cos waves
C4Real begin = itofix(360)*Random(100) / 100;
C4Real begin2 = itofix(360)*Random(100) / 100;
C4Real begin = itofix(360)*(int32_t)Random(100) / 100;
C4Real begin2 = itofix(360)*(int32_t)Random(100) / 100;
// parameters:
// the bigger the factor, the smaller the divergence from a circle

View File

@ -393,7 +393,7 @@ void C4ParticleValueProvider::RollRandom(const C4Particle *forParticle)
void C4ParticleValueProvider::RollRandomUnseeded()
{
float range = endValue - startValue;
float rnd = (float)(rand()) / (float)(RAND_MAX);
float rnd = (float)(SafeRandom(RAND_MAX)) / (float)(RAND_MAX);
currentValue = startValue + rnd * range;
}

View File

@ -40,7 +40,7 @@ void C4SVal::Set(int32_t std, int32_t rnd, int32_t min, int32_t max)
int32_t C4SVal::Evaluate()
{
return Clamp(Std+Random(2*Rnd+1)-Rnd,Min,Max);
return Clamp<int32_t>(Std+Random(2*Rnd+1)-Rnd,Min,Max);
}
void C4SVal::Default()

View File

@ -22,40 +22,41 @@
#include "control/C4Record.h"
int RandomCount = 0;
static unsigned int RandomHold = 0;
static pcg32 RandomRng;
pcg32 SafeRandom;
void FixedRandom(DWORD dwSeed)
void FixedRandom(uint64_t seed)
{
// for SafeRandom
srand((unsigned)time(NULL));
RandomHold = dwSeed;
SafeRandom.seed(seed);
RandomRng.seed(seed);
RandomCount = 0;
}
int Random(int iRange)
static void RecordRandom(uint32_t range, uint32_t val)
{
RandomCount++;
if (Config.General.DebugRec)
{
// next pseudorandom value
RandomCount++;
C4RCRandom rc;
rc.Cnt=RandomCount;
rc.Range=iRange;
if (iRange<=0)
rc.Val=0;
else
{
RandomHold = ((uint64_t)RandomHold * 16807) % 2147483647;
rc.Val = RandomHold % iRange;
}
rc.Range=range;
rc.Val=val;
AddDbgRec(RCT_Random, &rc, sizeof(rc));
return rc.Val;
}
else
{
RandomCount++;
if (iRange<=0) return 0;
RandomHold = ((uint64_t)RandomHold * 16807) % 2147483647;
return RandomHold % iRange;
}
}
uint32_t Random()
{
uint32_t result = RandomRng();
RecordRandom(UINT32_MAX, result);
return result;
}
uint32_t Random(uint32_t iRange)
{
uint32_t result = RandomRng(iRange);
RecordRandom(iRange, result);
return result;
}

View File

@ -20,24 +20,20 @@
#ifndef INC_C4Random
#define INC_C4Random
#include <pcg/pcg_random.hpp>
extern int RandomCount;
extern pcg32 SafeRandom;
void FixedRandom(DWORD dwSeed);
void FixedRandom(uint64_t dwSeed);
int Random(int iRange);
uint32_t Random(uint32_t iRange);
inline unsigned int SeededRandom(unsigned int iSeed, unsigned int iRange)
inline uint32_t SeededRandom(uint64_t iSeed, uint32_t iRange)
{
if (!iRange) return 0;
iSeed = iSeed * 214013L + 2531011L;
return (iSeed >> 16) % iRange;
}
inline int SafeRandom(int range)
{
if (!range) return 0;
return rand()%range;
pcg32 rng(iSeed);
return rng(iRange);
}
#endif // INC_C4Random

637
thirdparty/pcg/pcg_extras.hpp vendored 100644
View File

@ -0,0 +1,637 @@
/*
* PCG Random Number Generation for C++
*
* Copyright 2014 Melissa O'Neill <oneill@pcg-random.org>
*
* 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.
*
* For additional information about the PCG random number generation scheme,
* including its license and other licensing options, visit
*
* http://www.pcg-random.org
*/
/*
* This file provides support code that is useful for random-number generation
* but not specific to the PCG generation scheme, including:
* - 128-bit int support for platforms where it isn't available natively
* - bit twiddling operations
* - I/O of 128-bit and 8-bit integers
* - Handling the evilness of SeedSeq
* - Support for efficiently producing random numbers less than a given
* bound
*/
#ifndef PCG_EXTRAS_HPP_INCLUDED
#define PCG_EXTRAS_HPP_INCLUDED 1
#include <cinttypes>
#include <cstddef>
#include <cstdlib>
#include <cstring>
#include <cassert>
#include <limits>
#include <iostream>
#include <type_traits>
#include <utility>
#include <locale>
#include <iterator>
#include <utility>
#ifdef __GNUC__
#include <cxxabi.h>
#endif
/*
* Abstractions for compiler-specific directives
*/
#ifdef __GNUC__
#define PCG_NOINLINE __attribute__((noinline))
#else
#define PCG_NOINLINE
#endif
/*
* Some members of the PCG library use 128-bit math. When compiling on 64-bit
* platforms, both GCC and Clang provide 128-bit integer types that are ideal
* for the job.
*
* On 32-bit platforms (or with other compilers), we fall back to a C++
* class that provides 128-bit unsigned integers instead. It may seem
* like we're reinventing the wheel here, because libraries already exist
* that support large integers, but most existing libraries provide a very
* generic multiprecision code, but here we're operating at a fixed size.
* Also, most other libraries are fairly heavyweight. So we use a direct
* implementation. Sadly, it's much slower than hand-coded assembly or
* direct CPU support.
*
*/
#if __SIZEOF_INT128__
namespace pcg_extras {
typedef __uint128_t pcg128_t;
}
#define PCG_128BIT_CONSTANT(high,low) \
((pcg128_t(high) << 64) + low)
#else
#include "pcg_uint128.hpp"
namespace pcg_extras {
typedef pcg_extras::uint_x4<uint32_t,uint64_t> pcg128_t;
}
#define PCG_128BIT_CONSTANT(high,low) \
pcg128_t(high,low)
#define PCG_EMULATED_128BIT_MATH 1
#endif
namespace pcg_extras {
/*
* We often need to represent a "number of bits". When used normally, these
* numbers are never greater than 128, so an unsigned char is plenty.
* If you're using a nonstandard generator of a larger size, you can set
* PCG_BITCOUNT_T to have it define it as a larger size. (Some compilers
* might produce faster code if you set it to an unsigned int.)
*/
#ifndef PCG_BITCOUNT_T
typedef uint8_t bitcount_t;
#else
typedef PCG_BITCOUNT_T bitcount_t;
#endif
/*
* C++ requires us to be able to serialize RNG state by printing or reading
* it from a stream. Because we use 128-bit ints, we also need to be able
* ot print them, so here is code to do so.
*
* This code provides enough functionality to print 128-bit ints in decimal
* and zero-padded in hex. It's not a full-featured implementation.
*/
template <typename CharT, typename Traits>
std::basic_ostream<CharT,Traits>&
operator<<(std::basic_ostream<CharT,Traits>& out, pcg128_t value)
{
auto desired_base = out.flags() & out.basefield;
bool want_hex = desired_base == out.hex;
if (want_hex) {
uint64_t highpart = uint64_t(value >> 64);
uint64_t lowpart = uint64_t(value);
auto desired_width = out.width();
if (desired_width > 16) {
out.width(desired_width - 16);
}
if (highpart != 0 || desired_width > 16)
out << highpart;
CharT oldfill;
if (highpart != 0) {
out.width(16);
oldfill = out.fill('0');
}
auto oldflags = out.setf(decltype(desired_base){}, out.showbase);
out << lowpart;
out.setf(oldflags);
if (highpart != 0) {
out.fill(oldfill);
}
return out;
}
constexpr size_t MAX_CHARS_128BIT = 40;
char buffer[MAX_CHARS_128BIT];
char* pos = buffer+sizeof(buffer);
*(--pos) = '\0';
constexpr auto BASE = pcg128_t(10ULL);
do {
auto div = value / BASE;
auto mod = uint32_t(value - (div * BASE));
*(--pos) = '0' + mod;
value = div;
} while(value != pcg128_t(0ULL));
return out << pos;
}
template <typename CharT, typename Traits>
std::basic_istream<CharT,Traits>&
operator>>(std::basic_istream<CharT,Traits>& in, pcg128_t& value)
{
typename std::basic_istream<CharT,Traits>::sentry s(in);
if (!s)
return in;
constexpr auto BASE = pcg128_t(10ULL);
pcg128_t current(0ULL);
bool did_nothing = true;
bool overflow = false;
for(;;) {
CharT wide_ch = in.get();
if (!in.good())
break;
auto ch = in.narrow(wide_ch, '\0');
if (ch < '0' || ch > '9') {
in.unget();
break;
}
did_nothing = false;
pcg128_t digit(uint32_t(ch - '0'));
pcg128_t timesbase = current*BASE;
overflow = overflow || timesbase < current;
current = timesbase + digit;
overflow = overflow || current < digit;
}
if (did_nothing || overflow) {
in.setstate(std::ios::failbit);
if (overflow)
current = ~pcg128_t(0ULL);
}
value = current;
return in;
}
/*
* Likewise, if people use tiny rngs, we'll be serializing uint8_t.
* If we just used the provided IO operators, they'd read/write chars,
* not ints, so we need to define our own. We *can* redefine this operator
* here because we're in our own namespace.
*/
template <typename CharT, typename Traits>
std::basic_ostream<CharT,Traits>&
operator<<(std::basic_ostream<CharT,Traits>&out, uint8_t value)
{
return out << uint32_t(value);
}
template <typename CharT, typename Traits>
std::basic_istream<CharT,Traits>&
operator>>(std::basic_istream<CharT,Traits>& in, uint8_t target)
{
uint32_t value = 0xdecea5edU;
in >> value;
if (!in && value == 0xdecea5edU)
return in;
if (value > uint8_t(~0)) {
in.setstate(std::ios::failbit);
value = ~0U;
}
target = uint8_t(value);
return in;
}
/* Unfortunately, the above functions don't get found in preference to the
* built in ones, so we create some more specific overloads that will.
* Ugh.
*/
inline std::ostream& operator<<(std::ostream& out, uint8_t value)
{
return pcg_extras::operator<< <char>(out, value);
}
inline std::istream& operator>>(std::istream& in, uint8_t& value)
{
return pcg_extras::operator>> <char>(in, value);
}
/*
* Useful bitwise operations.
*/
/*
* XorShifts are invertable, but they are someting of a pain to invert.
* This function backs them out. It's used by the whacky "inside out"
* generator defined later.
*/
template <typename itype>
inline itype unxorshift(itype x, bitcount_t bits, bitcount_t shift)
{
if (2*shift >= bits) {
return x ^ (x >> shift);
}
itype lowmask1 = (itype(1U) << (bits - shift*2)) - 1;
itype highmask1 = ~lowmask1;
itype top1 = x;
itype bottom1 = x & lowmask1;
top1 ^= top1 >> shift;
top1 &= highmask1;
x = top1 | bottom1;
itype lowmask2 = (itype(1U) << (bits - shift)) - 1;
itype bottom2 = x & lowmask2;
bottom2 = unxorshift(bottom2, bits - shift, shift);
bottom2 &= lowmask1;
return top1 | bottom2;
}
/*
* Rotate left and right.
*
* In ideal world, compilers would spot idiomatic rotate code and convert it
* to a rotate instruction. Of course, opinions vary on what the correct
* idiom is and how to spot it. For clang, sometimes it generates better
* (but still crappy) code if you define PCG_USE_ZEROCHECK_ROTATE_IDIOM.
*/
template <typename itype>
inline itype rotl(itype value, bitcount_t rot)
{
constexpr bitcount_t bits = sizeof(itype) * 8;
constexpr bitcount_t mask = bits - 1;
#if PCG_USE_ZEROCHECK_ROTATE_IDIOM
return rot ? (value << rot) | (value >> (bits - rot)) : value;
#else
return (value << rot) | (value >> ((- rot) & mask));
#endif
}
template <typename itype>
inline itype rotr(itype value, bitcount_t rot)
{
constexpr bitcount_t bits = sizeof(itype) * 8;
constexpr bitcount_t mask = bits - 1;
#if PCG_USE_ZEROCHECK_ROTATE_IDIOM
return rot ? (value >> rot) | (value << (bits - rot)) : value;
#else
return (value >> rot) | (value << ((- rot) & mask));
#endif
}
/* Unfortunately, both Clang and GCC sometimes perform poorly when it comes
* to properly recognizing idiomatic rotate code, so for we also provide
* assembler directives (enabled with PCG_USE_INLINE_ASM). Boo, hiss.
* (I hope that these compilers get better so that this code can die.)
*
* These overloads will be preferred over the general template code above.
*/
#if PCG_USE_INLINE_ASM && __GNUC__ && (__x86_64__ || __i386__)
inline uint8_t rotr(uint8_t value, bitcount_t rot)
{
asm ("rorb %%cl, %0" : "=r" (value) : "0" (value), "c" (rot));
return value;
}
inline uint16_t rotr(uint16_t value, bitcount_t rot)
{
asm ("rorw %%cl, %0" : "=r" (value) : "0" (value), "c" (rot));
return value;
}
inline uint32_t rotr(uint32_t value, bitcount_t rot)
{
asm ("rorl %%cl, %0" : "=r" (value) : "0" (value), "c" (rot));
return value;
}
#if __x86_64__
inline uint64_t rotr(uint64_t value, bitcount_t rot)
{
asm ("rorq %%cl, %0" : "=r" (value) : "0" (value), "c" (rot));
return value;
}
#endif // __x86_64__
#endif // PCG_USE_INLINE_ASM
/*
* The C++ SeedSeq concept (modelled by seed_seq) can fill an array of
* 32-bit integers with seed data, but sometimes we want to produce
* larger or smaller integers.
*
* The following code handles this annoyance.
*
* uneven_copy will copy an array of 32-bit ints to an array of larger or
* smaller ints (actually, the code is general it only needing forward
* iterators). The copy is identical to the one that would be performed if
* we just did memcpy on a standard little-endian machine, but works
* regardless of the endian of the machine (or the weirdness of the ints
* involved).
*
* generate_to initializes an array of integers using a SeedSeq
* object. It is given the size as a static constant at compile time and
* tries to avoid memory allocation. If we're filling in 32-bit constants
* we just do it directly. If we need a separate buffer and it's small,
* we allocate it on the stack. Otherwise, we fall back to heap allocation.
* Ugh.
*
* generate_one produces a single value of some integral type using a
* SeedSeq object.
*/
/* uneven_copy helper, case where destination ints are less than 32 bit. */
template<class SrcIter, class DestIter>
SrcIter uneven_copy_impl(
SrcIter src_first, DestIter dest_first, DestIter dest_last,
std::true_type)
{
typedef typename std::iterator_traits<SrcIter>::value_type src_t;
typedef typename std::iterator_traits<DestIter>::value_type dest_t;
constexpr bitcount_t SRC_SIZE = sizeof(src_t);
constexpr bitcount_t DEST_SIZE = sizeof(dest_t);
constexpr bitcount_t DEST_BITS = DEST_SIZE * 8;
constexpr bitcount_t SCALE = SRC_SIZE / DEST_SIZE;
size_t count = 0;
src_t value;
while (dest_first != dest_last) {
if ((count++ % SCALE) == 0)
value = *src_first++; // Get more bits
else
value >>= DEST_BITS; // Move down bits
*dest_first++ = dest_t(value); // Truncates, ignores high bits.
}
return src_first;
}
/* uneven_copy helper, case where destination ints are more than 32 bit. */
template<class SrcIter, class DestIter>
SrcIter uneven_copy_impl(
SrcIter src_first, DestIter dest_first, DestIter dest_last,
std::false_type)
{
typedef typename std::iterator_traits<SrcIter>::value_type src_t;
typedef typename std::iterator_traits<DestIter>::value_type dest_t;
constexpr auto SRC_SIZE = sizeof(src_t);
constexpr auto SRC_BITS = SRC_SIZE * 8;
constexpr auto DEST_SIZE = sizeof(dest_t);
constexpr auto SCALE = (DEST_SIZE+SRC_SIZE-1) / SRC_SIZE;
while (dest_first != dest_last) {
dest_t value(0UL);
unsigned int shift = 0;
for (size_t i = 0; i < SCALE; ++i) {
value |= dest_t(*src_first++) << shift;
shift += SRC_BITS;
}
*dest_first++ = value;
}
return src_first;
}
/* uneven_copy, call the right code for larger vs. smaller */
template<class SrcIter, class DestIter>
inline SrcIter uneven_copy(SrcIter src_first,
DestIter dest_first, DestIter dest_last)
{
typedef typename std::iterator_traits<SrcIter>::value_type src_t;
typedef typename std::iterator_traits<DestIter>::value_type dest_t;
constexpr bool DEST_IS_SMALLER = sizeof(dest_t) < sizeof(src_t);
return uneven_copy_impl(src_first, dest_first, dest_last,
std::integral_constant<bool, DEST_IS_SMALLER>{});
}
/* generate_to, fill in a fixed-size array of integral type using a SeedSeq
* (actually works for any random-access iterator)
*/
template <size_t size, typename SeedSeq, typename DestIter>
inline void generate_to_impl(SeedSeq&& generator, DestIter dest,
std::true_type)
{
generator.generate(dest, dest+size);
}
template <size_t size, typename SeedSeq, typename DestIter>
void generate_to_impl(SeedSeq&& generator, DestIter dest,
std::false_type)
{
typedef typename std::iterator_traits<DestIter>::value_type dest_t;
constexpr auto DEST_SIZE = sizeof(dest_t);
constexpr auto GEN_SIZE = sizeof(uint32_t);
constexpr bool GEN_IS_SMALLER = GEN_SIZE < DEST_SIZE;
constexpr size_t FROM_ELEMS =
GEN_IS_SMALLER
? size * ((DEST_SIZE+GEN_SIZE-1) / GEN_SIZE)
: (size + (GEN_SIZE / DEST_SIZE) - 1)
/ ((GEN_SIZE / DEST_SIZE) + GEN_IS_SMALLER);
// this odd code ^^^^^^^^^^^^^^^^^ is work-around for
// a bug: http://llvm.org/bugs/show_bug.cgi?id=21287
if (FROM_ELEMS <= 1024) {
uint32_t buffer[FROM_ELEMS];
generator.generate(buffer, buffer+FROM_ELEMS);
uneven_copy(buffer, dest, dest+size);
} else {
uint32_t* buffer = (uint32_t*) malloc(GEN_SIZE * FROM_ELEMS);
generator.generate(buffer, buffer+FROM_ELEMS);
uneven_copy(buffer, dest, dest+size);
free(buffer);
}
}
template <size_t size, typename SeedSeq, typename DestIter>
inline void generate_to(SeedSeq&& generator, DestIter dest)
{
typedef typename std::iterator_traits<DestIter>::value_type dest_t;
constexpr bool IS_32BIT = sizeof(dest_t) == sizeof(uint32_t);
generate_to_impl<size>(std::forward<SeedSeq>(generator), dest,
std::integral_constant<bool, IS_32BIT>{});
}
/* generate_one, produce a value of integral type using a SeedSeq
* (optionally, we can have it produce more than one and pick which one
* we want)
*/
template <typename UInt, size_t i = 0UL, size_t N = i+1UL, typename SeedSeq>
inline UInt generate_one(SeedSeq&& generator)
{
UInt result[N];
generate_to<N>(std::forward<SeedSeq>(generator), result);
return result[i];
}
template <typename RngType>
auto bounded_rand(RngType& rng, typename RngType::result_type upper_bound)
-> typename RngType::result_type
{
typedef typename RngType::result_type rtype;
rtype threshold = (RngType::max() - RngType::min() + rtype(1) - upper_bound)
% upper_bound;
for (;;) {
rtype r = rng() - RngType::min();
if (r >= threshold)
return r % upper_bound;
}
}
template <typename Iter, typename RandType>
void shuffle(Iter from, Iter to, RandType&& rng)
{
typedef typename std::iterator_traits<Iter>::difference_type delta_t;
auto count = to - from;
while (count > 1) {
delta_t chosen(bounded_rand(rng, count));
--count;
--to;
using std::swap;
swap(*(from+chosen), *to);
}
}
/*
* Although std::seed_seq is useful, it isn't everything. Often we want to
* initialize a random-number generator some other way, such as from a random
* device.
*
* Technically, it does not meet the requirements of a SeedSequence because
* it lacks some of the rarely-used member functions (some of which would
* be impossible to provide). However the C++ standard is quite specific
* that actual engines only called the generate method, so it ought not to be
* a problem in practice.
*/
template <typename RngType>
class seed_seq_from {
private:
RngType rng_;
typedef uint_least32_t result_type;
public:
template<typename... Args>
seed_seq_from(Args&&... args) :
rng_(std::forward<Args>(args)...)
{
// Nothing (else) to do...
}
template<typename Iter>
void generate(Iter start, Iter finish)
{
for (auto i = start; i != finish; ++i)
*i = result_type(rng_());
}
constexpr size_t size() const
{
return (sizeof(typename RngType::result_type) > sizeof(result_type)
&& RngType::max() > ~size_t(0UL))
? ~size_t(0UL)
: size_t(RngType::max());
}
};
/*
* Sometimes you might want a distinct seed based on when the program
* was compiled. That way, a particular instance of the program will
* behave the same way, but when recompiled it'll produce a different
* value.
*/
template <typename IntType>
struct static_arbitrary_seed {
private:
static constexpr IntType fnv(IntType hash, const char* pos) {
return *pos == '\0'
? hash
: fnv((hash * IntType(16777619U)) ^ *pos, (pos+1));
}
public:
static constexpr IntType value = fnv(IntType(2166136261U ^ sizeof(IntType)),
__DATE__ __TIME__ __FILE__);
};
// Sometimes, when debugging or testing, it's handy to be able print the name
// of a (in human-readable form). This code allows the idiom:
//
// cout << printable_typename<my_foo_type_t>()
//
// to print out my_foo_type_t (or its concrete type if it is a synonym)
template <typename T>
struct printable_typename {};
template <typename T>
std::ostream& operator<<(std::ostream& out, printable_typename<T>) {
const char *implementation_typename = typeid(T).name();
#ifdef __GNUC__
int status;
const char* pretty_name =
abi::__cxa_demangle(implementation_typename, NULL, NULL, &status);
if (status == 0)
out << pretty_name;
free((void*) pretty_name);
if (status == 0)
return out;
#endif
out << implementation_typename;
return out;
}
} // namespace pcg_extras
#endif // PCG_EXTRAS_HPP_INCLUDED

1751
thirdparty/pcg/pcg_random.hpp vendored 100644

File diff suppressed because it is too large Load Diff

750
thirdparty/pcg/pcg_uint128.hpp vendored 100644
View File

@ -0,0 +1,750 @@
/*
* PCG Random Number Generation for C++
*
* Copyright 2014 Melissa O'Neill <oneill@pcg-random.org>
*
* 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.
*
* For additional information about the PCG random number generation scheme,
* including its license and other licensing options, visit
*
* http://www.pcg-random.org
*/
/*
* This code provides a a C++ class that can provide 128-bit (or higher)
* integers. To produce 2K-bit integers, it uses two K-bit integers,
* placed in a union that allowes the code to also see them as four K/2 bit
* integers (and access them either directly name, or by index).
*
* It may seem like we're reinventing the wheel here, because several
* libraries already exist that support large integers, but most existing
* libraries provide a very generic multiprecision code, but here we're
* operating at a fixed size. Also, most other libraries are fairly
* heavyweight. So we use a direct implementation. Sadly, it's much slower
* than hand-coded assembly or direct CPU support.
*/
#ifndef PCG_UINT128_HPP_INCLUDED
#define PCG_UINT128_HPP_INCLUDED 1
#include <cstdint>
#include <cstdio>
#include <cassert>
#include <climits>
#include <utility>
#include <initializer_list>
#include <type_traits>
/*
* We want to lay the type out the same way that a native type would be laid
* out, which means we must know the machine's endian, at compile time.
* This ugliness attempts to do so.
*/
#ifndef PCG_LITTLE_ENDIAN
#if defined(__BYTE_ORDER__)
#if __BYTE_ORDER__ == __ORDER_LITTLE_ENDIAN__
#define PCG_LITTLE_ENDIAN 1
#elif __BYTE_ORDER__ == __ORDER_BIG_ENDIAN__
#define PCG_LITTLE_ENDIAN 0
#else
#error __BYTE_ORDER__ does not match a standard endian, pick a side
#endif
#elif __LITTLE_ENDIAN__ || _LITTLE_ENDIAN
#define PCG_LITTLE_ENDIAN 1
#elif __BIG_ENDIAN__ || _BIG_ENDIAN
#define PCG_LITTLE_ENDIAN 0
#elif __x86_64 || __x86_64__ || __i386 || __i386__
#define PCG_LITTLE_ENDIAN 1
#elif __powerpc__ || __POWERPC__ || __ppc__ || __PPC__ \
|| __m68k__ || __mc68000__
#define PCG_LITTLE_ENDIAN 0
#else
#error Unable to determine target endianness
#endif
#endif
namespace pcg_extras {
// Recent versions of GCC have intrinsics we can use to quickly calculate
// the number of leading and trailing zeros in a number. If possible, we
// use them, otherwise we fall back to old-fashioned bit twiddling to figure
// them out.
#ifndef PCG_BITCOUNT_T
typedef uint8_t bitcount_t;
#else
typedef PCG_BITCOUNT_T bitcount_t;
#endif
/*
* Provide some useful helper functions
* * flog2 floor(log2(x))
* * trailingzeros number of trailing zero bits
*/
#ifdef __GNUC__ // Any GNU-compatible compiler supporting C++11 has
// some useful intrinsics we can use.
inline bitcount_t flog2(uint32_t v)
{
return 31 - __builtin_clz(v);
}
inline bitcount_t trailingzeros(uint32_t v)
{
return __builtin_ctz(v);
}
inline bitcount_t flog2(uint64_t v)
{
#if UINT64_MAX == ULONG_MAX
return 63 - __builtin_clzl(v);
#elif UINT64_MAX == ULLONG_MAX
return 63 - __builtin_clzll(v);
#else
#error Cannot find a function for uint64_t
#endif
}
inline bitcount_t trailingzeros(uint64_t v)
{
#if UINT64_MAX == ULONG_MAX
return __builtin_ctzl(v);
#elif UINT64_MAX == ULLONG_MAX
return __builtin_ctzll(v);
#else
#error Cannot find a function for uint64_t
#endif
}
#else // Otherwise, we fall back to bit twiddling
// implementations
inline bitcount_t flog2(uint32_t v)
{
// Based on code by Eric Cole and Mark Dickinson, which appears at
// https://graphics.stanford.edu/~seander/bithacks.html#IntegerLogDeBruijn
static const uint8_t multiplyDeBruijnBitPos[32] = {
0, 9, 1, 10, 13, 21, 2, 29, 11, 14, 16, 18, 22, 25, 3, 30,
8, 12, 20, 28, 15, 17, 24, 7, 19, 27, 23, 6, 26, 5, 4, 31
};
v |= v >> 1; // first round down to one less than a power of 2
v |= v >> 2;
v |= v >> 4;
v |= v >> 8;
v |= v >> 16;
return multiplyDeBruijnBitPos[(uint32_t)(v * 0x07C4ACDDU) >> 27];
}
inline bitcount_t trailingzeros(uint32_t v)
{
static const uint8_t multiplyDeBruijnBitPos[32] = {
0, 1, 28, 2, 29, 14, 24, 3, 30, 22, 20, 15, 25, 17, 4, 8,
31, 27, 13, 23, 21, 19, 16, 7, 26, 12, 18, 6, 11, 5, 10, 9
};
return multiplyDeBruijnBitPos[((uint32_t)((v & -v) * 0x077CB531U)) >> 27];
}
inline bitcount_t flog2(uint64_t v)
{
uint32_t high = v >> 32;
uint32_t low = uint32_t(v);
return high ? 32+flog2(high) : flog2(low);
}
inline bitcount_t trailingzeros(uint64_t v)
{
uint32_t high = v >> 32;
uint32_t low = uint32_t(v);
return low ? trailingzeros(low) : trailingzeros(high)+32;
}
#endif
template <typename UInt>
inline bitcount_t clog2(UInt v)
{
return flog2(v) + ((v & (-v)) != v);
}
template <typename UInt>
inline UInt addwithcarry(UInt x, UInt y, bool carryin, bool* carryout)
{
UInt half_result = y + carryin;
UInt result = x + half_result;
*carryout = (half_result < y) || (result < x);
return result;
}
template <typename UInt>
inline UInt subwithcarry(UInt x, UInt y, bool carryin, bool* carryout)
{
UInt half_result = y + carryin;
UInt result = x - half_result;
*carryout = (half_result < y) || (result > x);
return result;
}
template <typename UInt, typename UIntX2>
class uint_x4 {
// private:
public:
union {
#if PCG_LITTLE_ENDIAN
struct {
UInt v0, v1, v2, v3;
} w;
struct {
UIntX2 v01, v23;
} d;
#else
struct {
UInt v3, v2, v1, v0;
} w;
struct {
UIntX2 v23, v01;
} d;
#endif
// For the array access versions, the code that uses the array
// must handle endian itself. Yuck.
UInt wa[4];
UIntX2 da[2];
};
public:
uint_x4() = default;
constexpr uint_x4(UInt v3, UInt v2, UInt v1, UInt v0)
#if PCG_LITTLE_ENDIAN
: w{v0, v1, v2, v3}
#else
: w{v3, v2, v1, v0}
#endif
{
// Nothing (else) to do
}
constexpr uint_x4(UIntX2 v23, UIntX2 v01)
#if PCG_LITTLE_ENDIAN
: d{v01,v23}
#else
: d{v23,v01}
#endif
{
// Nothing (else) to do
}
template<class Integral,
typename std::enable_if<(std::is_integral<Integral>::value
&& sizeof(Integral) <= sizeof(UIntX2))
>::type* = nullptr>
constexpr uint_x4(Integral v01)
#if PCG_LITTLE_ENDIAN
: d{UIntX2(v01),0UL}
#else
: d{0UL,UIntX2(v01)}
#endif
{
// Nothing (else) to do
}
explicit constexpr operator uint64_t() const
{
return d.v01;
}
explicit constexpr operator uint32_t() const
{
return w.v0;
}
explicit constexpr operator int() const
{
return w.v0;
}
explicit constexpr operator uint16_t() const
{
return w.v0;
}
explicit constexpr operator uint8_t() const
{
return w.v0;
}
typedef typename std::conditional<std::is_same<uint64_t,
unsigned long>::value,
unsigned long long,
unsigned long>::type
uint_missing_t;
explicit constexpr operator uint_missing_t() const
{
return d.v01;
}
explicit constexpr operator bool() const
{
return d.v01 || d.v23;
}
template<typename U, typename V>
friend uint_x4<U,V> operator*(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend std::pair< uint_x4<U,V>,uint_x4<U,V> >
divmod(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend uint_x4<U,V> operator+(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend uint_x4<U,V> operator-(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend uint_x4<U,V> operator<<(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend uint_x4<U,V> operator>>(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend uint_x4<U,V> operator&(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend uint_x4<U,V> operator|(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend uint_x4<U,V> operator^(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend bool operator==(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend bool operator!=(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend bool operator<(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend bool operator<=(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend bool operator>(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend bool operator>=(const uint_x4<U,V>&, const uint_x4<U,V>&);
template<typename U, typename V>
friend uint_x4<U,V> operator~(const uint_x4<U,V>&);
template<typename U, typename V>
friend uint_x4<U,V> operator-(const uint_x4<U,V>&);
template<typename U, typename V>
friend bitcount_t flog2(const uint_x4<U,V>&);
template<typename U, typename V>
friend bitcount_t trailingzeros(const uint_x4<U,V>&);
uint_x4& operator*=(const uint_x4& rhs)
{
uint_x4 result = *this * rhs;
return *this = result;
}
uint_x4& operator/=(const uint_x4& rhs)
{
uint_x4 result = *this / rhs;
return *this = result;
}
uint_x4& operator%=(const uint_x4& rhs)
{
uint_x4 result = *this % rhs;
return *this = result;
}
uint_x4& operator+=(const uint_x4& rhs)
{
uint_x4 result = *this + rhs;
return *this = result;
}
uint_x4& operator-=(const uint_x4& rhs)
{
uint_x4 result = *this - rhs;
return *this = result;
}
uint_x4& operator&=(const uint_x4& rhs)
{
uint_x4 result = *this & rhs;
return *this = result;
}
uint_x4& operator|=(const uint_x4& rhs)
{
uint_x4 result = *this | rhs;
return *this = result;
}
uint_x4& operator^=(const uint_x4& rhs)
{
uint_x4 result = *this ^ rhs;
return *this = result;
}
uint_x4& operator>>=(bitcount_t shift)
{
uint_x4 result = *this >> shift;
return *this = result;
}
uint_x4& operator<<=(bitcount_t shift)
{
uint_x4 result = *this << shift;
return *this = result;
}
};
template<typename U, typename V>
bitcount_t flog2(const uint_x4<U,V>& v)
{
#if PCG_LITTLE_ENDIAN
for (uint8_t i = 4; i !=0; /* dec in loop */) {
--i;
#else
for (uint8_t i = 0; i < 4; ++i) {
#endif
if (v.wa[i] == 0)
continue;
return flog2(v.wa[i]) + (sizeof(U)*CHAR_BIT)*i;
}
abort();
}
template<typename U, typename V>
bitcount_t trailingzeros(const uint_x4<U,V>& v)
{
#if PCG_LITTLE_ENDIAN
for (uint8_t i = 0; i < 4; ++i) {
#else
for (uint8_t i = 4; i !=0; /* dec in loop */) {
--i;
#endif
if (v.wa[i] != 0)
return trailingzeros(v.wa[i]) + (sizeof(U)*CHAR_BIT)*i;
}
return (sizeof(U)*CHAR_BIT)*4;
}
template <typename UInt, typename UIntX2>
std::pair< uint_x4<UInt,UIntX2>, uint_x4<UInt,UIntX2> >
divmod(const uint_x4<UInt,UIntX2>& orig_dividend,
const uint_x4<UInt,UIntX2>& divisor)
{
// If the dividend is less than the divisor, the answer is always zero.
// This takes care of boundary cases like 0/x (which would otherwise be
// problematic because we can't take the log of zero. (The boundary case
// of division by zero is undefined.)
if (orig_dividend < divisor)
return { uint_x4<UInt,UIntX2>(0UL), orig_dividend };
auto dividend = orig_dividend;
auto log2_divisor = flog2(divisor);
auto log2_dividend = flog2(dividend);
// assert(log2_dividend >= log2_divisor);
bitcount_t logdiff = log2_dividend - log2_divisor;
constexpr uint_x4<UInt,UIntX2> ONE(1UL);
if (logdiff == 0)
return { ONE, dividend - divisor };
// Now we change the log difference to
// floor(log2(divisor)) - ceil(log2(dividend))
// to ensure that we *underestimate* the result.
logdiff -= 1;
uint_x4<UInt,UIntX2> quotient(0UL);
auto qfactor = ONE << logdiff;
auto factor = divisor << logdiff;
do {
dividend -= factor;
quotient += qfactor;
while (dividend < factor) {
factor >>= 1;
qfactor >>= 1;
}
} while (dividend >= divisor);
return { quotient, dividend };
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator/(const uint_x4<UInt,UIntX2>& dividend,
const uint_x4<UInt,UIntX2>& divisor)
{
return divmod(dividend, divisor).first;
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator%(const uint_x4<UInt,UIntX2>& dividend,
const uint_x4<UInt,UIntX2>& divisor)
{
return divmod(dividend, divisor).second;
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator*(const uint_x4<UInt,UIntX2>& a,
const uint_x4<UInt,UIntX2>& b)
{
uint_x4<UInt,UIntX2> r = {0U, 0U, 0U, 0U};
bool carryin = false;
bool carryout;
UIntX2 a0b0 = UIntX2(a.w.v0) * UIntX2(b.w.v0);
r.w.v0 = UInt(a0b0);
r.w.v1 = UInt(a0b0 >> 32);
UIntX2 a1b0 = UIntX2(a.w.v1) * UIntX2(b.w.v0);
r.w.v2 = UInt(a1b0 >> 32);
r.w.v1 = addwithcarry(r.w.v1, UInt(a1b0), carryin, &carryout);
carryin = carryout;
r.w.v2 = addwithcarry(r.w.v2, UInt(0U), carryin, &carryout);
carryin = carryout;
r.w.v3 = addwithcarry(r.w.v3, UInt(0U), carryin, &carryout);
UIntX2 a0b1 = UIntX2(a.w.v0) * UIntX2(b.w.v1);
carryin = false;
r.w.v2 = addwithcarry(r.w.v2, UInt(a0b1 >> 32), carryin, &carryout);
carryin = carryout;
r.w.v3 = addwithcarry(r.w.v3, UInt(0U), carryin, &carryout);
carryin = false;
r.w.v1 = addwithcarry(r.w.v1, UInt(a0b1), carryin, &carryout);
carryin = carryout;
r.w.v2 = addwithcarry(r.w.v2, UInt(0U), carryin, &carryout);
carryin = carryout;
r.w.v3 = addwithcarry(r.w.v3, UInt(0U), carryin, &carryout);
UIntX2 a1b1 = UIntX2(a.w.v1) * UIntX2(b.w.v1);
carryin = false;
r.w.v2 = addwithcarry(r.w.v2, UInt(a1b1), carryin, &carryout);
carryin = carryout;
r.w.v3 = addwithcarry(r.w.v3, UInt(a1b1 >> 32), carryin, &carryout);
r.d.v23 += a.d.v01 * b.d.v23 + a.d.v23 * b.d.v01;
return r;
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator+(const uint_x4<UInt,UIntX2>& a,
const uint_x4<UInt,UIntX2>& b)
{
uint_x4<UInt,UIntX2> r = {0U, 0U, 0U, 0U};
bool carryin = false;
bool carryout;
r.w.v0 = addwithcarry(a.w.v0, b.w.v0, carryin, &carryout);
carryin = carryout;
r.w.v1 = addwithcarry(a.w.v1, b.w.v1, carryin, &carryout);
carryin = carryout;
r.w.v2 = addwithcarry(a.w.v2, b.w.v2, carryin, &carryout);
carryin = carryout;
r.w.v3 = addwithcarry(a.w.v3, b.w.v3, carryin, &carryout);
return r;
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator-(const uint_x4<UInt,UIntX2>& a,
const uint_x4<UInt,UIntX2>& b)
{
uint_x4<UInt,UIntX2> r = {0U, 0U, 0U, 0U};
bool carryin = false;
bool carryout;
r.w.v0 = subwithcarry(a.w.v0, b.w.v0, carryin, &carryout);
carryin = carryout;
r.w.v1 = subwithcarry(a.w.v1, b.w.v1, carryin, &carryout);
carryin = carryout;
r.w.v2 = subwithcarry(a.w.v2, b.w.v2, carryin, &carryout);
carryin = carryout;
r.w.v3 = subwithcarry(a.w.v3, b.w.v3, carryin, &carryout);
return r;
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator&(const uint_x4<UInt,UIntX2>& a,
const uint_x4<UInt,UIntX2>& b)
{
return uint_x4<UInt,UIntX2>(a.d.v23 & b.d.v23, a.d.v01 & b.d.v01);
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator|(const uint_x4<UInt,UIntX2>& a,
const uint_x4<UInt,UIntX2>& b)
{
return uint_x4<UInt,UIntX2>(a.d.v23 | b.d.v23, a.d.v01 | b.d.v01);
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator^(const uint_x4<UInt,UIntX2>& a,
const uint_x4<UInt,UIntX2>& b)
{
return uint_x4<UInt,UIntX2>(a.d.v23 ^ b.d.v23, a.d.v01 ^ b.d.v01);
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator~(const uint_x4<UInt,UIntX2>& v)
{
return uint_x4<UInt,UIntX2>(~v.d.v23, ~v.d.v01);
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator-(const uint_x4<UInt,UIntX2>& v)
{
return uint_x4<UInt,UIntX2>(0UL,0UL) - v;
}
template <typename UInt, typename UIntX2>
bool operator==(const uint_x4<UInt,UIntX2>& a, const uint_x4<UInt,UIntX2>& b)
{
return (a.d.v01 == b.d.v01) && (a.d.v23 == b.d.v23);
}
template <typename UInt, typename UIntX2>
bool operator!=(const uint_x4<UInt,UIntX2>& a, const uint_x4<UInt,UIntX2>& b)
{
return !operator==(a,b);
}
template <typename UInt, typename UIntX2>
bool operator<(const uint_x4<UInt,UIntX2>& a, const uint_x4<UInt,UIntX2>& b)
{
return (a.d.v23 < b.d.v23)
|| ((a.d.v23 == b.d.v23) && (a.d.v01 < b.d.v01));
}
template <typename UInt, typename UIntX2>
bool operator>(const uint_x4<UInt,UIntX2>& a, const uint_x4<UInt,UIntX2>& b)
{
return operator<(b,a);
}
template <typename UInt, typename UIntX2>
bool operator<=(const uint_x4<UInt,UIntX2>& a, const uint_x4<UInt,UIntX2>& b)
{
return !(operator<(b,a));
}
template <typename UInt, typename UIntX2>
bool operator>=(const uint_x4<UInt,UIntX2>& a, const uint_x4<UInt,UIntX2>& b)
{
return !(operator<(a,b));
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator<<(const uint_x4<UInt,UIntX2>& v,
const bitcount_t shift)
{
uint_x4<UInt,UIntX2> r = {0U, 0U, 0U, 0U};
const bitcount_t bits = sizeof(UInt) * CHAR_BIT;
const bitcount_t bitmask = bits - 1;
const bitcount_t shiftdiv = shift / bits;
const bitcount_t shiftmod = shift & bitmask;
if (shiftmod) {
UInt carryover = 0;
#if PCG_LITTLE_ENDIAN
for (uint8_t out = shiftdiv, in = 0; out < 4; ++out, ++in) {
#else
for (uint8_t out = 4-shiftdiv, in = 4; out != 0; /* dec in loop */) {
--out, --in;
#endif
r.wa[out] = (v.wa[in] << shiftmod) | carryover;
carryover = (v.wa[in] >> (bits - shiftmod));
}
} else {
#if PCG_LITTLE_ENDIAN
for (uint8_t out = shiftdiv, in = 0; out < 4; ++out, ++in) {
#else
for (uint8_t out = 4-shiftdiv, in = 4; out != 0; /* dec in loop */) {
--out, --in;
#endif
r.wa[out] = v.wa[in];
}
}
return r;
}
template <typename UInt, typename UIntX2>
uint_x4<UInt,UIntX2> operator>>(const uint_x4<UInt,UIntX2>& v,
const bitcount_t shift)
{
uint_x4<UInt,UIntX2> r = {0U, 0U, 0U, 0U};
const bitcount_t bits = sizeof(UInt) * CHAR_BIT;
const bitcount_t bitmask = bits - 1;
const bitcount_t shiftdiv = shift / bits;
const bitcount_t shiftmod = shift & bitmask;
if (shiftmod) {
UInt carryover = 0;
#if PCG_LITTLE_ENDIAN
for (uint8_t out = 4-shiftdiv, in = 4; out != 0; /* dec in loop */) {
--out, --in;
#else
for (uint8_t out = shiftdiv, in = 0; out < 4; ++out, ++in) {
#endif
r.wa[out] = (v.wa[in] >> shiftmod) | carryover;
carryover = (v.wa[in] << (bits - shiftmod));
}
} else {
#if PCG_LITTLE_ENDIAN
for (uint8_t out = 4-shiftdiv, in = 4; out != 0; /* dec in loop */) {
--out, --in;
#else
for (uint8_t out = shiftdiv, in = 0; out < 4; ++out, ++in) {
#endif
r.wa[out] = v.wa[in];
}
}
return r;
}
} // namespace pcg_extras
#endif // PCG_UINT128_HPP_INCLUDED