［前：性能評価］
［ExRand目次］
共通プログラム
目次
define.h
bincode.h
xgf2.h
xgf2.c
define.h
/* define.h : coded by Isaku WADA(isaku@pb4.so-net.ne.jp) 2025 */ #if !defined(DEFINE_H__) #define DEFINE_H__ #if defined(_MSC_VER) #define _CRT_SECURE_NO_DEPRECATE 1 #define _CRT_NONSTDC_NO_WARNINGS 1 #pragma warning(disable:4619) /* warning number is not a valid */ #pragma warning(disable:4505) /* unreferenced local function removed */ #if _MSC_VER >= 1910 && _MSC_VER <= 1919 /* VC++2017 */ #pragma warning(disable:4310) /* cast truncates constant value */ #endif #if _MSC_VER == 1800 /* VC++2013 */ #pragma warning(disable:4752) /* found AVX, but AVX support not enabled */ #endif #if _MSC_VER < 1920 /* VC++2019 */ #pragma warning(disable:4307) /* integral constant overflow */ #endif #endif #if defined(__BORLANDC__) #pragma warn -use /* 'identifier' declared but never used */ #pragma warn -inl /* Functions containing * are not expanded inline */ #endif #include <stdlib.h> #define GOLDEN_GAMMA 0x9E3779B97F4A7C15U #define SPLITMIX_A 0xBF58476D1CE4E5B9U #define SPLITMIX_B 0x94D049BB133111EBU #if defined(_M_IX86) || defined(_MSC_VER) typedef __int64 sqword; typedef unsigned __int64 qword; #else typedef long long sqword; typedef unsigned long long qword; #endif #define U01(A) ((sqword)((qword)(A)>>11)*(1.0/9007199254740992.0)) #define ExInitInt(A,B) ExInit(A,(qword)(B)) #if defined(_M_IX86)||defined(__i386)||defined(_M_X64)||defined(__x86_64__) #if defined(_MSC_VER) && _MSC_VER >= 1400 /* VC++2005 */ #include <intrin.h> #define Rdtsc __rdtsc #define Cpuid(V,EAX) __cpuid(V,EAX) #if _MSC_FULL_VER >= 150030729 /* VC++2008 SP1 */ #define CpuidEcx0(V,EAX) __cpuidex(V,EAX,0) #endif #elif defined(_MSC_VER) && _MSC_VER >= 1300 /* VC++2002 */ #include <emmintrin.h> #elif defined(__GNUC__) #include <x86intrin.h> #include <cpuid.h> #define Rdtsc __rdtsc #define Cpuid(V,EAX) __cpuid(EAX,V[0],V[1],V[2],V[3]) #define CpuidEcx0(V,EAX) __cpuid_count(EAX,0,V[0],V[1],V[2],V[3]) #endif #if defined(DQU) #define Load128(A) _mm_loadu_si128(&(A)) #define Load256(A) _mm256_loadu_si256(&(A)) #define Load512(A) _mm512_loadu_si512(&(A)) #define Save128(A,B) _mm_storeu_si128(&(A),B) #define Save256(A,B) _mm256_storeu_si256(&(A),B) #define Save512(A,B) _mm512_storeu_si512(&(A),B) #elif defined(DQA) #define Load128(A) _mm_load_si128(&(A)) #define Load256(A) _mm256_load_si256(&(A)) #define Load512(A) _mm512_load_si512(&(A)) #define Save128(A,B) _mm_store_si128(&(A),B) #define Save256(A,B) _mm256_store_si256(&(A),B) #define Save512(A,B) _mm512_store_si512(&(A),B) #else #define Load128(A) (A) #define Load256(A) (A) #define Load512(A) (A) #define Save128(A,B) ((A)=(B)) #define Save256(A,B) ((A)=(B)) #define Save512(A,B) ((A)=(B)) #endif #endif #if defined(_MSC_VER) && _MSC_VER >= 1400 /* VC++2005 */ #define Bswap _byteswap_uint64 #else static __inline qword Bswap(qword r) { /* GCC might optimize this function using a BSWAP instruction. */ r=(r<<32)|(r>>32); r=((0x0000ffff0000ffffU&r)<<16)|((0xffff0000ffff0000U&r)>>16); r=((0x00ff00ff00ff00ffU&r)<< 8)|((0xff00ff00ff00ff00U&r)>> 8); return r; } #endif #if defined(__GNUC__) #define GCC_NOINLINE __attribute__ ((noinline)) #else #define GCC_NOINLINE #endif static __inline qword SplitMixFinal(qword z) { z=(z^(z>>30))*SPLITMIX_A; z=(z^(z>>27))*SPLITMIX_B; return z^(z>>31); } #if defined(_MSC_VER) #pragma auto_inline(off) #endif GCC_NOINLINE static qword SplitMixFinalN(qword z) { z=(z^(z>>30))*SPLITMIX_A; z=(z^(z>>27))*SPLITMIX_B; return z^(z>>31); } #if defined(_MSC_VER) #pragma auto_inline(on) #endif static __inline qword SplitMixFinalDummy(qword z) { return SplitMixFinalN(z); } #endif
bincode.h
/* bincode.h : coded by Isaku WADA(isaku@pb4.so-net.ne.jp) 2025 */ #if !defined(BINCODE_H__) && \ (defined(_M_IX86)||defined(__i386)||defined(_M_X64)||defined(__x86_64__)) #define BINCODE_H__ #if defined(__GNUC__) #define ALIGN(N) __attribute((aligned(N))) #if defined(DQU) #define DQU__ 1 #else #define DQU__ 0 #endif #elif defined(_MSC_VER) && _MSC_VER >= 1310 /* VC++2003 */ #define ALIGN(N) __declspec(align(N)) #if defined(DQU) #define DQU__ 1 #else #define DQU__ 0 #endif #else #define ALIGN(N) #define DQU__ 1 #endif #if defined(__clang__) #define BINCODE__(X) ALIGN(16) \ const unsigned char X[]__attribute__((section(".text"))) #elif defined(__GNUC__) #define BINCODE__(X) ALIGN(16) \ const unsigned char X[]__attribute__((section(".text#"))) #elif defined(_MSC_VER) #if defined(NO_BTEXT) || _MSC_VER <= 1600 /* VC2010 */ #pragma code_seg(".text") #define BINCODE__(X) __declspec(allocate(".text")) \ ALIGN(16) const unsigned char X[] #else #pragma section(".btext",execute,read) #if defined(MERGE_BTEXT) #pragma comment(linker,"/MERGE:.btext=.text") #endif #define BINCODE__(X) __declspec(allocate(".btext")) \ ALIGN(16) const unsigned char X[] #endif #else #define BINCODE__(X) const unsigned char X[] #endif #if defined(__x86_64__) && defined(__unix__) #define RDI__ 6 /* RCX -> RDI */ #define RSI__ 4 /* RDX -> RSI */ #else #define RDI__ 0 /* RCX */ #define RSI__ 0 /* RDX */ #endif #if ((defined(_M_X64)||defined(__x86_64__))&&!defined(__unix__)) #define CALLER_SAVED 1 #else #define CALLER_SAVED 0 #endif #define INT_BYTE(A) (A)&255,((A)>>8)&255,((A)>>16)&255,((A)>>24)&255 #define RET() 0xc3 /* https://yaya.lsv.jp/multibyte-nop/ */ #define NPAD1() 0x90 #define NPAD2() 0x66,0x90 #define NPAD3() 0x0f,0x1f,0 #define NPAD4() 0x0f,0x1f,0x40,0 #define NPAD5() 0x0f,0x1f,0x44,0,0 #define NPAD6() 0x66,0x0f,0x1f,0x44,0,0 #define NPAD7() 0x0f,0x1f,0x80,0,0,0,0 #define NPAD8() 0x0f,0x1f,0x84,0,0,0,0,0 #define NPAD9() 0x66,0x0f,0x1f,0x84,0,0,0,0,0 #define NPAD10() 0x66,0x66,0x0f,0x1f,0x84,0,0,0,0,0 #define NPAD11() 0x66,0x66,0x66,0x0f,0x1f,0x84,0,0,0,0,0 #define NPAD12() 0x0f,0x1f,0x40,0,0x0f,0x1f,0x84,0,0,0,0,0 #define NPAD13() 0x0f,0x1f,0x40,0,0x66,0x0f,0x1f,0x84,0,0,0,0,0 #define NPAD14() 0x0f,0x1f,0x40,0,0x66,0x66,0x0f,0x1f,0x84,0,0,0,0,0 #define NPAD15() 0x0f,0x1f,0x40,0,0x66,0x66,0x66,0x0f,0x1f,0x84,0,0,0,0,0 #if defined(SWAP_REG) enum { xmm5,xmm4,xmm3,xmm2,xmm1,xmm0,xmm7,xmm6 }; enum { ymm5,ymm4,ymm3,ymm2,ymm1,ymm0,ymm7,ymm6 }; enum { zmm5,zmm4,zmm3,zmm2,zmm1,zmm0,zmm7,zmm6 }; #else enum { xmm0,xmm1,xmm2,xmm3,xmm4,xmm5,xmm6,xmm7 }; enum { ymm0,ymm1,ymm2,ymm3,ymm4,ymm5,ymm6,ymm7 }; enum { zmm0,zmm1,zmm2,zmm3,zmm4,zmm5,zmm6,zmm7 }; #endif #if defined(_M_IX86) || defined(__i386) #define FIRST_ARG() 0x8b,0x4c,0x24,4, /* mov ecx,dword ptr[esp+4] */ #define SECOND_ARG() 0x8b,0x54,0x24,8, /* mov edx,dword ptr[esp+8] */ #else #define FIRST_ARG() #define SECOND_ARG() #endif #define VZEROUPPER() 0xc5,0xf8,0x77 #define MOVDQA(D,S) 0x66,0x0f,0x6f,0xc0+(D)*8+(S) #define XMM_PTR(R,A) 0x81+RDI__+(R)*8,INT_BYTE((A)*16) #define XMM_PTReax(R,A) 0x84+(R)*8,0x01+RDI__,INT_BYTE((A)*16) #define LOAD_XMM(R,A) 0x66+DQU__*141,0x0f,0x6f,XMM_PTR(R,A) #define LOAD_XMMeax(R,A) 0x66+DQU__*141,0x0f,0x6f,XMM_PTReax(R,A) #define LOAD_XMM64esp(R,A) 0xf3,0x0f,0x7e,0x44+(R)*8,0x24,(A) #define LOAD_XMMZ(R) 0x66+DQU__*141,0x0f,0x6f,0x01+(R)*8 #define LOAD_XMM64(R,A) 0xf3,0x0f,0x7e,XMM_PTR(R,A) #define SAVE_XMM(A,R) 0x66+DQU__*141,0x0f,0x7f,XMM_PTR(R,A) #define SAVE_XMMeax(A,R) 0x66+DQU__*141,0x0f,0x7f,XMM_PTReax(R,A) #define SAVE_XMMZ(R) 0x66+DQU__*141,0x0f,0x7f,0x01+RDI__+(R)*8 #define SAVE_XMM64esp(A,R) 0x66,0x0f,0xd6,0x44+(R)*8,0x24,(A) #define SAVE_XMM64(A,R) 0x66,0x0f,0xd6,XMM_PTR(R,A) #define PSLLQ(R,S) 0x66,0x0f,0x73,0xf0+(R),(S) #define PSRLQ(R,S) 0x66,0x0f,0x73,0xd0+(R),(S) #define PXOR(R0,R1) 0x66,0x0f,0xef,0xc0+(R0)*8+(R1) #define PXORp(R,A) 0x66,0x0f,0xef,0x81+(R)*8,INT_BYTE((A)*16) #define POR(R0,R1) 0x66,0x0f,0xeb,0xc0+(R0)*8+(R1) #define PADDQ(R0,R1) 0x66,0x0f,0xd4,0xc0+(R0)*8+(R1) #define PADDQp(R,A) 0x66,0x0f,0xd4,0x81+(R)*8,INT_BYTE((A)*16) #define PSUBQ(R0,R1) 0x66,0x0f,0xfb,0xc0+(R0)*8+(R1) #define PMULUDQ(R0,R1) 0x66,0x0f,0xf4,0xc0+(R0)*8+(R1) #define PSHUFB(R0,R1) 0x66,0x0f,0x38,0x00,0xc0+(R0)*8+(R1) #define PUNPCKLQDQ(R0,R1) 0x66,0x0f,0x6c,0xc0+(R0)*8+(R1) #define PMULLQ(D,R) \ MOVDQA(xmm4,R), MOVDQA(xmm5,D), /* 4+4= 8 8 */ \ PSRLQ(xmm4,32), PSRLQ(xmm5,32), /* 5+5=10 18 */ \ PMULUDQ(xmm4,D), PMULUDQ(xmm5,R),/* 4+4= 8 26 */ \ PSLLQ(xmm4,32), PSLLQ(xmm5,32), /* 5+5=10 36 */ \ PADDQ(xmm4,xmm5), PMULUDQ(D,R), /* 4+4= 8 44 */ \ PADDQ(D,xmm4) /* 4 48 */ #define VMOVDQA(D,S) 0xc5,0xfd,0x6f,0xc0+(D)*8+(S) #define YMM_PTR(R,A) 0x81+RDI__+(R)*8,INT_BYTE((A)*32) #define YMM_PTReax(R,A) 0x84+(R)*8,0x01+RDI__,INT_BYTE((A)*32) #define LOAD_YMM(R,A) 0xc5,0xfd+DQU__,0x6f,YMM_PTR(R,A) #define LOAD_YMMeax(R,A) 0xc5,0xfd+DQU__,0x6f,YMM_PTReax(R,A) #define SAVE_YMM(A,R) 0xc5,0xfd+DQU__,0x7f,YMM_PTR(R,A) #define SAVE_YMMeax(A,R) 0xc5,0xfd+DQU__,0x7f,YMM_PTReax(R,A) #define SAVE_YMMZ(R) 0xc5,0xfd+DQU__,0x7f,0x01+RDI__+(R)*8 #define VPSLLQ(R0,R1,S) 0xc5,0xfd-(R0)*8,0x73,0xf0+(R1),(S) #define VPSRLQ(R0,R1,S) 0xc5,0xfd-(R0)*8,0x73,0xd0+(R1),(S) #define VPXOR(R0,R1,R2) 0xc5,0xfd-(R1)*8,0xef,0xc0+(R0)*8+(R2) #define VPADDQ(R0,R1,R2) 0xc5,0xfd-(R1)*8,0xd4,0xc0+(R0)*8+(R2) #define VPADDQp(R0,R1,A) 0xc5,0xfd-(R1)*8,0xd4,YMM_PTR(R0,A) #define VPSUBQ(R0,R1,R2) 0xc5,0xfd-(R1)*8,0xfb,0xc0+(R0)*8+(R2) #define VPMULUDQ(R0,R1,R2) 0xc5,0xfd-(R1)*8,0xf4,0xc0+(R0)*8+(R2) #define VPSHUFB(R0,R1,R2) 0xc4,0xe2,0x7d-(R1)*8,0,0xc0+(R0)*8+(R2) #define VPBROADCASTQ(R,A) 0xc4,0xe2,0x7d,0x59,YMM_PTR(R,A) #define VPMULLQ(D,R) \ VPSRLQ(ymm4,R,32), VPSRLQ(ymm5,D,32), /* 5+5=10 10 */ \ VPMULUDQ(ymm4,ymm4,D), VPMULUDQ(ymm5,ymm5,R),/* 4+4= 8 18 */ \ VPSLLQ(ymm4,ymm4,32), VPSLLQ(ymm5,ymm5,32), /* 5+5=10 28 */ \ VPADDQ(ymm4,ymm4,ymm5),VPMULUDQ(D,D,R), /* 4+4= 8 36 */ \ VPADDQ(D,D,ymm4) /* 4 40 */ #define VMOVDQA64(D,S) 0x62,0xf1,0xfd,0x48,0x6f,0xc0+(D)*8+(S) #define ZMM_PTR(R,A) 0x81+RDI__+(R)*8,INT_BYTE((A)*64) #define ZMM_PTReax(R,A) 0x84+(R)*8,0x01+RDI__,INT_BYTE((A)*64) #define LOAD_ZMM(R,A) 0x62,0xf1,0xfd+DQU__,0x48,0x6f,ZMM_PTR(R,A) #define LOAD_ZMMeax(R,A) 0x62,0xf1,0xfd+DQU__,0x48,0x6f,ZMM_PTReax(R,A) #define SAVE_ZMM(A,R) 0x62,0xf1,0xfd+DQU__,0x48,0x7f,ZMM_PTR(R,A) #define SAVE_ZMMeax(A,R) 0x62,0xf1,0xfd+DQU__,0x48,0x7f,ZMM_PTReax(R,A) #define SAVE_ZMMZ(R) 0x62,0xf1,0xfd+DQU__,0x48,0x7f,0x01+RDI__+(R)*8 #define VPSLLQz(R0,R1,S) 0x62,0xf1,0xfd-(R0)*8,0x48,0x73,0xf0+(R1),(S) #define VPSRLQz(R0,R1,S) 0x62,0xf1,0xfd-(R0)*8,0x48,0x73,0xd0+(R1),(S) #define VPXORD(R0,R1,R2) 0x62,0xf1,0x7d-(R1)*8,0x48,0xef,0xc0+(R0)*8+(R2) #define VPADDQz(R0,R1,R2) 0x62,0xf1,0xfd-(R1)*8,0x48,0xd4,0xc0+(R0)*8+(R2) #define VPSUBQz(R0,R1,R2) 0x62,0xf1,0xfd-(R1)*8,0x48,0xfb,0xc0+(R0)*8+(R2) #define VPMULLQz(R0,R1,R2) 0x62,0xf2,0xfd-(R1)*8,0x48,0x40,0xc0+(R0)*8+(R2) #define VPADDQzp(R0,R1,A) 0x62,0xf1,0xfd-(R1)*8,0x48,0xd4,ZMM_PTR(R0,A) #define VPSHUFBz(R0,R1,R2) 0x62,0xf2,0x7d-(R1)*8,0x48,0,0xc0+(R0)*8+(R2) #define VPBROADCASTQz(R,A) 0x62,0xf2,0xfd,0x48,0x59,ZMM_PTR(R,A) typedef union m128i { ALIGN(16) qword q[2]; } m128i; typedef union m256i { ALIGN(32) qword q[4]; } m256i; typedef union m512i { ALIGN(64) qword q[8]; } m512i; #if !defined(Cpuid) || !defined(CpuidEcx0) static BINCODE__(Cpuid__)={ FIRST_ARG() SECOND_ARG() #if !defined(__unix__) || defined(__i386) 0x57, /* push edi */ 0x56, /* push esi */ 0x52, /* push edx */ 0x5e, /* pop esi */ 0x51, /* push ecx */ 0x5f, /* pop edi */ #endif 0x53, /* push ebx */ 0x31,0xc9, /* xor ecx,ecx */ 0x89,0xf0, /* mov eax,esi */ 0x0f,0xa2, /* cpuid */ 0x89,0x07, /* mov dword ptr[edi],eax */ 0x89,0x5f,0x04, /* mov dword ptr[edi+4],ebx */ 0x89,0x4f,0x08, /* mov dword ptr[edi+8],ecx */ 0x89,0x57,0x0c, /* mov dword ptr[edi+12],edx */ 0x5b, /* pop ebx */ #if !defined(__unix__) || defined(__i386) 0x5e, /* pop esi */ 0x5f, /* pop edi */ #endif 0xc3}; /* ret */ #if !defined(Cpuid) #define Cpuid ((void(*)(int*V,int EAX))(size_t)&Cpuid__) #endif #if !defined(CpuidEcx0) #define CpuidEcx0 ((void(*)(int*V,int EAX))(size_t)&Cpuid__) #endif #endif #if !defined(Bswap) && ( defined(_MSC_VER) || defined(__BORLANDC__) ) static BINCODE__(Bswap__)={ 0x8b,0x54,0x24,4, /* mov edx,dword ptr[esp+4] */ 0x8b,0x44,0x24,8, /* mov eax,dword ptr[esp+8] */ 0x0f,0xca, /* bswap edx */ 0x0f,0xc8, /* bswap eax */ 0xc3, /* ret */ }; #define Bswap ((qword(*)(qword))(size_t)Bswap__) #endif #if !defined(Rdtsc) #if defined(_M_X64) || defined(__x86_64__) static BINCODE__(Rdtsc__)={ 0x0f,0x31, /* rdtsc */ 0x48,0xc1,0xe2,0x20, /* shl rdx,32 */ 0x48,0x09,0xd0, /* or rax,rdx */ 0xc3}; /* ret */ #elif defined(_M_IX86) || defined(__i386) static BINCODE__(Rdtsc__)={0x0f,0x31,0xc3}; #endif #define Rdtsc ((qword(*)(void))(size_t)Rdtsc__) #endif #endif
xgf2.h
/* xgf2.h : coded by Isaku WADA(isaku@pb4.so-net.ne.jp) 2025 An Implementation of Xorshift RNG */ #if !defined(XGF2_H__) #define XGF2_H__ #if !defined(XID) #if defined(XIDA) #define XID 'A' #elif defined(XIDB) #define XID 'B' #elif defined(XIDC) #define XID 'C' #elif defined(XIDD) #define XID 'D' #elif defined(XIDE) #define XID 'E' #elif defined(XIDF) #define XID 'F' #elif defined(XIDG) #define XID 'G' #elif defined(XIDH) #define XID 'H' #elif defined(XIDI) #define XID 'I' #elif defined(XIDJ) #define XID 'J' #elif defined(XIDK) #define XID 'K' #elif defined(XIDL) #define XID 'L' #elif defined(XIDM) #define XID 'M' #elif defined(XIDN) #define XID 'N' #elif defined(XIDO) #define XID 'O' #elif defined(XIDP) #define XID 'P' #elif defined(XIDQ) #define XID 'Q' #elif defined(XIDR) #define XID 'R' #elif defined(XIDS) #define XID 'S' #elif defined(XIDT) #define XID 'T' #elif defined(XIDU) #define XID 'U' #elif defined(XIDV) #define XID 'V' #elif defined(XNSD3) #define XID 'A' #elif defined(XNSD4) #define XID 'C' #elif defined(XNSD5) #define XID 'K' #elif defined(XNSD6) #define XID 'M' #elif defined(XNSD7) #define XID 'N' #elif defined(XNSD8) #define XID 'S' #elif defined(XNSD9) #define XID 'T' #elif defined(XNSD10) #define XID 'U' #else #define XID 'U' #endif #endif #if XID == 'A' #define XNABC 3423902 /*'A'35%*/ #elif XID == 'B' #define XNABC 4353821 /*'B'29%*/ #elif XID == 'C' #define XNABC 4364111 /*'C'33%*/ #elif XID == 'D' #define XNABC 4374913 /*'D'35%*/ #elif XID == 'E' #define XNABC 4414511 /*'E'36%*/ #elif XID == 'F' #define XNABC 4514601 /*'F'22%*/ #elif XID == 'G' #define XNABC 4514613 /*'G'21%*/ #elif XID == 'H' #define XNABC 4596001 /*'H'5%*/ #elif XID == 'I' #define XNABC 4614203 /*'I'14%*/ #elif XID == 'J' #define XNABC 5363703 /*'J'38%*/ #elif XID == 'K' #define XNABC 5474317 /*'K'43%*/ #elif XID == 'L' #define XNABC 6374309 /*'L'34%*/ #elif XID == 'M' #define XNABC 6595505 /*'M'36%*/ #elif XID == 'N' #define XNABC 7383912 /*'N'23%*/ #elif XID == 'O' #define XNABC 7464507 /*'O'34%*/ #elif XID == 'P' #define XNABC 7475508 /*'P'5%*/ #elif XID == 'Q' #define XNABC 8353914 /*'Q'24%*/ #elif XID == 'R' #define XNABC 8364913 /*'R'8%*/ #elif XID == 'S' #define XNABC 8513613 /*'S'18%*/ #elif XID == 'T' #define XNABC 9343305 /*'T'20%*/ #elif XID == 'U' #define XNABC 10353720 /*'U'25%*/ #elif XID == 'V' #define XNABC 10493701 /*'V'19%*/ #else #error Bad XID #endif #define XNSD (XNABC/1000000) #define ENUM_ABCM(M) \ enum { a=XNABC/10000%100,b=XNABC/100%100,c=XNABC%100,m=M } #define JEXP (XNSD*64-8) #if defined(_MSC_VER) && _MSC_VER >= 1400 /* VC++2005 */ #define Bsr32(ptr,mask) \ ((int)_BitScanReverse((unsigned long*)(ptr),(unsigned long)(mask))) #elif defined(__GNUC__) && __GNUC__ >= 4 #include <x86intrin.h> static __inline int Bsr32(int*ptr,unsigned mask) { if (mask==0) return 0; *ptr=__bsrd((int)mask); return 1; } #else static __inline int Bsr32(int*ptr,unsigned mask) { union tagBsr32 { qword q; double d; } x; if (mask==0) return 0; x.d=(double)mask+0.5; x.q>>=52; x.q-=1023; *ptr=(int)x.q; return 1; } #endif #if !defined(XEXP) && defined(XEXP0) #define XEXP 0 #elif !defined(XEXP) && defined(XEXP1) #define XEXP 1 #elif !defined(XEXP) && defined(XEXP2) #define XEXP 2 #elif !defined(XEXP) && defined(XEXP3) #define XEXP 3 #elif !defined(XEXP) && defined(XEXP4) #define XEXP 4 #elif !defined(XEXP) && defined(XEXP5) #define XEXP 5 #elif !defined(XEXP) && defined(XEXP6) #define XEXP 6 #elif !defined(XEXP) && defined(XEXP7) #define XEXP 7 #elif !defined(XEXP) && defined(XEXP8) #define XEXP 8 #else #define XEXP 3 #endif #if defined(FAST3) #define IF_NOT_FAST3(A) #define IF_NOT_FAST2(A) #define IF_FAST1(A) #define IF_FAST0(A) #define XFAST 3 #elif defined(FAST2) #define IF_NOT_FAST3(A) A #define IF_NOT_FAST2(A) #define IF_FAST1(A) #define IF_FAST0(A) #define XFAST 2 #elif defined(FAST1) #define IF_NOT_FAST3(A) A #define IF_NOT_FAST2(A) #define IF_FAST1(A) A #define IF_FAST0(A) #define XFAST 1 #elif defined(FAST0) #define IF_NOT_FAST3(A) A #define IF_NOT_FAST2(A) A #define IF_FAST1(A) A #define IF_FAST0(A) A #define XFAST 0 #else #define IF_NOT_FAST3(A) A #define IF_NOT_FAST2(A) A #define IF_FAST1(A) A #define IF_FAST0(A) #define XFAST 9 #endif #define RND_EXP (XNSD*64) #define RND_NUM ((RND_EXP+32)/32) #define RND_NUM2 (((RND_EXP+1)*2+32-1)/32) typedef unsigned xgf2_t[RND_NUM2]; extern GCC_NOINLINE void InitJumpFunc(qword*p); extern GCC_NOINLINE void JumpHelp (qword*p,qword*o,const xgf2_t e,int mul,int back); extern GCC_NOINLINE void QxJumpHelp (double lag,qword*cnt,qword*xs,int num,int*idx); extern const xgf2_t NextChara; extern const xgf2_t BackChara; # if defined(FAST3) extern GCC_NOINLINE void XJumpHelp(qword*x,qword*o, const unsigned*m,int num,double lag); #else extern GCC_NOINLINE void XJumpHelp(qword*x,qword*aux,qword*o, const unsigned*m,int num,int qnum,double lag); #endif extern GCC_NOINLINE int InitByU64(qword*p,qword s); extern GCC_NOINLINE int InitByStr(unsigned*p,const char*key); extern void xgf2square(xgf2_t c,const xgf2_t a); extern void xgf2mult(xgf2_t c,const xgf2_t a,const xgf2_t b); extern void xgf2shiftleft(xgf2_t c,const xgf2_t a,const int s); extern int xgf2deg(const xgf2_t a); extern void xgf2rem(xgf2_t r,const xgf2_t a,const xgf2_t b); extern qword XorNextState1(qword*p); extern qword XorBackState1(qword*p); #endif
xgf2.c
/* xgf2.c : coded by Isaku WADA(isaku@pb4.so-net.ne.jp) 2025 An Implementation of Xorshift RNG */ #if !defined(XGF2_C__) #define XGF2_C__ #include "define.h" #include "xgf2.h" #include <math.h> #include <string.h> /* Efficient Jump Ahead for F2-Linear Random Number Generators http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/ARTICLES/jumpf2-printed.pdf */ #if XID == 'A' const xgf2_t NextChara=/*3'A'35%*/{0x1351001,0xf28721f2U,0x2e4b10f9, 0x4f90dd0e,0x12bc224,0x35,1}; const xgf2_t BackChara=/*3'A'35%*/{0x58000001,0x4887a901,0xe17613e4U, 0x3e11a4e8,0x9f09c29fU,0x115900,1}; static const unsigned InitJumpTbl[]/*3 184'A'35%*/={0x9902f954U,0xd374746dU, 0xbbe9b323U,0xb832d856U,0xa0681055U,0x3a629b4e}; #elif XID == 'B' const xgf2_t NextChara=/*4'B'29%*/{0x40000001,0xf9bc1415U,0x2504e07, 0xf240a755U,0xbb127a91U,0x5760b094,0x140000,0,1}; const xgf2_t BackChara=/*4'B'29%*/{1,0x5000,0x521a0dd4,0x12bc91ba, 0x55ca049f,0xc0e41481U,0x50507b3f,5,1}; static const unsigned InitJumpTbl[]/*4 248'B'29%*/={0x312277f6,0x9b047320U, 0xb6b33139U,0x29fa647c,0xc8ce1c1dU,0xd0357385U,0x98a43f67U,0x43f4d784}; #elif XID == 'C' const xgf2_t NextChara=/*4'C'33%*/{0xc0cc0001U,0xafc05090U,0x3c0fafcf, 0x33c0c3,0x50cc9390,0x9030fc63U,0xc00c900dU,0xc0,1}; const xgf2_t BackChara=/*4'C'33%*/{0x6000001,0x60126006,0x8c7e1813U, 0x13926615,0x86079800U,0xe7ebe079U,0x121407eb,0x6606,1}; static const unsigned InitJumpTbl[]/*4 248'C'33%*/={0xdb480769U,0x9bf3adefU, 0xbe13f37bU,0x8ce99b5fU,0x6d01368b,0xb4708c6fU,0x21f96047,0xc3ebfa22U}; #elif XID == 'D' const xgf2_t NextChara=/*4'D'35%*/{0xe84c8a01U,0xa080224aU,0xc0a44a6cU, 0x16d366c2,0x48ccb055,0x60246826,0xc0acc8afU,0x882c2,1}; const xgf2_t BackChara=/*4'D'35%*/{0x86822001U,0xea266a06U,0xc82c480dU, 0x541a6624,0x86cd96d1U,0x6ca44a06,0xa488020aU,0xa2642e,1}; static const unsigned InitJumpTbl[]/*4 248'D'35%*/={0x6de18507,0x4e03329f, 0xfd299e78U,0xfaa767ceU,0x4478ba6c,0xe49cb80,0xacfce590U,0x5c55ee3}; #elif XID == 'E' const xgf2_t NextChara=/*4'E'36%*/{0x40028001,0x82a04,0x8f42622aU, 0x376d57fb,0x4cacdfa6,0xc284486aU,0xc28ee0adU,0x28286,1}; const xgf2_t BackChara=/*4'E'36%*/{0xc2828001U,0x6a0ee286,0xac244287U, 0xcbf66a64U,0xbfd56dd8U,0xa88c85e3U,0x40a82000,0x28004,1}; static const unsigned InitJumpTbl[]/*4 248'E'36%*/={0xd2c3d8d5U,0xb1a8cd84U, 0x57064969,0x89dff015U,0x53451394,0xec2923,0xd4ac4f82U,0x881dc8aeU}; #elif XID == 'F' const xgf2_t NextChara=/*4'F'22%*/{0xa08001,0x22880080,0xa8aa88aaU,0xa00a, 0x20288080,0xa8a08a0aU,0xa80aaa83U,0xa0200a,1}; const xgf2_t BackChara=/*4'F'22%*/{0xa0080a01U,0x82aaa02aU,0xa0a20a2bU, 0x2022808,0xa00a0000U,0xaa22aa2aU,0x2002288,0x20a00,1}; static const unsigned InitJumpTbl[]/*4 248'F'22%*/={0x416b1098,0x3414918, 0x5526b3aa,0x323e1ec9,0xdb1ed708U,0xf476aac7U,0xb2eea960U,0x1d9c49f1}; #elif XID == 'G' const xgf2_t NextChara=/*4'G'21%*/{1,0x100800a2,0x3551f8,0x40000000, 0x434801b3,0x103d504a,0x3551f9,0,1}; const xgf2_t BackChara=/*4'G'21%*/{1,0x3f155800,0xa4157811U,0x9b002584U,5, 0x3f155800,0x8a002010U,0,1}; static const unsigned InitJumpTbl[]/*4 248'G'21%*/={0xdd78de56U,0x4ecb064f, 0x334f50fe,0xf18ea6d1U,0xb61b9893U,0x82b63304U,0xa760de35U,0x862d953cU}; #elif XID == 'H' const xgf2_t NextChara=/*4'H' 5%*/{1,0,0,0,0,0x213b4800,3,0,1}; const xgf2_t BackChara=/*4'H' 5%*/{1,0x80000000U,0x25b909,0,0,0,0,0,1}; static const unsigned InitJumpTbl[]/*4 248'H' 5%*/={0xd59019c,0x219728e0, 0xf6ec530eU,0xd848c63dU,0x3dd2b418,0xe3911798U,0xa3633d78U,0x3844d353}; #elif XID == 'I' const xgf2_t NextChara=/*4'I'14%*/{1,0x80000000U,0x82800002U,0x88a02, 0xa820a28aU,0x8088200,0x828882a9U,2,1}; const xgf2_t BackChara=/*4'I'14%*/{0x80000001U,0x2a822282,0x822021, 0xa28a082aU,0x80a22000U,0x80000282U,2,0,1}; static const unsigned InitJumpTbl[]/*4 248'I'14%*/={0x6eae969c,0x15938e5c, 0x6b54c395,0xff227c65U,0xc06be65bU,0xa16503c,0x5f37cfbb,0xa886790bU}; #elif XID == 'J' const xgf2_t NextChara=/*5'J'38%*/{1,0xf000000,0xc53ef033U,0xaedcc6c1U, 0xa7e43404U,0xf0ca9d25U,0xfc3af233U,0x59ee3502,0x55edc93f,3,1}; const xgf2_t BackChara=/*5'J'38%*/{0x80000001U,0xf9276f55U,0x8158ef35U, 0x989eb87eU,0x4972a61f,0x40584fcb,0x6c676ea,0x981ef947U,0x1e1,0,1}; static const unsigned InitJumpTbl[]/*5 312'J'38%*/={0x15b41534,0x5e8dedb5, 0xefd7fe6dU,0x75704abc,0x933c755,0x463432b8,0x39e734ce,0xfbe819d1U, 0x268a2a19,0xefe9c2a8U}; #elif XID == 'K' const xgf2_t NextChara=/*5'K'43%*/{0xf030001,0xc570005,0xf1f8c32dU, 0x39a62e4d,0x33055bf1,0x16fbcf96,0x960ea175U,0x2621efa,0xc5c124f,0xf0c59, 1}; const xgf2_t BackChara=/*5'K'43%*/{0x3461e001,0xe4907461U,0xbef08c81U, 0x5d0ae0d2,0xd3e7bed1U,0x1fb54198,0x64e8cb39,0x69863f1f,0x4001d461,0x181e1, 1}; static const unsigned InitJumpTbl[]/*5 312'K'43%*/={0x71afc427,0xc9d5bd73U, 0xb591b4f9U,0x762cc92c,0x666c6b4b,0x10cc9b1f,0xa365fca5U,0x909e93c3U, 0x4bb44d35,0x3282a64b}; #elif XID == 'L' const xgf2_t NextChara=/*6'L'34%*/{1,0xf0800000U,0x8a5f880aU,0xd782fa70U, 0x23637591,0x3eeed510,0xa3972270U,0x5970d8,0x207edd5b,0xd782b00dU,0xa54,0, 1}; const xgf2_t BackChara=/*6'L'34%*/{1,0x54a00000,0x601a83d6,0xb576fc09U, 0x361d3401,0x1c89d38a,0x1156eef8,0x135d8d88,0x1cbe83d7,0xa023f4a2U,0x21e,0, 1}; static const unsigned InitJumpTbl[]/*6 376'L'34%*/={0xeaff341dU,0x7bc73739, 0xa2f70af,0x1a7188dd,0x38b68644,0xf939221dU,0x946548c0U,0x8085a8bcU, 0xbf02e096U,0xb286195eU,0xcd940bd9U,0xdb1a6a8bU}; #elif XID == 'M' const xgf2_t NextChara=/*6'M'36%*/{0x6000001,0x7f866078,0x80000800U, 0x6e700608,0x7f8e977e,0x6000000,0x79866078,0x7fe66878,0x8e700e0fU, 0x11f91616,0x7f8e977f,0,1}; const xgf2_t BackChara=/*6'M'36%*/{1,0xfdd2e3fcU,0xd0d13f11U,0xe0e01ce2U, 0x3c2ccffd,0x3c0cc33c,0xc0,0xfdd2e3fcU,0x20c01cec,0x200002,0x3c0cc3fc,0xc0, 1}; static const unsigned InitJumpTbl[]/*6 376'M'36%*/={0x10b5bc00,0xa1903bd3U, 0x8a626856U,0x5e333c4f,0x935af335U,0xe14a8a4aU,0xb7ef7241U,0xdbfe7339U, 0x8ea9fbd6U,0xf2fa799bU,0x9d264e92U,0x85fcfdc5U}; #elif XID == 'N' const xgf2_t NextChara=/*7'N'23%*/{1,0x10050000,1,0xff000000U,0x6000cf, 0xde24ea14U,0xff003130U,0xa1ffcf,0x3430ea75,0x3031ee15,0x10050000, 0x10041004,0x11005,0,1}; const xgf2_t BackChara=/*7'N'23%*/{1,0x40110000,0x40104011,0x14010, 0x50ef1818,0x5cae1859,0xe7ff0a01U,0x191801ff,0x50ae48f6,0xe6000c00U,0x1ff, 0,0x14011,0,1}; static const unsigned InitJumpTbl[]/*7 440'N'23%*/={0x8b4eb637U,0x12b3132d, 0x4d9134c2,0x84ab614bU,0x73adaa74,0xc5a46640U,0x647c936a,0xbfd07fffU, 0xddc67bf4U,0xde64db76U,0x61092fba,0xf9fd5c7dU,0x74d73870,0xd763f4b9U}; #elif XID == 'O' const xgf2_t NextChara=/*7'O'34%*/{1,0x30000,0x11103300,0x55014554, 0xfafa2320U,0x3c00f300,0xfd333011U,0x8dde33ceU,0xda400dffU,0x8872e954U, 0x50fdfcd8,0x36062313,0x33003302,0,1}; const xgf2_t BackChara=/*7'O'34%*/{1,0x81980198U,0x9188c0d8U,0x367f7e15, 0x552e9c22,0xff6004b6U,0xe798f763U,0x1019997e,0x19e0079,0x988bebe, 0x55450154,0x1981110,0x18000,0,1}; static const unsigned InitJumpTbl[]/*7 440'O'34%*/={0xb763b4edU,0x8921f35, 0x46701258,0x553bcce9,0x9fc300c3U,0xf86ae88dU,0xb645a05eU,0x20e75cbc, 0x3f3d967,0xda92311aU,0xf8c28e38U,0xe726bb9,0x580f9628,0xb9b1fa32U}; #elif XID == 'P' const xgf2_t NextChara=/*7'P' 5%*/{1,0,0,0,0,0,0,0,0xfff00000U,0xc3cc,0,0,1, 0,1}; const xgf2_t BackChara=/*7'P' 5%*/{1,0,1,0,0x67860000,0x1ffe,0,0,0,0,0,0,0, 0,1}; static const unsigned InitJumpTbl[]/*7 440'P' 5%*/={0x2b655f01,0x6dff44c8, 0x6a3808a0,0xe53bcc32U,0xa43ec3a4U,0xeff6ec31U,0x4b3629a4,0xec0a92a8U, 0x349c6112,0xa5439f4eU,0x3677db17,0x478ee9a7,0x98aefe21U,0xa258bedcU}; #elif XID == 'Q' const xgf2_t NextChara=/*8'Q'24%*/{1,0,0,0xa00,0x54780ac0,0x6cbe8167, 0x6a307860,0x5c90bfb3,0xd91f5848U,0x181ee119,0x999e1418U,0x7fff39,0xaa0,0, 1,0,1}; const xgf2_t BackChara=/*8'Q'24%*/{1,0,1,0xaa00000,0x39fffc00,0x3050f333, 0x310ef030,0x2435f137,0x9bfa1274U,0xc3c18ad,0xcd02fa6cU,0x6a03c55,0xa00000, 0,0,0,1}; static const unsigned InitJumpTbl[]/*8 504'Q'24%*/={0x4defe857,0xb99097f7U, 0xf8a6b963U,0xee8705d6U,0x38900a39,0xdc733801U,0x53cc6e1f,0x79998d05, 0x6919e3b6,0x4d29596e,0xc0278ec7U,0x7e309b99,0x1dc82ebb,0x72b8ddc9, 0x14b99137,0x28547a56}; #elif XID == 'R' const xgf2_t NextChara=/*8'R' 8%*/{1,0,0,0,0,0,0,0x3c3c000,0xfcfff30fU,0xff, 0,0,0,0,1,0,1}; const xgf2_t BackChara=/*8'R' 8%*/{1,0,1,0,0,0,0xfe000000U,0xe19ffe7fU, 0x78781,0,0,0,0,0,0,0,1}; static const unsigned InitJumpTbl[]/*8 504'R' 8%*/={0xe4246f77U,0x6872bec7, 0x8a59eb8cU,0x46e109b6,0x22763c0e,0xb6d29598U,0x66340415,0xd3b9224eU, 0x3be806ef,0x941b7a5eU,0xb4de4812U,0x1eb2b5ee,0xa42c493fU,0x57d9e371, 0xbfe18282U,0x1c6149ce}; #elif XID == 'S' const xgf2_t NextChara=/*8'S'18%*/{1,0x80000000U,0x88000000U,0x2910014, 0x14000220,0xa0119114U,0xaa008800U,0x30200000,0xa2220000U,0x28002000, 0xa8912214U,0x26b13234,0xb633b116U,0x8111914,0x88008801U,0,1}; const xgf2_t BackChara=/*8'S'18%*/{1,0x220022,0x51311021,0xd11b98daU, 0x58991ac8,0x5089122a,0x80028,0x888a,0x818,0x2200aa,0x5113100a,0x8800050, 0x50011280,0x22,2,0,1}; static const unsigned InitJumpTbl[]/*8 504'S'18%*/={0xfb15fce0U,0x61688ea1, 0x54f24a39,0x962a693fU,0x742da0d3,0x5543ab8a,0x7adec671,0xf70bf101U, 0xb4c01cd3U,0x1fc47b27,0xc1708f8cU,0xc5017c20U,0xff10509aU,0x66d43210, 0xfb2d128aU,0x6d72ca9}; #elif XID == 'T' const xgf2_t NextChara=/*9'T'20%*/{1,0,0x1010000,0x505,0x10000,0xf5e0014, 0x1103223,0x1fe050b,0x5115242,0xf4b0114,0x33332c69,0x1fe05ea,0x31606e9, 0x1000105,0x11141115,0,1,0,1}; const xgf2_t BackChara=/*9'T'20%*/{1,0,1,0x51105110,0x41000101,0x2ec0d181, 0xaf40ff01U,0x2c699998,0x5101a5e1,0x84951140U,0xa140ff00U,0x88981101U, 0x5000f5e1,0x10000,0x41400000,0x10101,0,0,1}; static const unsigned InitJumpTbl[]/*9 568'T'20%*/={0x316493f9,0xc558ea9dU, 0x256857fb,0x7cfe2316,0xe7630b49U,0xbbf9001dU,0xb575f441U,0xf55f0edaU, 0x9f47b766U,0x2ee74ded,0x34887206,0xbf47e9e,0x81b5968cU,0xe796dbbaU, 0xe01b89b5U,0x40fcdc9f,0x6f2d95e3,0x970ab47aU}; #elif XID == 'U' const xgf2_t NextChara=/*10'U'25%*/{1,0,0,0,0xfe000000U,0xf80019ffU, 0xaab2e001U,0x19ce01e7,0x87e79fb0U,0x6ff767f6,0xe8179167U,0xab2acacaU, 0xe01e1829U,0xe1f981b1U,0x199981,0,0,0,1,0,1}; const xgf2_t BackChara=/*10'U'25%*/{1,0,1,0,0,0x3333000,0x1b033f0f, 0x2830f00f,0xa6a6a9abU,0xcd13d02eU,0xdfcddfedU,0x1bf3cfc2,0xcf00e730U, 0xe9aab,0xff30003fU,0xff,0,0,0,0,1}; static const unsigned InitJumpTbl[]/*10 632'U'25%*/={0x5a26d50d,0xd3120b7cU, 0xd6a5ccc4U,0x8b1e9065U,0x9c75de55U,0x8650bf8eU,0xdb8ee886U,0x39bba333, 0x14f81d9,0x48a54708,0x397c67a0,0x22d252ad,0xf6998c5eU,0x7315960d, 0xa3cfe859U,0xb751efdbU,0xb3e3c260U,0x7377147e,0xba83fa64U,0x670ec0e5}; #elif XID == 'V' const xgf2_t NextChara=/*10'V'19%*/{1,0,0,0x600000,0x780,0,0,0x1400, 0x7f54060,0x780,0x80006000U,0x660007,0x6007f8,0x2610780,0x40002f94, 0x812d6002U,0x2fa20262,0xad77b1b8U,0x7f547cb,0x780,1}; const xgf2_t BackChara=/*10'V'19%*/{0x3c00001,0xa7c55fc0U,0x3b1bdd6b, 0x8c808be8U,0x800d6902U,0x53e80004,0x3c10c80,0x3fc00c00,0xc000cc00U, 0xc0003,0x3c00000,0xc055fc0,0x500000,0,0,0x3c00000,0xc00,0,0,0,1}; static const unsigned InitJumpTbl[]/*10 632'V'19%*/={0x59e93248,0x68620b3e, 0xb1371a1,0xb6b0d3c9U,0xc9ec5e3dU,0x83343a3bU,0x7860fca5,0x9807aa28U, 0x76a469b7,0x2a568945,0x92dadc79U,0x22289e5,0xa31ff7b0U,0xe7d0b5d5U, 0xede63a4dU,0x436edf1e,0x5cd959d6,0xcc8b46ebU,0x9bea3b2bU,0x8947e176U}; #else #error Bad XID #endif static qword XorNextState(qword*p,const int s) { ENUM_ABCM(0); int i; qword x=p[0],t=x; for (i=0;i<(XNSD-1)*s;i+=s) p[i]=p[i+s]; x<<=a; t^=x; x=t; x>>=b; t^=x; x=p[(XNSD-1)*s]; t^=x; x>>=c; t^=x; return p[(XNSD-1)*s]=t; } static qword XorBackState(qword*p,const int s) { ENUM_ABCM(0); int i; qword t=p[(XNSD-1)*s],x=p[(XNSD-2)*s]; for (i=(XNSD-1)*s;i>0;i-=s) p[i]=p[i-s]; t^=x; x>>=c; t^=x; x=t; x>>=b; t^=x; x=t; x<<=a; t^=x; p[0]=t; return p[(XNSD-1)*s]; } static void SetNonZero(qword*p) { static char alt[]="RANDDEFAULTS"; int i; for (i=0;i<XNSD;i++) if (p[i]) break; if (i<XNSD) return; for (i=0;i<XNSD;i++) p[i]=(qword)alt[i]; } GCC_NOINLINE int InitByU64(qword*q,qword seed) { int i; for (i=0;i<XNSD;i++) q[i]=SplitMixFinalN(seed += GOLDEN_GAMMA); return 0; } GCC_NOINLINE int InitByStr(unsigned*p,const char*key) { const unsigned char*ptr=(const unsigned char*)(key?key:""); int i=1,j=0,y,slen=(int)strlen((const char*)ptr); int ilen=(slen+3)/4,k=(XNSD*2>ilen?XNSD*2:ilen); unsigned data,c; if (InitByU64((qword*)p,19650218)) return 1; for (;k;k--) { data=0; for (y=0;y<4;y++) { if ((c=ptr[4*j+y])==0) break; data|=c<<(y*8); } p[i]=(p[i]^((p[i-1]^(p[i-1]>>30))*1664525))+data+j; i++; j++; if (i>=XNSD*2) { p[0]=p[XNSD*2-1]; i=1; } if (j>=ilen) j=0; } for (k=XNSD*2-1;k;k--) { p[i]=(p[i]^((p[i-1]^(p[i-1]>>30))*1566083941))-i; if (++i>=XNSD*2) { p[0]=p[XNSD*2-1]; i=1; } } SetNonZero((qword*)p); return 0; } void xgf2square(xgf2_t c,const xgf2_t a) { unsigned ch,cl,aa; int i,j,i2,j2; for (i=0;i<RND_NUM2;i++) c[i]=0; for (j2=j=0;j<RND_NUM;j++,j2+=2) { aa=a[j]; if (aa==0) continue; ch=cl=0; for (i2=i=0;i<32;i++,i2+=2) { if (aa&(1U<<i)) { if (i2>=32) ch|=1U<<(i2-32); else cl|=1U<<i2; } } c[j2]^=cl; if (ch) c[j2+1]^=ch; } } void xgf2mult(xgf2_t c,const xgf2_t a,const xgf2_t b) { unsigned ch,cl,d,aa,bb; int i,j,k,ij; for (i=0;i<RND_NUM2;i++) c[i]=0; for (j=0;j<RND_NUM;j++) { bb=b[j]; if (bb==0) continue; for (i=0;i<RND_NUM;i++) { aa=a[i]; if (aa==0) continue; ij=i+j; cl=ch=0; for (k=0;k<32;k++) { if (aa&(1U<<k)) ch^=bb; d =(ch>>1)|(cl<<31); cl=(cl>>1)|(ch<<31); ch=d; } c[ij]^=cl; if (ch) c[ij+1]^=ch; } } } void xgf2shiftleft(xgf2_t c,const xgf2_t a,const int s) { int i,ib=s&31; int na1=RND_NUM2-((s-ib)>>5),iw=s>>5; unsigned x; for (i=0;i<RND_NUM2;i++) c[i]=0; if (ib==0) { for (i=0;i<na1;i++) c[i+iw]^=a[i]; return; } for (i=0;i<na1;i++) { x=a[i]; c[i+iw+1]^=x>>(32-ib); c[i+iw]^=x<<ib; } } int xgf2deg(const xgf2_t a) { int i,j=0; for (i=RND_NUM2-1;i>=0;i--) if (Bsr32(&j,a[i])) break; if (i==-1) return -1; return i*32+j; } void xgf2rem(xgf2_t r,const xgf2_t a,const xgf2_t b) { xgf2_t t; int i,dr=xgf2deg(a); for (i=0;i<RND_NUM2;i++) r[i]=a[i]; while (dr>=RND_EXP) { xgf2shiftleft(t,b,dr-RND_EXP); for (i=0;i<RND_NUM2;i++) r[i]^=t[i]; dr=xgf2deg(r); } } GCC_NOINLINE void InitJumpFunc(qword*p) { ENUM_ABCM(0); qword q[XNSD],t,x; int i,j; memset(q,0,sizeof q); for (j=RND_EXP-1;;j--) { if (InitJumpTbl[j>>5]&(1U<<(j&31))) for (i=0;i<XNSD;i++) q[i]^=p[i]; if (j==0) break; for (x=q[0],t=x,i=0;i<XNSD-1;i++) q[i]=q[i+1]; x<<=a; t^=x; x=t; x>>=b; t^=x; x=q[XNSD-1]; t^=x; x>>=c; t^=x; q[XNSD-1]=t; } memcpy(p,q,sizeof q); } GCC_NOINLINE void JumpHelp(qword*p,qword*o, const xgf2_t e,int num,int back) { ENUM_ABCM(0); int j,i,k,n=XNSD*num; memset(o,0,sizeof*o*n); for (j=RND_EXP-1;;j--) { if (e[j>>5]&(1U<<(j&31))) for (i=0;i<n;i++) o[i]^=p[i]; if (j==0) break; for (k=0;k<num;k++) if (back==0) XorNextState(o+k,num); else XorBackState(o+k,num); } for (i=0;i<n;i++) p[i]=o[i]; } #if defined(FAST3) GCC_NOINLINE void XJumpHelp(qword*x,qword*o, const unsigned*m,int num,double lag) { xgf2_t e,d,s; int i; e[0]=1; s[0]=2; for (i=1;i<RND_NUM2;i++) e[i]=s[i]=0; for (;;) { if (modf(lag/2.0,&lag)!=0) { xgf2mult(d,e,s); xgf2rem(e,d,m); } if (lag==0.0) break; xgf2square(d,s); xgf2rem(s,d,m); } JumpHelp(x,o,e,num,m==BackChara); } #else GCC_NOINLINE void XJumpHelp(qword*x,qword*aux,qword*o, const unsigned*m,int num,int qnum,double lag) { xgf2_t e,d,s; int i; qword add=GOLDEN_GAMMA*num; if (m==BackChara) add*=-1; e[0]=1; s[0]=2; for (i=1;i<RND_NUM2;i++) e[i]=s[i]=0; for (;;) { if (modf(lag/2.0,&lag)!=0) { xgf2mult(d,e,s); xgf2rem(e,d,m); for (i=0;i<qnum;i++) aux[i]+=add; } if (lag==0.0) break; xgf2square(d,s); xgf2rem(s,d,m); add+=add; } JumpHelp(x,o,e,num,m==BackChara); } #endif #define GG_MUL_INV 0xF1DE83E19937733DU #if 0 // Calculate and print GG_MUL_INV #include <stdio.h> #if defined(_MSC_VER) || defined(__BORLANDC__) typedef __int64 sqword; typedef unsigned __int64 qword; #else typedef long long sqword; typedef unsigned long long qword; #endif #define GOLDEN_GAMMA 0x9E3779B97F4A7C15U #define GG_MUL_INV 0xF1DE83E19937733DU int main(void) { qword r=1,g=GOLDEN_GAMMA; int i; for (i=0;i<62;i++) { r*=g; g*=g; } printf(" #define GG_MUL_INV 0x%08X%08XU\n", (unsigned)(r>>32),(unsigned)r); return 0; } #endif GCC_NOINLINE void QxJumpHelp(double lag,qword*cnt,qword*xs,int num,int*idx) { const unsigned*m=NextChara; int i; qword a=GOLDEN_GAMMA,b; int mod=(unsigned short)(idx?*idx:*cnt*GG_MUL_INV); if (lag<0) { lag=-lag; m=BackChara; a=0-a; } i=(int)(modf(lag/65536,&lag)*65536+0.25); if (idx) { if (m==NextChara) *idx+=i; else *idx-=i; } if (i) { if (num==1) b=GOLDEN_GAMMA*i; else { int mask=num-1,sum=i&mask; if (m==NextChara) sum+=(mod-1)&mask; else sum+=-mod&mask; b=GOLDEN_GAMMA*num*((qword)i/num+(sum>=num)); } if (m==NextChara) { if (i>65536-mod) XorNextState1(xs); for (i=0;i<num;i++) cnt[i]+=b; } else { if (i>mod) XorBackState1(xs); for (i=0;i<num;i++) cnt[i]-=b; } } if (lag) { xgf2_t e,d,s; qword y[XNSD]; a<<=16; e[0]=1; s[0]=2; for (i=1;i<RND_NUM2;i++) e[i]=s[i]=0; for (;;) { if (modf(lag/2.0,&lag)!=0) { xgf2mult(d,e,s); xgf2rem(e,d,m); for (i=0;i<num;i++) cnt[i]+=a; } if (lag==0.0) break; xgf2square(d,s); xgf2rem(s,d,m); a<<=1; } JumpHelp(xs,y,e,1,m==BackChara); } } qword XorNextState1(qword*p) { ENUM_ABCM(0); int i; qword x=p[0],t=x; for (i=0;i<XNSD-1;i++) p[i]=p[i+1]; x<<=a; t^=x; x=t; x>>=b; t^=x; x=p[XNSD-1]; t^=x; x>>=c; t^=x; return p[XNSD-1]=t; } qword XorBackState1(qword*p) { ENUM_ABCM(0); int i; qword t=p[XNSD-1],x=p[XNSD-2]; for (i=XNSD-1;i>0;i--) p[i]=p[i-1]; t^=x; x>>=c; t^=x; x=t; x>>=b; t^=x; x=t; x<<=a; t^=x; p[0]=t; return p[XNSD-1]; } #endif
［前：性能評価］
［ExRand目次］