Browse Source

* Cumulative patch including Bruno's work on large fixnums on 64 bit machines.

master
Richard Kreckel 20 years ago
parent
commit
6fdd87f5e7
  1. 1055
      ChangeLog
  2. 4
      doc/cln.tex
  3. 2
      include/cln/integer.h
  4. 2
      include/cln/numtheory.h
  5. 6
      include/cln/object.h
  6. 5
      include/cln/types.h
  7. 226
      src/base/cl_low.h
  8. 67
      src/base/low/cl_low_div.cc
  9. 2
      src/complex/input/cl_N_read.cc
  10. 2
      src/complex/transcendental/cl_C_expt_C.cc
  11. 13
      src/float/dfloat/elem/cl_DF_from_RA.cc
  12. 8
      src/float/dfloat/elem/cl_DF_scale_I.cc
  13. 2
      src/float/ffloat/conv/cl_RA_to_float.cc
  14. 2
      src/float/ffloat/elem/cl_FF_from_RA.cc
  15. 10
      src/float/ffloat/elem/cl_FF_scale_I.cc
  16. 4
      src/float/lfloat/elem/cl_LF_scale_I.cc
  17. 6
      src/float/output/cl_F_dprint.cc
  18. 2
      src/float/sfloat/elem/cl_SF_from_RA.cc
  19. 10
      src/float/sfloat/elem/cl_SF_scale_I.cc
  20. 6
      src/integer/algebraic/cl_I_sqrtp.cc
  21. 8
      src/integer/bitwise/cl_I_ash_I.cc
  22. 14
      src/integer/bitwise/cl_I_ilength.cc
  23. 8
      src/integer/bitwise/cl_I_log_aux.cc
  24. 2
      src/integer/bitwise/cl_I_logbitp_I.cc
  25. 10
      src/integer/bitwise/cl_I_logcount.cc
  26. 8
      src/integer/bitwise/cl_I_logtest.cc
  27. 150
      src/integer/cl_I.h
  28. 5
      src/integer/conv/cl_I_from_L2.cc
  29. 5
      src/integer/conv/cl_I_from_NDS.cc
  30. 2
      src/integer/conv/cl_I_from_Q.cc
  31. 5
      src/integer/conv/cl_I_from_UL2.cc
  32. 2
      src/integer/conv/cl_I_from_UQ.cc
  33. 101
      src/integer/conv/cl_I_to_L.cc
  34. 2
      src/integer/conv/cl_I_to_Q.cc
  35. 7
      src/integer/conv/cl_I_to_UL.cc
  36. 4
      src/integer/conv/cl_I_to_UQ.cc
  37. 2
      src/integer/conv/cl_I_to_digits.cc
  38. 69
      src/integer/elem/cl_I_div.cc
  39. 42
      src/integer/elem/cl_I_minus.cc
  40. 25
      src/integer/elem/cl_I_mul.cc
  41. 42
      src/integer/elem/cl_I_plus.cc
  42. 20
      src/integer/elem/cl_I_square.cc
  43. 10
      src/integer/elem/cl_I_uminus.cc
  44. 8
      src/integer/gcd/cl_I_gcd.cc
  45. 6
      src/integer/gcd/cl_low_gcd.cc
  46. 2
      src/integer/hash/cl_I_hashcode.cc
  47. 2
      src/integer/input/cl_I_read.cc
  48. 18
      src/integer/misc/cl_I_eqhashcode.cc
  49. 9
      src/integer/misc/cl_I_ord2.cc
  50. 6
      src/integer/misc/cl_I_power2p.cc
  51. 82
      src/integer/misc/combin/cl_I_doublefactorial.cc
  52. 54
      src/integer/misc/combin/cl_I_factorial.cc
  53. 24
      src/modinteger/cl_MI_fix16.h
  54. 24
      src/modinteger/cl_MI_fix29.h
  55. 26
      src/modinteger/cl_MI_fix32.h
  56. 2
      src/modinteger/cl_MI_std.h
  57. 2
      src/numtheory/cl_nt_cornacchia4.cc
  58. 8
      src/numtheory/cl_nt_jacobi.cc
  59. 10
      src/numtheory/cl_nt_jacobi_low.cc
  60. 2
      src/rational/input/cl_RA_read.cc
  61. 2
      src/real/input/cl_R_read.cc
  62. 20
      src/vector/cl_GV_I.cc
  63. 2
      tests/timefact.cc

1055
ChangeLog
File diff suppressed because it is too large
View File

4
doc/cln.tex

@ -1919,7 +1919,7 @@ If @code{x} = 2^(n-1), it returns n. Else it returns 0.
@subsection Number theoretic functions
@table @code
@item uint32 gcd (uint32 a, uint32 b)
@item uint32 gcd (unsigned long a, unsigned long b)
@cindex @code{gcd ()}
@itemx cl_I gcd (const cl_I& a, const cl_I& b)
This function returns the greatest common divisor of @code{a} and @code{b},
@ -1948,7 +1948,7 @@ normalized to be >= 0.
rational number, this function returns true and sets *l = log(a,b), else
it returns false.
@item int jacobi (sint32 a, sint32 b)
@item int jacobi (signed long a, signed long b)
@cindex @code{jacobi()}
@itemx int jacobi (const cl_I& a, const cl_I& b)
Returns the Jacobi symbol

2
include/cln/integer.h

@ -409,7 +409,7 @@ struct cl_I_div_t {
// > a,b: zwei Integers
// < ergebnis: (gcd a b), ein Integer >=0
extern const cl_I gcd (const cl_I& a, const cl_I& b);
extern uint32 gcd (uint32 a, uint32 b);
extern uintV gcd (uintV a, uintV b);
// Liefert den ggT zweier Integers samt Beifaktoren.
// g = xgcd(a,b,&u,&v)

2
include/cln/numtheory.h

@ -15,7 +15,7 @@ namespace cln {
// ( --- )
// ( b )
// a, b must be integers, b > 0, b odd. The result is 0 iff gcd(a,b) > 1.
extern int jacobi (sint32 a, sint32 b);
extern int jacobi (sintV a, sintV b);
extern int jacobi (const cl_I& a, const cl_I& b);
// isprobprime(n), n integer > 0,

6
include/cln/object.h

@ -92,11 +92,7 @@ inline cl_boolean cl_immediate_p (cl_uint word)
#define cl_tag_len 3
#endif
#define cl_tag_shift 0
#if (cl_pointer_size == 64)
#define cl_value_shift 32
#else
#define cl_value_shift (cl_tag_len+cl_tag_shift)
#endif
#define cl_value_shift (cl_tag_len+cl_tag_shift)
#define cl_value_len (cl_pointer_size - cl_value_shift)
#define cl_tag_mask (((1UL << cl_tag_len) - 1) << cl_tag_shift)
#define cl_value_mask (((1UL << cl_value_len) - 1) << cl_value_shift)

5
include/cln/types.h

@ -101,6 +101,11 @@
typedef long sintP;
typedef unsigned long uintP;
// Integer type used for the value of a fixnum.
#define intVsize long_bitsize
typedef long sintV;
typedef unsigned long uintV;
// Numbers in the heap are stored as "digit" sequences.
// A digit is an unsigned int with intDsize bits.
// intDsize should be 8 or 16 or 32 or 64 and it should match mp_limb_t,

226
src/base/cl_low.h

@ -236,13 +236,13 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
})
#elif defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
#define mulu32(x,y,hi_zuweisung,lo_zuweisung) \
({ var register uint64 _hi; \
var register uint64 _lo; \
__asm__("umul %2,%3,%1\n\trd %y,%0" \
: "=r" (_hi), "=r" (_lo) \
({ var register uint64 _prod; \
__asm__("umul %1,%2,%0" \
: "=r" (_prod) \
: "r" ((uint32)(x)), "r" ((uint32)(y)) \
); \
hi_zuweisung (uint32)_hi; lo_zuweisung (uint32)_lo; \
hi_zuweisung (uint32)(_prod>>32); \
lo_zuweisung (uint32)(_prod); \
})
#elif defined(__GNUC__) && defined(__sparc__) && !defined(NO_ASM)
#define mulu32(x,y,hi_zuweisung,lo_zuweisung) \
@ -314,7 +314,17 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
// mulu32_w(arg1,arg2)
// > arg1, arg2 : zwei 32-Bit-Zahlen
// < result : eine 64-Bit-Zahl
#if defined(__GNUC__)
#if defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
// Prefer the umul instruction over the mulx instruction (overkill).
#define mulu32_w(x,y) \
({ var register uint64 _prod; \
__asm__("umul %1,%2,%0" \
: "=r" (_prod) \
: "r" ((uint32)(x)), "r" ((uint32)(y)) \
); \
_prod; \
})
#elif defined(__GNUC__)
#define mulu32_w(x,y) ((uint64)(uint32)(x) * (uint64)(uint32)(y))
#else
extern "C" uint64 mulu32_w (uint32 arg1, uint32 arg2);
@ -356,7 +366,7 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
var register uint64 _lo; \
__asm__("mulq %2" \
: "=d" /* %rdx */ (_hi), "=a" /* %rax */ (_lo) \
: "g" ((uint64)(x)), "1" /* %rax */ ((uint64)(y)) \
: "rm" ((uint64)(x)), "1" /* %rax */ ((uint64)(y)) \
); \
hi_zuweisung _hi; lo_zuweisung _lo; \
})
@ -468,9 +478,9 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
})
#elif defined(__GNUC__) && defined(__arm__) && 0 // see comment cl_asm_arm.cc
#define divu_3216_1616(x,y,q_zuweisung,r_zuweisung) \
{ var uint32 _q = divu_3216_1616_(x,y); /* extern in Assembler */ \
var register uint32 _r __asm__("%r1"/*"%a2"*/); \
q_zuweisung _q; r_zuweisung _r; \
{ var uint32 __q = divu_3216_1616_(x,y); /* extern in Assembler */ \
var register uint32 __r __asm__("%r1"/*"%a2"*/); \
q_zuweisung __q; r_zuweisung __r; \
}
#elif defined(__GNUC__) && !defined(__arm__)
#define divu_3216_1616(x,y,q_zuweisung,r_zuweisung) \
@ -768,13 +778,54 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
// < uint32 q: floor(x/y)
// < uint32 r: x mod y
// < x = q*y+r
#if defined(__GNUC__)
#if defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
// Prefer the udiv and umul instructions over the udivx and mulx instructions
// (overkill).
#define divu_6432_3232_w(x,y,q_zuweisung,r_zuweisung) \
({var uint64 __x = (x); \
var uint32 __xhi = high32(__x); \
var uint32 __xlo = low32(__x); \
var uint32 __y = (y); \
var uint64 __q; \
var uint64 __r; \
__asm__ __volatile__ ( \
"wr %2,%%g0,%%y\n\t" \
"udiv %3,%4,%0\n\t" \
"umul %0,%4,%1" \
"sub %3,%1,%1" \
: "=&r" (__q), "=&r" (__r) \
: "r" (__xhi), "r" (__xlo), "r" (__y)); \
q_zuweisung (uint32)__q; \
r_zuweisung (uint32)__r; \
})
#elif defined(__GNUC__) && (defined(__alpha__) || defined(__ia64__) || defined(__mips64__) || defined(__sparc64__))
// On __alpha__, computing the remainder by multiplication is just two
// instructions, compared to the __remqu (libc) function call for the %
// operator.
// On __ia64__, computing the remainder by multiplication is just four
// instructions, compared to the __umoddi3 (libgcc) function call for the %
// operator.
// On __mips64__, computing the remainder by multiplication is just two
// instructions, compared to the __umoddi3 (libgcc) function call for the %
// operator.
// On __sparc64__, computing the remainder by multiplication uses a 32-bit
// multiplication instruction, compared to a 64-bit multiplication when the %
// operator is used.
#define divu_6432_3232_w(x,y,q_zuweisung,r_zuweisung) \
({var uint64 __x = (x); \
var uint32 __y = (y); \
var uint32 __q = floor(__x,(uint64)__y); \
q_zuweisung __q; r_zuweisung (uint32)__x - __q * __y; \
})
#elif defined(__GNUC__) && defined(__x86_64__)
// On __x86_64__, gcc 4.0 performs both quotient and remainder computation
// in a single instruction.
#define divu_6432_3232_w(x,y,q_zuweisung,r_zuweisung) \
({var uint64 __x = (x); \
var uint32 __y = (y); \
var uint32 __q = floor(__x,(uint64)__y); \
q_zuweisung __q; r_zuweisung __x % (uint64)__y; \
})
#else
#define divu_6432_3232_w(x,y,q_zuweisung,r_zuweisung) \
{ var uint64 __x = (x); \
@ -782,22 +833,125 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
}
#endif
// Dividiert eine 64-Bit-Zahl durch eine 32-Bit-Zahl und
// liefert einen 64-Bit-Quotienten und einen 32-Bit-Rest.
// divu_6432_6432(x,y,q=,r=);
// > uint64 x: Zähler
// > uint32 y: Nenner
// > Es sei bekannt, daß y>0.
// < uint64 q: floor(x/y)
// < uint32 r: x mod y
// < x = q*y+r
#if defined(__GNUC__) && (defined(__alpha__) || defined(__ia64__) || defined(__mips64__) || defined(__sparc64__))
// On __alpha__, computing the remainder by multiplication is just two
// instructions, compared to the __remqu (libc) function call for the %
// operator.
// On __ia64__, computing the remainder by multiplication is just four
// instructions, compared to the __umoddi3 (libgcc) function call for the %
// operator.
// On __mips64__, computing the remainder by multiplication is just two
// instructions, compared to the __umoddi3 (libgcc) function call for the %
// operator.
// On __sparc64__, computing the remainder by multiplication uses a 32-bit
// multiplication instruction, compared to a 64-bit multiplication when the %
// operator is used.
#define divu_6432_6432(x,y,q_zuweisung,r_zuweisung) \
({var uint64 _x = (x); \
var uint32 _y = (y); \
var uint64 _q; \
q_zuweisung _q = floor(_x,(uint64)_y); \
r_zuweisung low32(_x) - low32(_q) * _y; \
})
#elif defined(__GNUC__) && defined(__x86_64__)
// On __x86_64__, gcc 4.0 performs both quotient and remainder computation
// in a single instruction.
#define divu_6432_6432(x,y,q_zuweisung,r_zuweisung) \
({var uint64 _x = (x); \
var uint32 _y = (y); \
q_zuweisung floor(_x,(uint64)_y); \
r_zuweisung _x % (uint64)_y; \
})
#else
// Methode: (beta = 2^32)
// x = x1*beta+x0 schreiben.
// Division mit Rest: x1 = q1*y + r1, wobei 0 <= x1 < beta <= beta*y.
// Also 0 <= q1 < beta, 0 <= r1 < y.
// Division mit Rest: (r1*beta+x0) = q0*y + r0, wobei 0 <= r1*beta+x0 < beta*y.
// Also 0 <= q0 < beta, 0 <= r0 < y
// und x = x1*beta+x0 = (q1*beta+q0)*y + r0.
// Setze q := q1*beta+q0 und r := r0.
#if defined(__GNUC__)
#define divu_6432_6432(x,y,q_zuweisung,r_zuweisung) \
({var uint64 _x = (x); \
var uint32 _y = (y); \
var uint32 _q1; \
var uint32 _q0; \
var uint32 _r1; \
divu_6432_3232(0,high32(_x),_y, _q1 = , _r1 = ); \
divu_6432_3232(_r1,low32(_x),_y, _q0 = , r_zuweisung); \
q_zuweisung highlow64(_q1,_q0); \
})
#else
#define divu_6432_6432(x,y,q_zuweisung,r_zuweisung) \
{var uint64 _x = (x); \
var uint32 _y = (y); \
var uint32 _q1; \
var uint32 _q0; \
var uint32 _r1; \
divu_6432_3232(0,high32(_x),_y, _q1 = , _r1 = ); \
divu_6432_3232(_r1,low32(_x),_y, _q0 = , r_zuweisung); \
q_zuweisung highlow64(_q1,_q0); \
}
#endif
#endif
// Dividiert eine 64-Bit-Zahl durch eine 64-Bit-Zahl und
// liefert einen 64-Bit-Quotienten und einen 64-Bit-Rest.
// divu_6464_6464(x,y,q=,r=);
// > uint64 x: Zähler
// > uint64 y: Nenner
// Es sei bekannt, daß y>0.
// > Es sei bekannt, daß y>0.
// < uint64 q: floor(x/y)
// < uint64 r: x mod y
// < x = q*y+r
#if defined(__alpha__) || 1
#if defined(__GNUC__) && (defined(__alpha__) || defined(__ia64__) || defined(__mips64__) || defined(__sparc64__))
// On __alpha__, computing the remainder by multiplication is just two
// instructions, compared to the __remqu (libc) function call for the %
// operator.
// On __ia64__, computing the remainder by multiplication is just four
// instructions, compared to the __umoddi3 (libgcc) function call for the %
// operator.
// On __mips64__, computing the remainder by multiplication is just two
// instructions, compared to the __umoddi3 (libgcc) function call for the %
// operator.
// On __sparc64__, it doesn't matter.
#define divu_6464_6464(x,y,q_zuweisung,r_zuweisung) \
{ var uint64 __x = (x); \
var uint64 __y = (y); \
q_zuweisung (__x / __y); \
r_zuweisung (__x % __y); \
}
({var uint64 _x = (x); \
var uint64 _y = (y); \
var uint64 _q; \
q_zuweisung _q = floor(_x,_y); \
r_zuweisung _x - _q * _y; \
})
#elif defined(__GNUC__) && (defined(__sparc64__) || defined(__x86_64__))
// On __sparc64__, it doesn't matter.
// On __x86_64__, gcc 4.0 performs both quotient and remainder computation
// in a single instruction.
#define divu_6464_6464(x,y,q_zuweisung,r_zuweisung) \
({var uint64 _x = (x); \
var uint64 _y = (y); \
q_zuweisung floor(_x,_y); \
r_zuweisung _x % _y; \
})
#else
// For unknown CPUs, we don't know whether gcc's __udivdi3 function plus a
// multiplication is slower or faster than our own divu_6464_6464_ routine.
// Anyway, call our own routine.
extern "C" uint64 divu_6464_6464_ (uint64 x, uint64 y); // -> Quotient q
extern "C" uint64 divu_64_rest; // -> Rest r
#define divu_6464_6464(x,y,q_zuweisung,r_zuweisung) \
{ q_zuweisung divu_6464_6464_(x,y); r_zuweisung divu_64_rest; }
#define NEED_VAR_divu_64_rest
#define NEED_FUNCTION_divu_6464_6464_
#endif
// Dividiert eine 128-Bit-Zahl durch eine 64-Bit-Zahl und
@ -831,6 +985,7 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
#else
#define divu_12864_6464(xhi,xlo,y,q_zuweisung,r_zuweisung) \
{ q_zuweisung divu_12864_6464_(xhi,xlo,y); r_zuweisung divu_64_rest; }
#define NEED_VAR_divu_64_rest
#define NEED_FUNCTION_divu_12864_6464_
#endif
@ -1230,10 +1385,11 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
/* _x32 hat höchstens 4 Bits. */\
if (_x32 >= bit(2)) { _x32 = _x32>>2; _bitsize += 2; } \
/* _x32 hat höchstens 2 Bits. */\
if (_x32 >= bit(1)) { /* _x32 = _x32>>1; */ _bitsize += 1; } \
if (_x32 >= bit(1)) { /* _x32 = _x32>>1; */ _bitsize += 1; } \
/* _x32 hat höchstens 1 Bit. Dieses Bit muß gesetzt sein. */\
size_zuweisung _bitsize; \
}
#define GENERIC_INTEGERLENGTH32
#endif
// Bits einer 64-Bit-Zahl zählen:
@ -1241,9 +1397,10 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
// setzt size auf die höchste in digit vorkommende Bitnummer.
// > digit: ein uint64 >0
// < size: >0, <=64, mit 2^(size-1) <= digit < 2^size
#ifdef GENERIC_INTEGERLENGTH32
#define integerlength64(digit,size_zuweisung) \
{ var uintC _bitsize = 1; \
var uint64 _x64 = (uint64)(digit); \
var uint64 _x64 = (uint64)(digit); \
/* _x64 hat höchstens 64 Bits. */\
if (_x64 >= bit(32)) { _x64 = _x64>>32; _bitsize += 32; } \
/* _x64 hat höchstens 32 Bits. */\
@ -1255,10 +1412,23 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
/* _x64 hat höchstens 4 Bits. */\
if (_x64 >= bit(2)) { _x64 = _x64>>2; _bitsize += 2; } \
/* _x64 hat höchstens 2 Bits. */\
if (_x64 >= bit(1)) { /* _x64 = _x64>>1; */ _bitsize += 1; } \
if (_x64 >= bit(1)) { /* _x64 = _x64>>1; */ _bitsize += 1; } \
/* _x64 hat höchstens 1 Bit. Dieses Bit muß gesetzt sein. */\
size_zuweisung _bitsize; \
}
#else
#define integerlength64(digit,size_zuweisung) \
{ var uint64 _x64 = (digit); \
var uintC _bitsize64 = 0; \
var uint32 _x32_from_integerlength64; \
if (_x64 >= (1ULL << 32)) { \
_x32_from_integerlength64 = _x64>>32; _bitsize64 += 32; \
} else { \
_x32_from_integerlength64 = _x64; \
} \
integerlength32(_x32_from_integerlength64, size_zuweisung _bitsize64 + ); \
}
#endif
// Hintere Nullbits eines 32-Bit-Wortes zählen:
// ord2_32(digit,count=);
@ -1288,11 +1458,23 @@ inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
// Sei n = ord2(x). Dann ist logxor(x,x-1) = 2^n + (2^n-1) = 2^(n+1)-1.
// Also (ord2 x) = (1- (integer-length (logxor x (1- x)))) .
#define ord2_32(digit,count_zuweisung) \
{ var uint32 _digit = digit ^ (digit - 1); \
{ var uint32 _digit = (digit) ^ ((digit) - 1); \
integerlength32(_digit,count_zuweisung -1 + ) \
}
#endif
// Hintere Nullbits eines 64-Bit-Wortes zählen:
// ord2_64(digit,count=);
// setzt size auf die kleinste in digit vorkommende Bitnummer.
// > digit: ein uint64 >0
// < count: >=0, <64, mit 2^count | digit, digit/2^count ungerade
// Sei n = ord2(x). Dann ist logxor(x,x-1) = 2^n + (2^n-1) = 2^(n+1)-1.
// Also (ord2 x) = (1- (integer-length (logxor x (1- x)))) .
#define ord2_64(digit,count_zuweisung) \
{ var uint64 _digit = (digit) ^ ((digit) - 1); \
integerlength64(_digit,count_zuweisung -1 + ) \
}
// Bits eines Wortes zählen.
// logcount_NN();

67
src/base/low/cl_low_div.cc

@ -210,8 +210,73 @@ uint32 divu_6432_3232_(uint32 xhi, uint32 xlo, uint32 y)
uint64 divu_64_rest;
#endif
#ifdef NEED_FUNCTION_divu_6464_6464_
namespace cln {
uint64 divu_6464_6464_(uint64 x, uint64 y)
// Methode: (beta = 2^n = 2^32, n = 32)
// Falls y < beta, handelt es sich um eine 64-durch-32-Bit-Division.
// Falls y >= beta:
// Quotient q = floor(x/y) < beta (da 0 <= x < beta^2, y >= beta).
// y habe genau n+k Bits (1 <= k <= n), d.h. 2^(n+k-1) <= y < 2^(n+k).
// Schreibe x = 2^k*x1 + x0 mit x1 := floor(x/2^k)
// und y = 2^k*y1 + y0 mit y1 := floor(y/2^k)
// und bilde den Näherungs-Quotienten floor(x1/y1)
// oder (noch besser) floor(x1/(y1+1)).
// Wegen 0 <= x1 < 2^(2n) und 0 < 2^(n-1) <= y1 < 2^n
// und x1/(y1+1) <= x/y < x1/(y1+1) + 2
// (denn x1/(y1+1) = (x1*2^k)/((y1+1)*2^k) <= (x1*2^k)/y <= x/y
// und x/y - x1/(y1+1) = (x+x*y1-x1*y)/(y*(y1+1))
// = (x+x0*y1-x1*y0)/(y*(y1+1)) <= (x+x0*y1)/(y*(y1+1))
// <= x/(y*(y1+1)) + x0/y
// <= 2^(2n)/(2^(n+k-1)*(2^(n-1)+1)) + 2^k/2^(n+k-1)
// = 2^(n-k+1)/(2^(n-1)+1) + 2^(1-n) <= 2^n/(2^(n-1)+1) + 2^(1-n) < 2 )
// gilt floor(x1/(y1+1)) <= floor(x/y) <= floor(x1/(y1+1)) + 2 .
// Man bildet also q:=floor(x1/(y1+1)) (ein Shift um n Bit oder
// eine (2n)-durch-n-Bit-Division, mit Ergebnis q <= floor(x/y) < beta)
// und x-q*y und muss hiervon noch höchstens 2 mal y abziehen und q
// incrementieren, um den Quotienten q = floor(x/y) und den Rest
// x-floor(x/y)*y der Division zu bekommen.
{
if (y <= (uint64)(((uint64)1<<32)-1))
{ var uint32 q1;
var uint32 q0;
var uint32 r1;
divu_6432_3232(0,high32(x),y, q1 = , r1 = );
divu_6432_3232(r1,low32(x),y, q0 = , divu_64_rest = );
return highlow64(q1,q0);
}
else
{ var uint64 x1 = x; // x1 := x
var uint64 y1 = y; // y1 := y
var uint32 q;
do { x1 = floor(x1,2); y1 = floor(y1,2); } // k erhöhen
while (!(y1 <= (uint64)(((uint64)1<<32)-1))); // bis y1 < beta
{ var uint32 y2 = low32(y1)+1; // y1+1 bilden
if (y2==0)
{ q = high32(x1); } // y1+1=beta -> ein Shift
else
{ divu_6432_3232(high32(x1),low32(x1),y2,q=,); } // Division von x1 durch y1+1
}
// q = floor(x1/(y1+1))
// x-q*y bilden (eine 32-mal-64-Bit-Multiplikation ohne Überlauf):
x -= highlow64_0(mulu32_64(q,high32(y))); // q * high32(y) * beta
// gefahrlos, da q*high32(y) <= q*y/beta <= x/beta < beta
x -= mulu32_64(q,low32(y)); // q * low32(y)
// gefahrlos, da q*high32(y)*beta + q*low32(y) = q*y <= x
// Noch höchstens 2 mal y abziehen:
if (x >= y)
{ q += 1; x -= y;
if (x >= y)
{ q += 1; x -= y;
} }
divu_64_rest = x;
return (uint64)q;
}
}
} // namespace cln
#endif
#ifdef NEED_FUNCTION_divu_12864_6464_
uint64 divu_64_rest;
namespace cln {
uint64 divu_12864_6464_(uint64 xhi, uint64 xlo, uint64 y)
// Methode:

2
src/complex/input/cl_N_read.cc

@ -97,7 +97,7 @@ const cl_N read_complex (const cl_read_flags& flags, const char * string, const
fprint(std::cerr, "\n");
cl_abort();
}
rational_base = FN_to_UL(base); ptr = base_end_ptr;
rational_base = FN_to_UV(base); ptr = base_end_ptr;
break;
}
ptr++;

2
src/complex/transcendental/cl_C_expt_C.cc

@ -166,7 +166,7 @@ const cl_N expt (const cl_N& x, const cl_N& y)
}
}
if (fixnump(m) && fixnump(n)) { // |m| und n klein?
var uintL _n = FN_to_UL(n);
var uintV _n = FN_to_UV(n);
if ((_n & (_n-1)) == 0) { // n Zweierpotenz?
Mutable(cl_N,x);
until ((_n = _n >> 1) == 0)

13
src/float/dfloat/elem/cl_DF_from_RA.cc

@ -66,10 +66,15 @@ const cl_DF cl_RA_to_DF (const cl_RA& x)
var cl_I_div_t q_r = cl_divide(zaehler,nenner);
var cl_I& q = q_r.quotient;
var cl_I& r = q_r.remainder;
// 2^53 <= q < 2^55, also ist q Bignum mit ceiling(55/intDsize) Digits.
var const uintD* ptr = BN_MSDptr(q);
#if (cl_word_size==64)
#if (cl_value_len>=55)
// 2^53 <= q < 2^55: q is fixnum.
var uint64 mant = FN_to_UV(q);
#else
// 2^53 <= q < 2^55: q is bignum with one Digit.
var const uintD* ptr = BN_MSDptr(q);
var uint64 mant = get_max64_Dptr(55,ptr);
#endif
if (mant >= bit(DF_mant_len+2))
// 2^54 <= q < 2^55, schiebe um 2 Bits nach rechts
{ var uint64 rounding_bits = mant & (bit(2)-1);
@ -106,7 +111,9 @@ const cl_DF cl_RA_to_DF (const cl_RA& x)
ab:
// Fertig.
return encode_DF(sign,lendiff,mant);
#else
#else // (cl_word_size<64)
// 2^53 <= q < 2^55: q is bignum with two Digits.
var const uintD* ptr = BN_MSDptr(q);
var uint32 manthi = get_max32_Dptr(23,ptr);
var uint32 mantlo = get_32_Dptr(ptr mspop ceiling(23,intDsize));
if (manthi >= bit(DF_mant_len-32+2))

8
src/float/dfloat/elem/cl_DF_scale_I.cc

@ -34,9 +34,9 @@ const cl_DF scale_float (const cl_DF& x, const cl_I& delta)
#endif
if (!minusp(delta))
// delta>=0
{ var uintL udelta;
{ var uintV udelta;
if (fixnump(delta)
&& ((udelta = FN_to_L(delta)) <= (uintL)(DF_exp_high-DF_exp_low))
&& ((udelta = FN_to_V(delta)) <= (uintV)(DF_exp_high-DF_exp_low))
)
{ exp = exp+udelta;
#if (cl_word_size==64)
@ -52,8 +52,8 @@ const cl_DF scale_float (const cl_DF& x, const cl_I& delta)
// delta<0
{ var uintL udelta;
if (fixnump(delta)
&& ((udelta = -FN_to_L(delta)) <= (uintL)(DF_exp_high-DF_exp_low))
&& ((cl_value_len+1<intLsize) || !(udelta==0))
&& ((udelta = -FN_to_V(delta)) <= (uintV)(DF_exp_high-DF_exp_low))
&& ((cl_value_len+1<intVsize) || !(udelta==0))
)
{ exp = exp-udelta;
#if (cl_word_size==64)

2
src/float/ffloat/conv/cl_RA_to_float.cc

@ -58,7 +58,7 @@ float float_approx (const cl_RA& x)
var cl_I& r = q_r.remainder;
// 2^24 <= q < 2^26, also ist q Fixnum oder Bignum mit bn_minlength Digits.
var uint32 mant = ((FF_mant_len+3 < cl_value_len)
? FN_to_UL(q)
? FN_to_UV(q)
: cl_I_to_UL(q)
);
if (mant >= bit(FF_mant_len+2))

2
src/float/ffloat/elem/cl_FF_from_RA.cc

@ -68,7 +68,7 @@ const cl_FF cl_RA_to_FF (const cl_RA& x)
var cl_I& r = q_r.remainder;
// 2^24 <= q < 2^26, also ist q Fixnum oder Bignum mit bn_minlength Digits.
var uint32 mant = ((FF_mant_len+3 < cl_value_len)
? FN_to_UL(q)
? FN_to_UV(q)
: cl_I_to_UL(q)
);
if (mant >= bit(FF_mant_len+2))

10
src/float/ffloat/elem/cl_FF_scale_I.cc

@ -28,9 +28,9 @@ const cl_FF scale_float (const cl_FF& x, const cl_I& delta)
FF_decode(x, { return x; }, sign=,exp=,mant=);
if (!minusp(delta))
// delta>=0
{ var uintL udelta;
{ var uintV udelta;
if (fixnump(delta)
&& ((udelta = FN_to_L(delta)) <= (uintL)(FF_exp_high-FF_exp_low))
&& ((udelta = FN_to_V(delta)) <= (uintV)(FF_exp_high-FF_exp_low))
)
{ exp = exp+udelta;
return encode_FF(sign,exp,mant);
@ -40,10 +40,10 @@ const cl_FF scale_float (const cl_FF& x, const cl_I& delta)
}
else
// delta<0
{ var uintL udelta;
{ var uintV udelta;
if (fixnump(delta)
&& ((udelta = -FN_to_L(delta)) <= (uintL)(FF_exp_high-FF_exp_low))
&& ((cl_value_len+1<intLsize) || !(udelta==0))
&& ((udelta = -FN_to_V(delta)) <= (uintV)(FF_exp_high-FF_exp_low))
&& ((cl_value_len+1<intVsize) || !(udelta==0))
)
{ exp = exp-udelta;
return encode_FF(sign,exp,mant);

4
src/float/lfloat/elem/cl_LF_scale_I.cc

@ -26,11 +26,11 @@ const cl_LF scale_float (const cl_LF& x, const cl_I& delta)
if (eq(delta,0)) { return x; } // delta=0 -> x als Ergebnis
var uintL uexp = TheLfloat(x)->expo;
if (uexp==0) { return x; }
var uintL udelta;
var uintV udelta;
// |delta| muß <= LF_exp_high-LF_exp_low < 2^32 sein. Wie bei I_to_UL:
if (fixnump(delta)) {
// Fixnum
var sint32 sdelta = FN_to_L(delta);
var sintV sdelta = FN_to_V(delta);
if (sdelta >= 0)
{ udelta = sdelta; goto pos; }
else

6
src/float/output/cl_F_dprint.cc

@ -390,13 +390,13 @@ void print_float (std::ostream& stream, const cl_print_float_flags& flags, const
// erst Null und Punkt, dann -expo Nullen, dann alle Ziffern
fprintchar(stream,'0');
fprintchar(stream,'.');
for (uintL i = -FN_to_L(expo); i > 0; i--)
for (uintV i = -FN_to_V(expo); i > 0; i--)
fprintchar(stream,'0');
fprint(stream,mantstring);
expo = 0; // auszugebender Exponent ist 0
} else {
// "fixed-point notation" mit expo > 0 oder "scientific notation"
var uintL scale = (flag ? FN_to_L(expo) : 1);
var uintV scale = (flag ? FN_to_V(expo) : 1);
// Der Dezimalpunkt wird um scale Stellen nach rechts geschoben,
// d.h. es gibt scale Vorkommastellen. scale > 0.
if (scale < mantlen) {
@ -412,7 +412,7 @@ void print_float (std::ostream& stream, const cl_print_float_flags& flags, const
// scale>=mantlen -> es bleibt nichts für die Nachkommastellen.
// alle Ziffern, dann scale-mantlen Nullen, dann Punkt und Null
fprint(stream,mantstring);
for (uintL i = mantlen; i < scale; i++)
for (uintV i = mantlen; i < scale; i++)
fprintchar(stream,'0');
fprintchar(stream,'.');
fprintchar(stream,'0');

2
src/float/sfloat/elem/cl_SF_from_RA.cc

@ -66,7 +66,7 @@ const cl_SF cl_RA_to_SF (const cl_RA& x)
var cl_I& q = q_r.quotient;
var cl_I& r = q_r.remainder;
// 2^17 <= q < 2^19, also ist q Fixnum.
var uint32 mant = FN_to_UL(q);
var uint32 mant = FN_to_UV(q);
if (mant >= bit(SF_mant_len+2))
// 2^18 <= q < 2^19, schiebe um 2 Bits nach rechts
{ var uintL rounding_bits = mant & (bit(2)-1);

10
src/float/sfloat/elem/cl_SF_scale_I.cc

@ -28,9 +28,9 @@ const cl_SF scale_float (const cl_SF& x, const cl_I& delta)
SF_decode(x, { return x; }, sign=,exp=,mant=);
if (!minusp(delta))
// delta>=0
{ var uintL udelta;
{ var uintV udelta;
if (fixnump(delta)
&& ((udelta = FN_to_L(delta)) <= (uintL)(SF_exp_high-SF_exp_low))
&& ((udelta = FN_to_V(delta)) <= (uintV)(SF_exp_high-SF_exp_low))
)
{ exp = exp+udelta;
return encode_SF(sign,exp,mant);
@ -40,10 +40,10 @@ const cl_SF scale_float (const cl_SF& x, const cl_I& delta)
}
else
// delta<0
{ var uintL udelta;
{ var uintV udelta;
if (fixnump(delta)
&& ((udelta = -FN_to_L(delta)) <= (uintL)(SF_exp_high-SF_exp_low))
&& ((cl_value_len+1<intLsize) || !(udelta==0))
&& ((udelta = -FN_to_V(delta)) <= (uintV)(SF_exp_high-SF_exp_low))
&& ((cl_value_len+1<intVsize) || !(udelta==0))
)
{ exp = exp-udelta;
return encode_SF(sign,exp,mant);

6
src/integer/algebraic/cl_I_sqrtp.cc

@ -63,17 +63,17 @@ cl_boolean sqrtp (const cl_I& x, cl_I* w)
}
// Check mod 63.
{ var cl_I_div_t div63 = cl_divide(x,L_to_FN(63));
if (!squares_mod_63[FN_to_UL(div63.remainder)])
if (!squares_mod_63[FN_to_UV(div63.remainder)])
{ return cl_false; } // not a square mod 63 -> not a square
}
// Check mod 65.
{ var cl_I_div_t div65 = cl_divide(x,L_to_FN(65));
if (!squares_mod_65[FN_to_UL(div65.remainder)])
if (!squares_mod_65[FN_to_UV(div65.remainder)])
{ return cl_false; } // not a square mod 65 -> not a square
}
// Check mod 11.
{ var cl_I_div_t div11 = cl_divide(x,L_to_FN(11));
if (!squares_mod_11[FN_to_UL(div11.remainder)])
if (!squares_mod_11[FN_to_UV(div11.remainder)])
{ return cl_false; } // not a square mod 11 -> not a square
}
// Check with full precision.

8
src/integer/bitwise/cl_I_ash_I.cc

@ -35,7 +35,7 @@ const cl_I ash (const cl_I& x, const cl_I& y)
if (!minusp(y)) {
// y>=0
var uintL i; // i = y mod intDsize, >=0, <intDsize
var uintL k; // k = y div intDsize, >=0, <2^intCsize
var cl_uint k; // k = y div intDsize, >=0, <2^intCsize
if (bignump(y)) {
#if (log2_intDsize+intCsize <= cl_value_len-1)
// y >= 2^(cl_value_len-1) >= intDsize*2^intCsize
@ -78,7 +78,7 @@ const cl_I ash (const cl_I& x, const cl_I& y)
i = lspref(arrayLSDptr(bn->data,len),0) % intDsize;
#endif
} else {
var uintL y_ = FN_to_L(y); // Wert von y, >=0, <intDsize*2^intCsize
var uintV y_ = FN_to_V(y); // Wert von y, >=0, <intDsize*2^intCsize
i = y_%intDsize;
k = floor(y_,intDsize);
}
@ -111,7 +111,7 @@ const cl_I ash (const cl_I& x, const cl_I& y)
} else {
// y<0
var uintL i; // i = (-y) mod intDsize, >=0, <intDsize
var uintL k; // k = (-y) div intDsize, >=0, <2^intCsize
var cl_uint k; // k = (-y) div intDsize, >=0, <2^intCsize
if (bignump(y)) {
#if (log2_intDsize+intCsize <= cl_value_len-1)
// -y-1 >= 2^(cl_value_len-1) >= intDsize*2^intCsize
@ -158,7 +158,7 @@ const cl_I ash (const cl_I& x, const cl_I& y)
goto sign;
#endif
} else {
var uintL y_ = -FN_to_L(y); // Wert von -y, >0, <intDsize*2^intCsize
var uintV y_ = -FN_to_V(y); // Wert von -y, >0, <intDsize*2^intCsize
i = y_%intDsize;
k = floor(y_,intDsize);
}

14
src/integer/bitwise/cl_I_ilength.cc

@ -18,10 +18,16 @@ uintL integer_length (const cl_I& x)
{
if (fixnump(x))
{ var uintL bitcount = 0;
var uint32 x_ = FN_to_L(x); // x als 32-Bit-Zahl
if (FN_L_minusp(x,(sint32)x_)) { x_ = ~ x_; } // falls <0, komplementieren
if (!(x_==0)) { integerlength32(x_,bitcount=); }
return bitcount; // 0 <= bitcount < 32.
var uintV x_ = FN_to_V(x); // x als intVsize-Bit-Zahl
if (FN_V_minusp(x,(sintV)x_)) { x_ = ~ x_; } // falls <0, komplementieren
if (!(x_==0)) {
#if (intVsize>32)
integerlength64(x_,bitcount=);
#else
integerlength32(x_,bitcount=);
#endif
}
return bitcount; // 0 <= bitcount < intVsize.
}
else
{ var const uintD* MSDptr;

8
src/integer/bitwise/cl_I_log_aux.cc

@ -23,8 +23,8 @@ uintD* I_to_DS_n_aux (const cl_I& obj, uintC n, uintD* destptr)
#if (intDsize==64) // && (FN_maxlength==1)
lsprefnext(destptr) = FN_to_Q(obj);
#else // (intDsize<=32)
var uint32 wert = FN_to_L(obj);
#define FN_maxlength_a (intLsize/intDsize)
var uintV wert = FN_to_V(obj);
#define FN_maxlength_a (intVsize/intDsize)
#define FN_maxlength_b (FN_maxlength<=FN_maxlength_a ? FN_maxlength : FN_maxlength_a)
// FN_maxlength Digits ablegen. Davon kann man FN_maxlength_b Digits aus wert nehmen.
#if (FN_maxlength_b > 1)
@ -34,9 +34,9 @@ uintD* I_to_DS_n_aux (const cl_I& obj, uintC n, uintD* destptr)
#endif
lsprefnext(destptr) = (uintD)wert;
#if (FN_maxlength > FN_maxlength_b)
// Es ist oint_data_len = intLsize, brauche
// Es ist cl_value_len-1 = intVsize, brauche
// noch FN_maxlength-FN_maxlength_b = 1 Digit.
lsprefnext(destptr) = (sintD)sign_of(FN_to_L(obj));
lsprefnext(destptr) = (sintD)sign_of(FN_to_V(obj));
#endif
#endif
n -= FN_maxlength;

2
src/integer/bitwise/cl_I_logbitp_I.cc

@ -25,7 +25,7 @@ cl_boolean logbitp (const cl_I& x, const cl_I& y)
// Sonst x=intDsize*k+i, Teste Bit i vom Worte Nr. k+1 (von oben herab).
if (!minusp(x)) // x>=0 ?
{ if (fixnump(x))
{ var uintL x_ = FN_to_L(x);
{ var uintV x_ = FN_to_V(x);
var uintC ylen;
var const uintD* yLSDptr;
I_to_NDS_nocopy(y, ,ylen=,yLSDptr=,cl_true, { return cl_false; } ); // DS zu y

10
src/integer/bitwise/cl_I_logcount.cc

@ -19,9 +19,15 @@ namespace cln {
uintL logcount (const cl_I& x)
{
if (fixnump(x))
{ var uint32 x32 = FN_to_L(x); // x als 32-Bit-Zahl
if (FN_L_minusp(x,(sint32)x32)) { x32 = ~ x32; } // falls <0, komplementieren
{ var uintV x32 = FN_to_V(x); // x als intDsize-Bit-Zahl
if (FN_V_minusp(x,(sintV)x32)) { x32 = ~ x32; } // falls <0, komplementieren
#if (intVsize>32)
#define x64 x32
logcount_64(); // Bits von x32 zählen
#undef x64
#else
logcount_32(); // Bits von x32 zählen
#endif
return x32;
}
else

8
src/integer/bitwise/cl_I_logtest.cc

@ -33,10 +33,10 @@ cl_boolean logtest (const cl_I& x, const cl_I& y)
}
else
// x Fixnum, y Bignum, also ist x echt kürzer
{ if (FN_L_minusp(x,FN_to_L(x))) return cl_true; // x<0 -> ja.
{ if (FN_V_minusp(x,FN_to_V(x))) return cl_true; // x<0 -> ja.
// x>=0. Kombiniere x mit den pFN_maxlength letzten Digits von y.
{var const uintD* yLSDptr;
var uintL x_ = FN_to_L(x);
var uintV x_ = FN_to_V(x);
BN_to_NDS_nocopy(y, ,,yLSDptr=);
#if (pFN_maxlength > 1)
doconsttimes(pFN_maxlength-1,
@ -50,10 +50,10 @@ cl_boolean logtest (const cl_I& x, const cl_I& y)
else
if (fixnump(y))
// x Bignum, y Fixnum, analog wie oben, nur x und y vertauscht
{ if (FN_L_minusp(y,FN_to_L(y))) return cl_true; // y<0 -> ja.
{ if (FN_V_minusp(y,FN_to_V(y))) return cl_true; // y<0 -> ja.
// y>=0. Kombiniere y mit den pFN_maxlength letzten Digits von x.
{var const uintD* xLSDptr;
var uintL y_ = FN_to_L(y);
var uintV y_ = FN_to_V(y);
BN_to_NDS_nocopy(x, ,,xLSDptr=);
#if (pFN_maxlength > 1)
doconsttimes(pFN_maxlength-1,

150
src/integer/cl_I.h

@ -112,32 +112,32 @@ inline cl_boolean eq (const cl_I& x, sint32 y)
// Umwandlungsroutinen Integer <--> Longword:
// Wandelt Fixnum >=0 in Unsigned Longword um.
// FN_to_UL(obj)
// FN_to_UV(obj)
// > obj: ein Fixnum >=0
// < ergebnis: der Wert des Fixnum als 32-Bit-Zahl.
inline uint32 FN_to_UL (const cl_I& x)
// < ergebnis: der Wert des Fixnum als intVsize-Bit-Zahl.
inline uintV FN_to_UV (const cl_I& x)
{
// This assumes cl_value_shift + cl_value_len == cl_pointer_size.
return (cl_uint)(x.word) >> cl_value_shift;
}
// Wandelt Fixnum in Longword um.
// FN_to_L(obj)
// FN_to_V(obj)
// > obj: ein Fixnum
// < ergebnis: der Wert des Fixnum als 32-Bit-Zahl.
inline sint32 FN_to_L (const cl_I& x)
// < ergebnis: der Wert des Fixnum als intVsize-Bit-Zahl.
inline sintV FN_to_V (const cl_I& x)
{
// This assumes cl_value_shift + cl_value_len == cl_pointer_size.
return (cl_sint)(x.word) >> cl_value_shift;
}
// FN_L_zerop(x,x_) stellt fest, ob x = 0 ist.
// Dabei ist x ein Fixnum und x_ = FN_to_L(x).
#define FN_L_zerop(x,x_) (x_==0)
// FN_V_zerop(x,x_) stellt fest, ob x = 0 ist.
// Dabei ist x ein Fixnum und x_ = FN_to_V(x).
#define FN_V_zerop(x,x_) (x_==0)
// FN_L_minusp(x,x_) stellt fest, ob x < 0 ist.
// Dabei ist x ein Fixnum und x_ = FN_to_L(x).
#define FN_L_minusp(x,x_) (x_<0)
// FN_V_minusp(x,x_) stellt fest, ob x < 0 ist.
// Dabei ist x ein Fixnum und x_ = FN_to_V(x).
#define FN_V_minusp(x,x_) (x_<0)
#ifdef intQsize
@ -212,11 +212,42 @@ inline sint64 FN_to_Q (const cl_I& x)
}
#endif
#ifdef intQsize
// Wandelt Quadword in Integer um.
// Q_to_I(wert)
// > wert: Wert des Integers, ein signed 64-Bit-Integer.
// < ergebnis: Integer mit diesem Wert.
extern cl_private_thing cl_I_constructor_from_Q (sint64 wert);
inline const cl_I Q_to_I (sint64 wert)
{
return cl_I(cl_I_constructor_from_Q(wert));
}
// Wandelt Unsigned Quadword in Integer >=0 um.
// UQ_to_I(wert)
// > wert: Wert des Integers, ein unsigned 64-Bit-Integer.
// < ergebnis: Integer mit diesem Wert.
extern cl_private_thing cl_I_constructor_from_UQ (uint64 wert);
inline const cl_I UQ_to_I (uint64 wert)
{
return cl_I(cl_I_constructor_from_UQ(wert));
}
#endif
// Wandelt Doppel-Longword in Integer um.
// L2_to_I(wert_hi,wert_lo)
// > wert_hi|wert_lo: Wert des Integers, ein signed 64-Bit-Integer.
// < ergebnis: Integer mit diesem Wert.
#if (cl_word_size==64)
inline cl_private_thing cl_I_constructor_from_L2 (sint32 wert_hi, uint32 wert_lo)
{
return cl_I_constructor_from_Q(((sint64)wert_hi<<32) | (sint64)wert_lo);
}
#else
extern cl_private_thing cl_I_constructor_from_L2 (sint32 wert_hi, uint32 wert_lo);
#endif
inline const cl_I L2_to_I (sint32 wert_hi, uint32 wert_lo)
{
return cl_I(cl_I_constructor_from_L2(wert_hi,wert_lo));
@ -226,34 +257,37 @@ inline sint64 FN_to_Q (const cl_I& x)
// UL2_to_I(wert_hi,wert_lo)
// > wert_hi|wert_lo: Wert des Integers, ein unsigned 64-Bit-Integer.
// < ergebnis: Integer mit diesem Wert.
#if (cl_word_size==64)
inline cl_private_thing cl_I_constructor_from_UL2 (uint32 wert_hi, uint32 wert_lo)
{
return cl_I_constructor_from_UQ(((uint64)wert_hi<<32) | (uint64)wert_lo);
}
#else
extern cl_private_thing cl_I_constructor_from_UL2 (uint32 wert_hi, uint32 wert_lo);
#endif
inline const cl_I UL2_to_I (uint32 wert_hi, uint32 wert_lo)
{
return cl_I(cl_I_constructor_from_UL2(wert_hi,wert_lo));
}
#ifdef intQsize
// Wandelt Quadword in Integer um.
// Q_to_I(wert)
// > wert: Wert des Integers, ein signed 64-Bit-Integer.
// Wandelt sintV in Integer um.
// V_to_I(wert)
// > wert: Wert des Integers, ein sintV.
// < ergebnis: Integer mit diesem Wert.
extern cl_private_thing cl_I_constructor_from_Q (sint64 wert);
inline const cl_I Q_to_I (sint64 wert)
{
return cl_I(cl_I_constructor_from_Q(wert));
}
#if (intVsize<=32)
#define V_to_I(wert) L_to_I(wert)
#else
#define V_to_I(wert) Q_to_I(wert)
#endif
// Wandelt Unsigned Quadword in Integer >=0 um.
// UQ_to_I(wert)
// > wert: Wert des Integers, ein unsigned 64-Bit-Integer.
// Wandelt uintV in Integer >=0 um.
// UV_to_I(wert)
// > wert: Wert des Integers, ein uintV.
// < ergebnis: Integer mit diesem Wert.
extern cl_private_thing cl_I_constructor_from_UQ (uint64 wert);
inline const cl_I UQ_to_I (uint64 wert)
{
return cl_I(cl_I_constructor_from_UQ(wert));
}
#if (intVsize<=32)
#define UV_to_I(wert) UL_to_I(wert)
#else
#define UV_to_I(wert) UQ_to_I(wert)
#endif
// Wandelt uintD in Integer >=0 um.
@ -284,22 +318,30 @@ inline const cl_I minus (uintL x, uintL y)
#if (intDsize<=32)
// Holt die nächsten pFN_maxlength Digits in ein uint32.
inline uint32 pFN_maxlength_digits_at (const uintD* ptr)
// Holt die nächsten pFN_maxlength Digits in ein uintV.
inline uintV pFN_maxlength_digits_at (const uintD* ptr)
{
#if (pFN_maxlength==1)
return (uint32)lspref(ptr,0);
return (uintV)lspref(ptr,0);
#elif (pFN_maxlength==2)
return ((uint32)lspref(ptr,1)<<intDsize) | (uint32)lspref(ptr,0);
return ((uintV)lspref(ptr,1)<<intDsize) | (uintV)lspref(ptr,0);
#elif (pFN_maxlength==3)
return ((((uint32)lspref(ptr,2)<<intDsize) | (uint32)lspref(ptr,1))<<intDsize) | (uint32)lspref(ptr,0);
return ((((uintV)lspref(ptr,2)<<intDsize) | (uintV)lspref(ptr,1))<<intDsize) | (uintV)lspref(ptr,0);
#elif (pFN_maxlength==4)
return ((((((uint32)lspref(ptr,3)<<intDsize) | (uint32)lspref(ptr,2))<<intDsize) | (uint32)lspref(ptr,1))<<intDsize) | (uint32)lspref(ptr,0);
return ((((((uintV)lspref(ptr,3)<<intDsize) | (uintV)lspref(ptr,2))<<intDsize) | (uintV)lspref(ptr,1))<<intDsize) | (uintV)lspref(ptr,0);
#elif (pFN_maxlength==5)
return ((((((((uintV)lspref(ptr,4)<<intDsize) | (uintV)lspref(ptr,3))<<intDsize) | (uintV)lspref(ptr,2))<<intDsize) | (uintV)lspref(ptr,1))<<intDsize) | (uintV)lspref(ptr,0);
#elif (pFN_maxlength==6)
return ((((((((((uintV)lspref(ptr,5)<<intDsize) | (uintV)lspref(ptr,4))<<intDsize) | (uintV)lspref(ptr,3))<<intDsize) | (uintV)lspref(ptr,2))<<intDsize) | (uintV)lspref(ptr,1))<<intDsize) | (uintV)lspref(ptr,0);
#elif (pFN_maxlength==7)
return ((((((((((((uintV)lspref(ptr,6)<<intDsize) | (uintV)lspref(ptr,5))<<intDsize) | (uintV)lspref(ptr,4))<<intDsize) | (uintV)lspref(ptr,3))<<intDsize) | (uintV)lspref(ptr,2))<<intDsize) | (uintV)lspref(ptr,1))<<intDsize) | (uintV)lspref(ptr,0);
#elif (pFN_maxlength==8)
return ((((((((((((((uintV)lspref(ptr,7)<<intDsize) | (uintV)lspref(ptr,6))<<intDsize) | (uintV)lspref(ptr,5))<<intDsize) | (uintV)lspref(ptr,4))<<intDsize) | (uintV)lspref(ptr,3))<<intDsize) | (uintV)lspref(ptr,2))<<intDsize) | (uintV)lspref(ptr,1))<<intDsize) | (uintV)lspref(ptr,0);
#endif
}
// Schreibt ein uint32 in die nächsten pFN_maxlength Digits.
inline void set_pFN_maxlength_digits_at (uintD* ptr, uint32 wert)
// Schreibt ein uintV in die nächsten pFN_maxlength Digits.
inline void set_pFN_maxlength_digits_at (uintD* ptr, uintV wert)
{
#if (pFN_maxlength==1)
lspref(ptr,0) = (uintD)wert;
@ -315,6 +357,36 @@ inline void set_pFN_maxlength_digits_at (uintD* ptr, uint32 wert)
lspref(ptr,2) = (uintD)(wert>>(2*intDsize));
lspref(ptr,1) = (uintD)(wert>>intDsize);
lspref(ptr,0) = (uintD)(wert);
#elif (pFN_maxlength==5)
lspref(ptr,4) = (uintD)(wert>>(4*intDsize));
lspref(ptr,3) = (uintD)(wert>>(3*intDsize));
lspref(ptr,2) = (uintD)(wert>>(2*intDsize));
lspref(ptr,1) = (uintD)(wert>>intDsize);
lspref(ptr,0) = (uintD)(wert);
#elif (pFN_maxlength==6)
lspref(ptr,5) = (uintD)(wert>>(5*intDsize));
lspref(ptr,4) = (uintD)(wert>>(4*intDsize));
lspref(ptr,3) = (uintD)(wert>>(3*intDsize));
lspref(ptr,2) = (uintD)(wert>>(2*intDsize));
lspref(ptr,1) = (uintD)(wert>>intDsize);
lspref(ptr,0) = (uintD)(wert);
#elif (pFN_maxlength==7)
lspref(ptr,6) = (uintD)(wert>>(6*intDsize));
lspref(ptr,5) = (uintD)(wert>>(5*intDsize));
lspref(ptr,4) = (uintD)(wert>>(4*intDsize));
lspref(ptr,3) = (uintD)(wert>>(3*intDsize));
lspref(ptr,2) = (uintD)(wert>>(2*intDsize));
lspref(ptr,1) = (uintD)(wert>>intDsize);
lspref(ptr,0) = (uintD)(wert);
#elif (pFN_maxlength==8)
lspref(ptr,7) = (uintD)(wert>>(7*intDsize));
lspref(ptr,6) = (uintD)(wert>>(6*intDsize));
lspref(ptr,5) = (uintD)(wert>>(5*intDsize));
lspref(ptr,4) = (uintD)(wert>>(4*intDsize));
lspref(ptr,3) = (uintD)(wert>>(3*intDsize));
lspref(ptr,2) = (uintD)(wert>>(2*intDsize));
lspref(ptr,1) = (uintD)(wert>>intDsize);
lspref(ptr,0) = (uintD)(wert);
#endif
}

5
src/integer/conv/cl_I_from_L2.cc

@ -10,6 +10,9 @@
// Implementation.
#include "cln/number.h"
#if (cl_word_size < 64)
#include "cl_DS.h"
namespace cln {
@ -133,3 +136,5 @@ cl_private_thing cl_I_constructor_from_L2 (sint32 wert_hi, uint32 wert_lo)
}
} // namespace cln
#endif

5
src/integer/conv/cl_I_from_NDS.cc

@ -47,6 +47,7 @@ const cl_I NDS_to_I (const uintD* MSDptr, uintC len)
// 5 Digits
len_5:
{ wert = get_sint4D_Dptr(MSDptr mspop 5); }
return L_to_FN(wert);
#else // (cl_word_size==64)
var sint64 wert;
#if (intDsize==32)
@ -60,13 +61,13 @@ const cl_I NDS_to_I (const uintD* MSDptr, uintC len)
{ wert = ((sint64)(sintD)mspref(MSDptr,0) << intDsize) | (uint64)(uintD)mspref(MSDptr,1); }
#endif
#if (intDsize==64)
if (TRUE)
if (FALSE)
// 1 Digit
len_1:
{ wert = (sintD)mspref(MSDptr,0); }
#endif
return cl_I_from_word(cl_combine(cl_FN_tag,wert));
#endif
return L_to_FN(wert);
}
#if (cl_value_len > (bn_minlength-1)*intDsize)
if (len == bn_minlength)

2
src/integer/conv/cl_I_from_Q.cc

@ -22,7 +22,7 @@ cl_private_thing cl_I_constructor_from_Q (sint64 wert)
var uint64 test = wert & (sint64)minus_bit(cl_value_len-1);
// test enthält die Bits, die nicht in den Fixnum-Wert >= 0 reinpassen.
if ((test == 0) || (test == (uint64)(sint64)minus_bit(cl_value_len-1)))
return (cl_private_thing)(cl_combine(cl_FN_tag,(sint32)wert));
return (cl_private_thing)(cl_combine(cl_FN_tag,wert));
// Bignum erzeugen:
// (dessen Länge bn_minlength <= n <= ceiling(32/intDsize) erfüllt)
if (wert >= 0) {

5
src/integer/conv/cl_I_from_UL2.cc

@ -10,6 +10,9 @@
// Implementation.
#include "cln/number.h"
#if (cl_word_size < 64)
#include "cl_DS.h"
namespace cln {
@ -104,3 +107,5 @@ cl_private_thing cl_I_constructor_from_UL2 (uint32 wert_hi, uint32 wert_lo)
}
} // namespace cln
#endif

2
src/integer/conv/cl_I_from_UQ.cc

@ -21,7 +21,7 @@ cl_private_thing cl_I_constructor_from_UQ (uint64 wert)
{
if ((wert & (sint64)minus_bit(cl_value_len-1)) == 0)
// Bits, die nicht in den Fixnum-Wert >= 0 reinpassen.
return (cl_private_thing)(cl_combine(cl_FN_tag,(uint32)wert));
return (cl_private_thing)(cl_combine(cl_FN_tag,wert));
// Bignum erzeugen:
// (dessen Länge bn_minlength <= n <= ceiling((32+1)/intDsize) erfüllt)
#define UQ_maxlength ceiling(64+1,intDsize)

101
src/integer/conv/cl_I_to_L.cc

@ -20,61 +20,66 @@ namespace cln {
sint32 cl_I_to_L (const cl_I& obj)
{
if (fixnump(obj))
if (fixnump(obj)) {
// Fixnum
return FN_to_L(obj);
{ // Bignum
var cl_heap_bignum* bn = TheBignum(obj);
var uintC len = bn->length;
if ((sintD)mspref(arrayMSDptr(bn->data,len),0) >= 0) {
// Bignum > 0
#define IF_LENGTH(i) \
if (bn_minlength <= i) /* genau i Digits überhaupt möglich? */\
if (len == i) /* genau i Digits? */ \
/* 2^((i-1)*intDsize-1) <= obj < 2^(i*intDsize-1) */ \
if ( (i*intDsize > 32) \
&& ( ((i-1)*intDsize >= 32) \
|| (mspref(arrayMSDptr(bn->data,len),0) >= (uintD)bitc(31-(i-1)*intDsize)) \
) ) \
goto bad; \
else
IF_LENGTH(1)
return get_uint1D_Dptr(arrayLSDptr(bn->data,1));
IF_LENGTH(2)
return get_uint2D_Dptr(arrayLSDptr(bn->data,2));
IF_LENGTH(3)
return get_uint3D_Dptr(arrayLSDptr(bn->data,3));
IF_LENGTH(4)
return get_uint4D_Dptr(arrayLSDptr(bn->data,4));
#undef IF_LENGTH
} else {
// Bignum < 0
#define IF_LENGTH(i) \
if (bn_minlength <= i) /* genau i Digits überhaupt möglich? */\
if (len == i) /* genau i Digits? */ \
/* - 2^(i*intDsize-1) <= obj < - 2^((i-1)*intDsize-1) */ \
if ( (i*intDsize > 32) \
&& ( ((i-1)*intDsize >= 32) \
|| (mspref(arrayMSDptr(bn->data,len),0) < (uintD)(-bitc(31-(i-1)*intDsize))) \
) ) \
goto bad; \
else
IF_LENGTH(1)
return get_sint1D_Dptr(arrayLSDptr(bn->data,1));
IF_LENGTH(2)
return get_sint2D_Dptr(arrayLSDptr(bn->data,2));
IF_LENGTH(3)
return get_sint3D_Dptr(arrayLSDptr(bn->data,3));
IF_LENGTH(4)
return get_sint4D_Dptr(arrayLSDptr(bn->data,4));
#undef IF_LENGTH
var sintV wert = FN_to_V(obj);
#if (intVsize>32)
if ((sintV)(sint32)wert != wert)
goto bad;
#endif
return (sint32)wert;
} else { // Bignum
var cl_heap_bignum* bn = TheBignum(obj);
var uintC len = bn->length;
if ((sintD)mspref(arrayMSDptr(bn->data,len),0) >= 0) {
// Bignum > 0
#define IF_LENGTH(i) \
if (bn_minlength <= i) /* genau i Digits überhaupt möglich? */\
if (len == i) /* genau i Digits? */ \
/* 2^((i-1)*intDsize-1) <= obj < 2^(i*intDsize-1) */ \
if ( (i*intDsize > 32) \
&& ( ((i-1)*intDsize >= 32) \
|| (mspref(arrayMSDptr(bn->data,len),0) >= (uintD)bitc(31-(i-1)*intDsize)) \
) ) \
goto bad; \
else
IF_LENGTH(1)
return get_uint1D_Dptr(arrayLSDptr(bn->data,1));
IF_LENGTH(2)
return get_uint2D_Dptr(arrayLSDptr(bn->data,2));
IF_LENGTH(3)
return get_uint3D_Dptr(arrayLSDptr(bn->data,3));
IF_LENGTH(4)
return get_uint4D_Dptr(arrayLSDptr(bn->data,4));
#undef IF_LENGTH
} else {
// Bignum < 0
#define IF_LENGTH(i) \
if (bn_minlength <= i) /* genau i Digits überhaupt möglich? */\
if (len == i) /* genau i Digits? */ \
/* - 2^(i*intDsize-1) <= obj < - 2^((i-1)*intDsize-1) */ \
if ( (i*intDsize > 32) \
&& ( ((i-1)*intDsize >= 32) \
|| (mspref(arrayMSDptr(bn->data,len),0) < (uintD)(-bitc(31-(i-1)*intDsize))) \
) ) \
goto bad; \
else
IF_LENGTH(1)
return get_sint1D_Dptr(arrayLSDptr(bn->data,1));
IF_LENGTH(2)
return get_sint2D_Dptr(arrayLSDptr(bn->data,2));
IF_LENGTH(3)
return get_sint3D_Dptr(arrayLSDptr(bn->data,3));
IF_LENGTH(4)
return get_sint4D_Dptr(arrayLSDptr(bn->data,4));
#undef IF_LENGTH
}
}
bad: // unpassendes Objekt
fprint(std::cerr, "Not a 32-bit integer: ");
fprint(std::cerr, obj);
fprint(std::cerr, "\n");
cl_abort();
}
}
} // namespace cln

2
src/integer/conv/cl_I_to_Q.cc

@ -23,7 +23,7 @@ sint64 cl_I_to_Q (const cl_I& obj)
{
if (fixnump(obj))
// Fixnum
return (sint64)(sint32)FN_to_L(obj);
return (sint64)(sintV)FN_to_V(obj);
{ // Bignum
var cl_heap_bignum* bn = TheBignum(obj);
var uintC len = bn->length;

7
src/integer/conv/cl_I_to_UL.cc

@ -22,9 +22,12 @@ uint32 cl_I_to_UL (const cl_I& obj)
{
if (fixnump(obj)) {
// Fixnum
var sint32 wert = FN_to_L(obj);
var sintV wert = FN_to_V(obj);
if (wert >= 0)
return (uint32)wert;
#if (intVsize>32)
if (wert < bit(32))
#endif
return (uint32)wert;
goto bad;
} else { // Bignum
var cl_heap_bignum* bn = TheBignum(obj);

4
src/integer/conv/cl_I_to_UQ.cc

@ -23,9 +23,9 @@ uint64 cl_I_to_UQ (const cl_I& obj)
{
if (fixnump(obj)) {
// Fixnum
var sint32 wert = FN_to_L(obj);
var sintV wert = FN_to_V(obj);
if (wert >= 0)
return (uint64)(uint32)wert;
return (uint64)(uintV)wert;
goto bad;
} else { // Bignum
var cl_heap_bignum* bn = TheBignum(obj);

2
src/integer/conv/cl_I_to_digits.cc

@ -54,7 +54,7 @@ static inline void I_to_digits_noshrink (const cl_I& X, uintD base, uintL erg_le
if (count > 0)
{ var uintB* ptr = erg->MSBptr;
do { *--ptr = '0'; } while (--count > 0);
erg->MSBptr = ptr; erg->len = erg->len;
erg->MSBptr = ptr; erg->len = erg_len;
}
}

69
src/integer/elem/cl_I_div.cc

@ -23,31 +23,60 @@ namespace cln {
// x Fixnum >=0
{ if (fixnump(y))
// auch y Fixnum >=0
{ var uint32 x_ = FN_to_UL(x);
var uint32 y_ = FN_to_UL(y);
{ var uintV x_ = FN_to_UV(x);
var uintV y_ = FN_to_UV(y);
if (y_==0) { cl_error_division_by_0(); }
elif (x_ < y_)
// Trivialfall: q=0, r=x
goto trivial;
elif (y_ < bit(16))
// 32-durch-16-Bit-Division
{ var uint32 q;
var uint16 r;
divu_3216_3216(x_,y_,q=,r=);
return cl_I_div_t(
/* result.quotient = */ UL_to_I(q),
/* result.remainder = */ L_to_FN((uintL)r)
);
}
else
// volle 32-durch-32-Bit-Division
{ var uint32 q;
var uint32 r;
divu_3232_3232(x_,y_,q=,r=);
return cl_I_div_t(
/* result.quotient = */ UL_to_I(q),
/* result.remainder = */ UL_to_I(r)
);
{
#if (intVsize>32)
if (x_ >= bit(32))
{ if (y_ < bit(32))
// 64-durch-32-Bit-Division
{ var uint64 q;
var uint32 r;
divu_6432_6432(x_,y_,q=,r=);
return cl_I_div_t(
/* result.quotient = */ UQ_to_I(q),
/* result.remainder = */ UL_to_I(r)
);
}
else
// volle 64-durch-64-Bit-Division
{ var uint64 q;
var uint64 r;
divu_6464_6464(x_,y_,q=,r=);
return cl_I_div_t(
/* result.quotient = */ UQ_to_I(q),
/* result.remainder = */ UQ_to_I(r)
);
}
}
else
#endif
{ if (y_ < bit(16))
// 32-durch-16-Bit-Division
{ var uint32 q;
var uint16 r;
divu_3216_3216(x_,y_,q=,r=);
return cl_I_div_t(
/* result.quotient = */ UL_to_I(q),
/* result.remainder = */ L_to_FN((uintL)r)
);
}
else
// volle 32-durch-32-Bit-Division
{ var uint32 q;
var uint32 r;
divu_3232_3232(x_,y_,q=,r=);
return cl_I_div_t(
/* result.quotient = */ UL_to_I(q),
/* result.remainder = */ UL_to_I(r)
);
}
}
}
}
else

42
src/integer/elem/cl_I_minus.cc

@ -31,15 +31,15 @@ const cl_I operator- (const cl_I& x, const cl_I& y)
{ // x ist Fixnum
if (fixnump(y))
{ // x,y sind Fixnums
#if (cl_value_len < intLsize)
return L_to_I( FN_to_L(x) - FN_to_L(y) ); // als 32-Bit-Zahlen subtrahieren
#if (cl_value_len < intVsize)
return V_to_I( FN_to_V(x) - FN_to_V(y) ); // als intVsize-Bit-Zahlen subtrahieren
#elif (cl_word_size==64)
return Q_to_I( FN_to_Q(x) - FN_to_Q(y) ); // als 64-Bit-Zahlen subtrahieren
#else // (cl_value_len == intLsize)
var sint32 xhi = sign_of(FN_to_L(x));
var uint32 xlo = FN_to_L(x);
var sint32 yhi = sign_of(FN_to_L(y));
var uint32 ylo = FN_to_L(y);
#elif (intVsize==32) // && (cl_value_len == intVsize)
var sint32 xhi = sign_of(FN_to_V(x));
var uint32 xlo = FN_to_V(x);
var sint32 yhi = sign_of(FN_to_V(y));
var uint32 ylo = FN_to_V(y);
xhi -= yhi;
if (xlo < ylo) { xhi -= 1; }
xlo -= ylo;
@ -49,11 +49,11 @@ const cl_I operator- (const cl_I& x, const cl_I& y)
else
{ // x ist Fixnum, y ist Bignum, also y länger
#if (intDsize==64)
var sint64 x_ = FN_to_L(x); // Wert von x
var sint64 x_ = FN_to_V(x); // Wert von x
#else
var sint32 x_ = FN_to_L(x); // Wert von x
var sintV x_ = FN_to_V(x); // Wert von x
#endif
if (FN_L_zerop(x,x_)) { return -y; } // bei x=0 Ergebnis (- y)
if (FN_V_zerop(x,x_)) { return -y; } // bei x=0 Ergebnis (- y)
CL_ALLOCA_STACK;
BN_to_NDS_1(y, MSDptr=,len=,LSDptr=); // NDS zu y bilden.
// vorsorglich 1 Digit mehr belegen:
@ -70,21 +70,21 @@ const cl_I operator- (const cl_I& x, const cl_I& y)
var uint64 y_new = y_+(uint64)x_;
lspref(LSDptr,0) = y_new;
#else
var uint32 y_ = pFN_maxlength_digits_at(LSDptr);
var uint32 y_new = y_+(uint32)x_;
var uintV y_ = pFN_maxlength_digits_at(LSDptr);
var uintV y_new = y_+(uintV)x_;
set_pFN_maxlength_digits_at(LSDptr,y_new);
#endif
var uintD* midptr = LSDptr lspop pFN_maxlength;
if (y_new < y_)
{ // Carry.
if (!FN_L_minusp(x,x_)) // kürzerer Summand war positiv
if (!FN_V_minusp(x,x_)) // kürzerer Summand war positiv
// Dann ist ein positiver Übertrag weiterzutragen
// (Beispiel: 0002FFFC + 0007 = 00030003)
{ DS_1_plus(midptr,len-pFN_maxlength); }
}
else
{ // Kein Carry.
if (FN_L_minusp(x,x_)) // kürzerer Summand war negativ
if (FN_V_minusp(x,x_)) // kürzerer Summand war negativ
// Dann ist ein negativer Übertrag weiterzutragen
// (Beispiel: 00020003 + FFF5 = 0001FFF8)
{ DS_minus1_plus(midptr,len-pFN_maxlength); }
@ -97,11 +97,11 @@ const cl_I operator- (const cl_I& x, const cl_I& y)
if (fixnump(y))
{ // x ist Bignum, y ist Fixnum, also x länger
#if (intDsize==64)
var sint64 y_ = FN_to_L(y); // Wert von y
var sint64 y_ = FN_to_V(y); // Wert von y
#else
var sint32 y_ = FN_to_L(y); // Wert von y
var sintV y_ = FN_to_V(y); // Wert von y
#endif
if (FN_L_zerop(y,y_)) { return x; } // bei y=0 Ergebnis x
if (FN_V_zerop(y,y_)) { return x; } // bei y=0 Ergebnis x
CL_ALLOCA_STACK;
BN_to_NDS_1(x, MSDptr=,len=,LSDptr=); // NDS zu x bilden.
// len>=bn_minlength. len>pFN_maxlength erzwingen:
@ -116,21 +116,21 @@ const cl_I operator- (const cl_I& x, const cl_I& y)
var uint64 x_new = x_-(uint64)y_;
lspref(LSDptr,0) = x_new;
#else
var uint32 x_ = pFN_maxlength_digits_at(LSDptr);
var uint32 x_new = x_-(uint32)y_;
var uintV x_ = pFN_maxlength_digits_at(LSDptr);
var uintV x_new = x_-(uintV)y_;
set_pFN_maxlength_digits_at(LSDptr,x_new);
#endif
var uintD* midptr = LSDptr lspop pFN_maxlength;
if (x_new > x_)
{ // Carry.
if (!FN_L_minusp(y,y_)) // kürzerer Summand war positiv
if (!FN_V_minusp(y,y_)) // kürzerer Summand war positiv
// Dann ist ein negativer Übertrag weiterzutragen
// (Beispiel: 00030003 - 0007 = 0002FFFC)
{ DS_minus1_plus(midptr,len-pFN_maxlength); }
}
else
{ // Kein Carry.
if (FN_L_minusp(y,y_)) // kürzerer Summand war negativ
if (FN_V_minusp(y,y_)) // kürzerer Summand war negativ
// Dann ist ein positiver Übertrag weiterzutragen
// (Beispiel: 0002FFF8 - FFF5 = 00030003)
{ DS_1_plus(midptr,len-pFN_maxlength); }

25
src/integer/elem/cl_I_mul.cc

@ -26,15 +26,22 @@ const cl_I operator* (const cl_I& x, const cl_I& y)
if (zerop(y))
{ return 0; }
if (fixnump(x) && fixnump(y))
{ var sint32 x_ = FN_to_L(x);
var sint32 y_ = FN_to_L(y);
// Werte direkt multiplizieren:
var uint32 hi;
var uint32 lo;
mulu32((uint32)x_,(uint32)y_,hi=,lo=); // erst unsigned multiplizieren
if (x_ < 0) { hi -= (uint32)y_; } // dann Korrektur für Vorzeichen
if (y_ < 0) { hi -= (uint32)x_; } // (vgl. DS_DS_mul_DS)
return L2_to_I(hi,lo);
{ var sintV x_ = FN_to_V(x);
var sintV y_ = FN_to_V(y);
#if (cl_value_len > 32)
// nur falls x und y Integers mit höchstens 32 Bit sind:
if (((uintV)((sintV)sign_of(x_) ^ x_) < bit(31))
&& ((uintV)((sintV)sign_of(y_) ^ y_) < bit(31)))
#endif
{
// Werte direkt multiplizieren:
var uint32 hi;
var uint32 lo;
mulu32((uint32)x_,(uint32)y_,hi=,lo=); // erst unsigned multiplizieren
if (x_ < 0) { hi -= (uint32)y_; } // dann Korrektur für Vorzeichen
if (y_ < 0) { hi -= (uint32)x_; } // (vgl. DS_DS_mul_DS)
return L2_to_I(hi,lo);
}
}
CL_ALLOCA_STACK;
var const uintD* xMSDptr;

42
src/integer/elem/cl_I_plus.cc

@ -31,15 +31,15 @@ const cl_I operator+ (const cl_I& x, const cl_I& y)
{ // x ist Fixnum
if (fixnump(y))
{ // x,y sind Fixnums
#if (cl_value_len < intLsize)
return L_to_I( FN_to_L(x) + FN_to_L(y) ); // als 32-Bit-Zahlen addieren
#if (cl_value_len < intVsize)
return V_to_I( FN_to_V(x) + FN_to_V(y) ); // als intVsize-Bit-Zahlen addieren
#elif (cl_word_size==64)
return Q_to_I( FN_to_Q(x) + FN_to_Q(y) ); // als 64-Bit-Zahlen addieren
#else // (cl_value_len == intLsize)
var sint32 xhi = sign_of(FN_to_L(x));
var uint32 xlo = FN_to_L(x);
var sint32 yhi = sign_of(FN_to_L(y));
var uint32 ylo = FN_to_L(y);
#elif (intVsize==32) // && (cl_value_len == intVsize)
var sint32 xhi = sign_of(FN_to_V(x));
var uint32 xlo = FN_to_V(x);
var sint32 yhi = sign_of(FN_to_V(y));
var uint32 ylo = FN_to_V(y);
xhi += yhi;
xlo += ylo;
if (xlo < ylo) { xhi += 1; }
@ -49,11 +49,11 @@ const cl_I operator+ (const cl_I& x, const cl_I& y)
else
{ // x ist Fixnum, y ist Bignum, also y länger
#if (intDsize==64)
var sint64 x_ = FN_to_L(x); // Wert von x
var sint64 x_ = FN_to_V(x); // Wert von x
#else
var sint32 x_ = FN_to_L(x); // Wert von x
var sintV x_ = FN_to_V(x); // Wert von x
#endif
if (FN_L_zerop(x,x_)) { return y; } // bei x=0 Ergebnis y
if (FN_V_zerop(x,x_)) { return y; } // bei x=0 Ergebnis y
CL_ALLOCA_STACK;
BN_to_NDS_1(y, MSDptr=,len=,LSDptr=); // NDS zu y bilden.
// len>=bn_minlength. len>pFN_maxlength erzwingen:
@ -68,21 +68,21 @@ const cl_I operator+ (const cl_I& x, const cl_I& y)
var uint64 y_new = y_+(uint64)x_;
lspref(LSDptr,0) = y_new;
#else
var uint32 y_ = pFN_maxlength_digits_at(LSDptr);
var uint32 y_new = y_+(uint32)x_;
var uintV y_ = pFN_maxlength_digits_at(LSDptr);
var uintV y_new = y_+(uintV)x_;
set_pFN_maxlength_digits_at(LSDptr,y_new);
#endif
var uintD* midptr = LSDptr lspop pFN_maxlength;
if (y_new < y_)
{ // Carry.
if (!FN_L_minusp(x,x_)) // kürzerer Summand war positiv
if (!FN_V_minusp(x,x_)) // kürzerer Summand war positiv
// Dann ist ein positiver Übertrag weiterzutragen
// (Beispiel: 0002FFFC + 0007 = 00030003)
{ DS_1_plus(midptr,len-pFN_maxlength); }
}
else
{ // Kein Carry.
if (FN_L_minusp(x,x_)) // kürzerer Summand war negativ
if (FN_V_minusp(x,x_)) // kürzerer Summand war negativ
// Dann ist ein negativer Übertrag weiterzutragen
// (Beispiel: 00020003 + FFF5 = 0001FFF8)
{ DS_minus1_plus(midptr,len-pFN_maxlength); }
@ -95,11 +95,11 @@ const cl_I operator+ (const cl_I& x, const cl_I& y)
if (fixnump(y))
{ // x ist Bignum, y ist Fixnum, also x länger
#if (intDsize==64)
var sint64 y_ = FN_to_L(y); // Wert von y
var sint64 y_ = FN_to_V(y); // Wert von y
#else
var sint32 y_ = FN_to_L(y); // Wert von y
var sintV y_ = FN_to_V(y); // Wert von y
#endif
if (FN_L_zerop(y,y_)) { return x; } // bei y=0 Ergebnis x
if (FN_V_zerop(y,y_)) { return x; } // bei y=0 Ergebnis x
CL_ALLOCA_STACK;
BN_to_NDS_1(x, MSDptr=,len=,LSDptr=); // NDS zu x bilden.
// len>=bn_minlength. len>pFN_maxlength erzwingen:
@ -114,21 +114,21 @@ const cl_I operator+ (const cl_I& x, const cl_I& y)
var uint64 x_new = x_+(uint64)y_;
lspref(LSDptr,0) = x_new;
#else
var uint32 x_ = pFN_maxlength_digits_at(LSDptr);
var uint32 x_new = x_+(uint32)y_;
var uintV x_ = pFN_maxlength_digits_at(LSDptr);
var uintV x_new = x_+(uintV)y_;
set_pFN_maxlength_digits_at(LSDptr,x_new);
#endif
var uintD* midptr = LSDptr lspop pFN_maxlength;
if (x_new < x_)
{ // Carry.
if (!FN_L_minusp(y,y_)) // kürzerer Summand war positiv
if (!FN_V_minusp(y,y_)) // kürzerer Summand war positiv
// Dann ist ein positiver Übertrag weiterzutragen
// (Beispiel: 0002FFFC + 0007 = 00030003)
{ DS_1_plus(midptr,len-pFN_maxlength); }
}
else
{ // Kein Carry.
if (FN_L_minusp(y,y_)) // kürzerer Summand war negativ
if (FN_V_minusp(y,y_)) // kürzerer Summand war negativ
// Dann ist ein negativer Übertrag weiterzutragen
// (Beispiel: 00020003 + FFF5 = 0001FFF8)
{ DS_minus1_plus(midptr,len-pFN_maxlength); }

20
src/integer/elem/cl_I_square.cc

@ -21,13 +21,19 @@ const cl_I square (const cl_I& x)
// x Fixnum -> direkt multiplizieren
// sonst: zu DS machen, multiplizieren.
if (fixnump(x))
{ var sint32 x_ = FN_to_L(x);
// Werte direkt multiplizieren:
var uint32 hi;
var uint32 lo;
mulu32((uint32)x_,(uint32)x_,hi=,lo=); // erst unsigned multiplizieren
if (x_ < 0) { hi -= 2*(uint32)x_; } // dann Korrektur für Vorzeichen
return L2_to_I(hi,lo);
{ var sintV x_ = FN_to_V(x);
#if (cl_value_len > 32)
// nur falls x ein Integer mit höchstens 32 Bit ist:
if ((uintV)((sintV)sign_of(x_) ^ x_) < bit(31))
#endif
{
// Werte direkt multiplizieren:
var uint32 hi;
var uint32 lo;
mulu32((uint32)x_,(uint32)x_,hi=,lo=); // erst unsigned multiplizieren
if (x_ < 0) { hi -= 2*(uint32)x_; } // dann Korrektur für Vorzeichen
return L2_to_I(hi,lo);
}
}
CL_ALLOCA_STACK;
var const uintD* xMSDptr;

10
src/integer/elem/cl_I_uminus.cc

@ -17,13 +17,13 @@ namespace cln {
const cl_I operator- (const cl_I& x)
{
if (fixnump(x)) {
#if (cl_value_len < intLsize)
return L_to_I(- FN_to_L(x));
#if (cl_value_len < intVsize)
return V_to_I(- FN_to_V(x));
#elif (cl_word_size==64)
return Q_to_I(- FN_to_Q(x));
#else // (cl_value_len == intLsize)
var sint32 xhi = sign_of(FN_to_L(x));
var uint32 xlo = FN_to_L(x);
#elif (intVsize==32) // && (cl_value_len == intVsize)
var sint32 xhi = sign_of(FN_to_V(x));
var uint32 xlo = FN_to_V(x);
return L2_to_I(-xhi-(xlo!=0),-xlo);
#endif
} else {

8
src/integer/gcd/cl_I_gcd.cc

@ -154,7 +154,7 @@ namespace cln {
var cl_I& b = abs_b;
if (fixnump(a) && fixnump(b)) // ggT zweier Fixnums >0
{ // bleibt Fixnum, da (gcd a b) <= (min a b)
return L_to_FN(gcd(FN_to_UL(a),FN_to_UL(b)));
return V_to_FN(gcd(FN_to_UV(a),FN_to_UV(b)));
}
{ var cl_signean vergleich = compare(a,b);
if (vergleich == 0) { return a; } // a=b -> fertig
@ -296,11 +296,11 @@ namespace cln {
if (eq(b,0)) { return abs(a); } // b=0 -> (abs a)
if (eq(a,0)) { return abs(b); } // a=0 -> (abs b)
if (fixnump(a) && fixnump(b)) // ggT zweier Fixnums /=0
{ var sintL a_ = FN_to_L(a);
{ var sintV a_ = FN_to_V(a);
if (a_ < 0) { a_ = -a_; }
var sintL b_ = FN_to_L(b);
var sintV b_ = FN_to_V(b);
if (b_ < 0) { b_ = -b_; }
return UL_to_I(gcd((uint32)a_,(uint32)b_));
return UV_to_I(gcd((uintV)a_,(uintV)b_));
}
CL_ALLOCA_STACK;
var uintD* a_MSDptr;

6
src/integer/gcd/cl_low_gcd.cc

@ -15,7 +15,7 @@ namespace cln {
// gcd(a,b)
// > a,b: zwei Integers
// < ergebnis: (gcd a b), ein Integer >=0
uint32 gcd (uint32 a, uint32 b)
uintV gcd (uintV a, uintV b)
// binäre Methode:
// (gcd a b) :==
// (prog ((j 0))
@ -44,7 +44,7 @@ namespace cln {
#ifdef DUMMER_GGT // so macht's ein Mathematiker:
if (a==0) { return b; }
if (b==0) { return a; }
var uint32 bit_j = bit(0);
var uintV bit_j = bit(0);
loop
{ // a,b >0
if (!((a & bit_j) ==0))
@ -54,7 +54,7 @@ namespace cln {
bit_j = bit_j<<1;
}
#else // Trick von B. Degel:
var uint32 bit_j = (a | b); // endet mit einer 1 und j Nullen
var uintV bit_j = (a | b); // endet mit einer 1 und j Nullen
bit_j = bit_j ^ (bit_j - 1); // Maske = bit(j) | bit(j-1) | ... | bit(0)
if (!((a & bit_j) ==0))
{ if (!((b & bit_j) ==0))

2
src/integer/hash/cl_I_hashcode.cc

@ -18,7 +18,7 @@ unsigned long hashcode (const cl_I& x)
// integers, but it's better than completely ignoring some limbs.
if (fixnump(x)) {
#if (cl_value_len <= intLsize)
code += FN_to_L(x);
code += FN_to_V(x);
#elif (cl_word_size==64)
code += FN_to_Q(x);
code ^= (code >> 32);

2
src/integer/input/cl_I_read.cc

@ -86,7 +86,7 @@ const cl_I read_integer (const cl_read_flags& flags, const char * string, const
fprint(std::cerr, "\n");
cl_abort();
}
rational_base = FN_to_UL(base); ptr = base_end_ptr;
rational_base = FN_to_UV(base); ptr = base_end_ptr;
break;
}
ptr++;

18
src/integer/misc/cl_I_eqhashcode.cc

@ -17,19 +17,25 @@ namespace cln {
inline uint32 equal_hashcode (const cl_FN& x)
{
var cl_signean sign;
var uint32 x32 = FN_to_L(x); // x als 32-Bit-Zahl
if (FN_L_minusp(x,(sint32)x32)) {
x32 = -x32;
var uintV x_ = FN_to_V(x); // x als intVsize-Bit-Zahl
if (FN_V_minusp(x,(sintV)x_)) {
x_ = -x_;
sign = -1;
} else {
sign = 0;
if (x32 == 0)
if (x_ == 0)
return 0;
}
var uintL s;
integerlength32(x32, s = 32 - );
var uint32 msd = x32 << s;
#if (intVsize > 32)
integerlength64(x_, s = 64 - );
var uint32 msd = (x_ << s) >> 32;
var sintL exp = 64-s;
#else
integerlength32(x_, s = 32 - );
var uint32 msd = x_ << s;
var sintL exp = 32-s;
#endif
return equal_hashcode_low(msd,exp,sign);
}

9
src/integer/misc/cl_I_ord2.cc

@ -31,8 +31,13 @@ namespace cln {
uintL ord2 (const cl_I& x) // x /= 0
{
if (fixnump(x))
{ var uint32 x_ = FN_to_L(x); // x als 32-Bit-Zahl
ord2_32(x_,return);
{ var uintV x_ = FN_to_V(x); // x als intVsize-Bit-Zahl
// This assumes cl_value_len <= intVsize.
#if (intVsize>32)
ord2_64(x_,return);
#else
ord2_32(x_,return);
#endif
}
else
{ var uintL bitcount = 0;

6
src/integer/misc/cl_I_power2p.cc

@ -21,9 +21,13 @@ uintL power2p (const cl_I& x) // x > 0
// Methode 3: Wenn das erste Digit /=0 eine Zweierpotenz ist und alle weiteren
// Digits Null sind.
if (fixnump(x))
{ var uintL x_ = FN_to_UL(x);
{ var uintV x_ = FN_to_UV(x);
if (!((x_ & (x_-1)) == 0)) return 0; // keine Zweierpotenz
#if (intVsize>32)
integerlength64(x_,return); // Zweierpotenz: n = integer_length(x)
#else
integerlength32(x_,return); // Zweierpotenz: n = integer_length(x)
#endif
}
else
{ var const uintD* MSDptr;

82
src/integer/misc/combin/cl_I_doublefactorial.cc

@ -29,42 +29,70 @@ namespace cln {
// two in factorial(m).
static cl_I const doublefakul_table [] = {
L_to_FN(1),
L_to_FN(1UL),
L_to_FN(1UL*2),
L_to_FN(1UL*3),
1,
1UL,
1UL*2,
1UL*3,
#if (cl_value_len>=5)
L_to_FN(1UL*2*4),
L_to_FN(1UL*3*5),
1UL*2*4,
1UL*3*5,
#if (cl_value_len>=7)
L_to_FN(1UL*2*4*6),
1UL*2*4*6,
#if (cl_value_len>=8)
L_to_FN(1UL*3*5*7),
1UL*3*5*7,
#if (cl_value_len>=10)
L_to_FN(1UL*2*4*6*8),
1UL*2*4*6*8,
#if (cl_value_len>=11)
L_to_FN(1UL*3*5*7*9),
1UL*3*5*7*9,
#if (cl_value_len>=13)
L_to_FN(1UL*2*4*6*8*10),
1UL*2*4*6*8*10,
#if (cl_value_len>=15)
L_to_FN(1UL*3*5*7*9*11),
1UL*3*5*7*9*11,
#if (cl_value_len>=17)
L_to_FN(1UL*2*4*6*8*10*12),
1UL*2*4*6*8*10*12,
#if (cl_value_len>=19)
L_to_FN(1UL*3*5*7*9*11*13),
1UL*3*5*7*9*11*13,
#if (cl_value_len>=21)
L_to_FN(1UL*2*4*6*8*10*12*14),
1UL*2*4*6*8*10*12*14,
#if (cl_value_len>=22)
L_to_FN(1UL*3*5*7*9*11*13*15),
1UL*3*5*7*9*11*13*15,
#if (cl_value_len>=25)
L_to_FN(1UL*2*4*6*8*10*12*14*16),
1UL*2*4*6*8*10*12*14*16,
#if (cl_value_len>=27)
L_to_FN(1UL*3*5*7*9*11*13*15*17),
1UL*3*5*7*9*11*13*15*17,
#if (cl_value_len>=29)
L_to_FN(1UL*2*4*6*8*10*12*14*16*18),
1UL*2*4*6*8*10*12*14*16*18,
#if (cl_value_len>=31)
L_to_FN(1UL*3*5*7*9*11*13*15*17*19),
1UL*3*5*7*9*11*13*15*17*19,
#if (cl_value_len>=33)
1UL*2*4*6*8*10*12*14*16*18*20,
#if (cl_value_len>=35)
1UL*3*5*7*9*11*13*15*17*19*21,
#if (cl_value_len>=38)
1UL*2*4*6*8*10*12*14*16*18*20*22,
#if (cl_value_len>=40)
1UL*3*5*7*9*11*13*15*17*19*21*23,
#if (cl_value_len>=42)
1UL*2*4*6*8*10*12*14*16*18*20*22*24,
#if (cl_value_len>=44)
1UL*3*5*7*9*11*13*15*17*19*21*23*25,
#if (cl_value_len>=47)
1UL*2*4*6*8*10*12*14*16*18*20*22*24*26,
#if (cl_value_len>=49)
1UL*3*5*7*9*11*13*15*17*19*21*23*25*27,
#if (cl_value_len>=52)
1UL*2*4*6*8*10*12*14*16*18*20*22*24*26*28,
#if (cl_value_len>=54)
1UL*3*5*7*9*11*13*15*17*19*21*23*25*27*29,
#if (cl_value_len>=57)
1UL*2*4*6*8*10*12*14*16*18*20*22*24*26*28*30,
#if (cl_value_len>=59)
1UL*3*5*7*9*11*13*15*17*19*21*23*25*27*29*31,
#if (cl_value_len>=62)
1UL*2*4*6*8*10*12*14*16*18*20*22*24*26*28*30*32,
#if (cl_value_len>=64)
1UL*3*5*7*9*11*13*15*17*19*21*23*25*27*29*31*33,
#if (cl_value_len>=67)
...
#endif
#endif
@ -82,6 +110,20 @@ static cl_I const doublefakul_table [] = {
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
};
const cl_I doublefactorial (uintL n) // assume n >= 0 small

54
src/integer/misc/combin/cl_I_factorial.cc

@ -30,31 +30,47 @@ namespace cln {
// vermeidet, daß oft große Zahlen mit ganz kleinen Zahlen multipliziert
// werden.
static cl_I const fakul_table [] = {
L_to_FN(1),
L_to_FN(1UL),
L_to_FN(1UL*2),
static uintV const fakul_table [] = {
1,
1UL,
1UL*2,
#if (cl_value_len>=4)
L_to_FN(1UL*2*3),
1UL*2*3,
#if (cl_value_len>=6)
L_to_FN(1UL*2*3*4),
1UL*2*3*4,
#if (cl_value_len>=8)
L_to_FN(1UL*2*3*4*5),
1UL*2*3*4*5,
#if (cl_value_len>=11)
L_to_FN(1UL*2*3*4*5*6),
1UL*2*3*4*5*6,
#if (cl_value_len>=14)
L_to_FN(1UL*2*3*4*5*6*7),
1UL*2*3*4*5*6*7,
#if (cl_value_len>=17)
L_to_FN(1UL*2*3*4*5*6*7*8),
1UL*2*3*4*5*6*7*8,
#if (cl_value_len>=20)
L_to_FN(1UL*2*3*4*5*6*7*8*9),
1UL*2*3*4*5*6*7*8*9,
#if (cl_value_len>=23)
L_to_FN(1UL*2*3*4*5*6*7*8*9*10),
1UL*2*3*4*5*6*7*8*9*10,
#if (cl_value_len>=27)
L_to_FN(1UL*2*3*4*5*6*7*8*9*10*11),
1UL*2*3*4*5*6*7*8*9*10*11,
#if (cl_value_len>=30)
L_to_FN(1UL*2*3*4*5*6*7*8*9*10*11*12),
1UL*2*3*4*5*6*7*8*9*10*11*12,
#if (cl_value_len>=34)
1UL*2*3*4*5*6*7*8*9*10*11*12*13,
#if (cl_value_len>=38)
1UL*2*3*4*5*6*7*8*9*10*11*12*13*14,
#if (cl_value_len>=42)
1UL*2*3*4*5*6*7*8*9*10*11*12*13*14*15,
#if (cl_value_len>=46)
1UL*2*3*4*5*6*7*8*9*10*11*12*13*14*15*16,
#if (cl_value_len>=50)
1UL*2*3*4*5*6*7*8*9*10*11*12*13*14*15*16*17,
#if (cl_value_len>=54)
1UL*2*3*4*5*6*7*8*9*10*11*12*13*14*15*16*17*18,
#if (cl_value_len>=58)
1UL*2*3*4*5*6*7*8*9*10*11*12*13*14*15*16*17*18*19,
#if (cl_value_len>=63)
1UL*2*3*4*5*6*7*8*9*10*11*12*13*14*15*16*17*18*19*20,
#if (cl_value_len>=67)
...
#endif
#endif
@ -67,12 +83,20 @@ static cl_I const fakul_table [] = {
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
#endif
};
const cl_I factorial (uintL n) // assume n >= 0 small
{
if (n < sizeof(fakul_table)/sizeof(cl_I))
{ return fakul_table[n]; }
{ return UV_to_I(fakul_table[n]); }
else
{ var cl_I prod = 1; // bisheriges Produkt := 1
var uintL k = 1;

24
src/modinteger/cl_MI_fix16.h

@ -4,41 +4,41 @@ namespace cln {
static const _cl_MI fix16_plus (cl_heap_modint_ring* R, const _cl_MI& x, const _cl_MI& y)
{
var uint32 zr = FN_to_UL(x.rep) + FN_to_UL(y.rep);
if (zr >= FN_to_UL(R->modulus)) { zr = zr - FN_to_UL(R->modulus); }
var uint32 zr = FN_to_UV(x.rep) + FN_to_UV(y.rep);
if (zr >= FN_to_UV(R->modulus)) { zr = zr - FN_to_UV(R->modulus); }
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix16_minus (cl_heap_modint_ring* R, const _cl_MI& x, const _cl_MI& y)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 yr = FN_to_UL(y.rep);
var uint32 xr = FN_to_UV(x.rep);
var uint32 yr = FN_to_UV(y.rep);
var sint32 zr = xr - yr;
if (zr < 0) { zr = zr + FN_to_UL(R->modulus); }
if (zr < 0) { zr = zr + FN_to_UV(R->modulus); }
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix16_uminus (cl_heap_modint_ring* R, const _cl_MI& x)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 zr = (xr==0 ? 0 : FN_to_UL(R->modulus)-xr);
var uint32 xr = FN_to_UV(x.rep);
var uint32 zr = (xr==0 ? 0 : FN_to_UV(R->modulus)-xr);
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix16_mul (cl_heap_modint_ring* R, const _cl_MI& x, const _cl_MI& y)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 yr = FN_to_UL(y.rep);
var uint32 xr = FN_to_UV(x.rep);
var uint32 yr = FN_to_UV(y.rep);
var uint32 zr = mulu16(xr,yr);
divu_3216_1616(zr,FN_to_UL(R->modulus),,zr=);
divu_3216_1616(zr,FN_to_UV(R->modulus),,zr=);
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix16_square (cl_heap_modint_ring* R, const _cl_MI& x)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 xr = FN_to_UV(x.rep);
var uint32 zr = mulu16(xr,xr);
divu_3216_1616(zr,FN_to_UL(R->modulus),,zr=);
divu_3216_1616(zr,FN_to_UV(R->modulus),,zr=);
return _cl_MI(R, L_to_FN(zr));
}

24
src/modinteger/cl_MI_fix29.h

@ -5,47 +5,47 @@ namespace cln {
static const _cl_MI fix29_plus (cl_heap_modint_ring* R, const _cl_MI& x, const _cl_MI& y)
{
var uint32 zr = FN_to_UL(x.rep) + FN_to_UL(y.rep);
if (zr >= FN_to_UL(R->modulus)) { zr = zr - FN_to_UL(R->modulus); }
var uint32 zr = FN_to_UV(x.rep) + FN_to_UV(y.rep);
if (zr >= FN_to_UV(R->modulus)) { zr = zr - FN_to_UV(R->modulus); }
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix29_minus (cl_heap_modint_ring* R, const _cl_MI& x, const _cl_MI& y)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 yr = FN_to_UL(y.rep);
var uint32 xr = FN_to_UV(x.rep);
var uint32 yr = FN_to_UV(y.rep);
var sint32 zr = xr - yr;
if (zr < 0) { zr = zr + FN_to_UL(R->modulus); }
if (zr < 0) { zr = zr + FN_to_UV(R->modulus); }
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix29_uminus (cl_heap_modint_ring* R, const _cl_MI& x)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 zr = (xr==0 ? 0 : FN_to_UL(R->modulus)-xr);
var uint32 xr = FN_to_UV(x.rep);
var uint32 zr = (xr==0 ? 0 : FN_to_UV(R->modulus)-xr);
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix29_mul (cl_heap_modint_ring* R, const _cl_MI& x, const _cl_MI& y)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 yr = FN_to_UL(y.rep);
var uint32 xr = FN_to_UV(x.rep);
var uint32 yr = FN_to_UV(y.rep);
var uint32 zrhi;
var uint32 zrlo;
mulu32(xr,yr,zrhi=,zrlo=);
var uint32 zr;
divu_6432_3232(zrhi,zrlo,FN_to_UL(R->modulus),,zr=);
divu_6432_3232(zrhi,zrlo,FN_to_UV(R->modulus),,zr=);
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix29_square (cl_heap_modint_ring* R, const _cl_MI& x)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 xr = FN_to_UV(x.rep);
var uint32 zrhi;
var uint32 zrlo;
mulu32(xr,xr,zrhi=,zrlo=);
var uint32 zr;
divu_6432_3232(zrhi,zrlo,FN_to_UL(R->modulus),,zr=);
divu_6432_3232(zrhi,zrlo,FN_to_UV(R->modulus),,zr=);
return _cl_MI(R, L_to_FN(zr));
}

26
src/modinteger/cl_MI_fix32.h

@ -5,49 +5,49 @@ namespace cln {
static const _cl_MI fix32_plus (cl_heap_modint_ring* R, const _cl_MI& x, const _cl_MI& y)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 yr = FN_to_UL(y.rep);
var uint32 xr = FN_to_UV(x.rep);
var uint32 yr = FN_to_UV(y.rep);
var uint32 zr = xr + yr;
var uint32 m = FN_to_UL(R->modulus);
var uint32 m = FN_to_UV(R->modulus);
if ((zr < xr) || (zr >= m)) { zr = zr - m; }
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix32_minus (cl_heap_modint_ring* R, const _cl_MI& x, const _cl_MI& y)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 yr = FN_to_UL(y.rep);
var sint32 zr = (xr >= yr ? xr - yr : xr - yr + FN_to_UL(R->modulus));
var uint32 xr = FN_to_UV(x.rep);
var uint32 yr = FN_to_UV(y.rep);
var sint32 zr = (xr >= yr ? xr - yr : xr - yr + FN_to_UV(R->modulus));
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix32_uminus (cl_heap_modint_ring* R, const _cl_MI& x)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 zr = (xr==0 ? 0 : FN_to_UL(R->modulus)-xr);
var uint32 xr = FN_to_UV(x.rep);
var uint32 zr = (xr==0 ? 0 : FN_to_UV(R->modulus)-xr);
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix32_mul (cl_heap_modint_ring* R, const _cl_MI& x, const _cl_MI& y)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 yr = FN_to_UL(y.rep);
var uint32 xr = FN_to_UV(x.rep);
var uint32 yr = FN_to_UV(y.rep);
var uint32 zrhi;
var uint32 zrlo;
mulu32(xr,yr,zrhi=,zrlo=);
var uint32 zr;
divu_6432_3232(zrhi,zrlo,FN_to_UL(R->modulus),,zr=);
divu_6432_3232(zrhi,zrlo,FN_to_UV(R->modulus),,zr=);
return _cl_MI(R, L_to_FN(zr));
}
static const _cl_MI fix32_square (cl_heap_modint_ring* R, const _cl_MI& x)
{
var uint32 xr = FN_to_UL(x.rep);
var uint32 xr = FN_to_UV(x.rep);
var uint32 zrhi;
var uint32 zrlo;
mulu32(xr,xr,zrhi=,zrlo=);
var uint32 zr;
divu_6432_3232(zrhi,zrlo,FN_to_UL(R->modulus),,zr=);
divu_6432_3232(zrhi,zrlo,FN_to_UV(R->modulus),,zr=);
return _cl_MI(R, L_to_FN(zr));
}

2
src/modinteger/cl_MI_std.h

@ -175,7 +175,7 @@ static const _cl_MI std_expt_pos (cl_heap_modint_ring* R, const _cl_MI& x, const
// n has nn bits.
if (nn <= 8) {
// k = 1, normal Left-Right Binary algorithm.
var uintL _n = FN_to_UL(n);
var uintL _n = FN_to_UV(n);
var _cl_MI a = x;
for (var int i = nn-2; i >= 0; i--) {
a = R->_square(a);

2
src/numtheory/cl_nt_cornacchia4.cc

@ -57,7 +57,7 @@ const cornacchia_t cornacchia4 (const cl_I& d, const cl_I& p)
if (d==7) return cornacchia_t(1, 1,1);
return cornacchia_t(0);
}
switch (FN_to_L(logand(d,7))) {
switch (FN_to_V(logand(d,7))) {
case 0: case 4: {
// d == 0 mod 4
var cornacchia_t s = cornacchia1(d>>2,p);

8
src/numtheory/cl_nt_jacobi.cc

@ -29,7 +29,7 @@ int jacobi (const cl_I& a, const cl_I& b)
a = mod(a,b);
// If a and b are fixnums, choose faster routine.
if (fixnump(b))
return jacobi(FN_to_L(a),FN_to_L(b));
return jacobi(FN_to_V(a),FN_to_V(b));
var int v = 1;
for (;;) {
// (a/b) * v is invariant.
@ -43,7 +43,7 @@ int jacobi (const cl_I& a, const cl_I& b)
// a > b/2, so (a/b) = (-1/b) * ((b-a)/b),
// and (-1/b) = -1 if b==3 mod 4.
a = b-a;
if (FN_to_L(logand(b,3)) == 3)
if (FN_to_V(logand(b,3)) == 3)
v = -v;
continue;
}
@ -51,14 +51,14 @@ int jacobi (const cl_I& a, const cl_I& b)
// b>1 and a=2a', so (a/b) = (2/b) * (a'/b),
// and (2/b) = -1 if b==3,5 mod 8.
a = a>>1;
switch (FN_to_L(logand(b,7))) {
switch (FN_to_V(logand(b,7))) {
case 3: case 5: v = -v; break;
}
continue;
}
// a and b odd, 0 < a < b/2 < b, so apply quadratic reciprocity
// law (a/b) = (-1)^((a-1)/2)((b-1)/2) * (b/a).
if (FN_to_L(logand(logand(a,b),3)) == 3)
if (FN_to_V(logand(logand(a,b),3)) == 3)
v = -v;
swap(cl_I, a,b);
// Now a > 2*b, set a := a mod b.

10
src/numtheory/cl_nt_jacobi_low.cc

@ -15,7 +15,7 @@
namespace cln {
// Assume 0 <= a < b.
inline int jacobi_aux (uint32 a, uint32 b)
inline int jacobi_aux (uintV a, uintV b)
{
var int v = 1;
for (;;) {
@ -52,7 +52,7 @@ inline int jacobi_aux (uint32 a, uint32 b)
// law (a/b) = (-1)^((a-1)/2)((b-1)/2) * (b/a).
if ((a & b & 3) == 3)
v = -v;
swap(sint32, a,b);
swap(uintV, a,b);
// Now a > 2*b, set a := a mod b.
if ((a >> 3) >= b)
a = a % b;
@ -61,7 +61,7 @@ inline int jacobi_aux (uint32 a, uint32 b)
}
}
int jacobi (sint32 a, sint32 b)
int jacobi (sintV a, sintV b)
{
// Check b > 0, b odd.
if (!(b > 0))
@ -70,9 +70,9 @@ int jacobi (sint32 a, sint32 b)
cl_abort();
// Ensure 0 <= a < b.
if (a >= 0)
a = (uint32)a % (uint32)b;
a = (uintV)a % (uintV)b;
else
a = b-1-((uint32)(~a) % (uint32)b);
a = b-1-((uintV)(~a) % (uintV)b);
return jacobi_aux(a,b);
}

2
src/rational/input/cl_RA_read.cc

@ -88,7 +88,7 @@ const cl_RA read_rational (const cl_read_flags& flags, const char * string, cons
fprint(std::cerr, "\n");
cl_abort();
}
rational_base = FN_to_UL(base); ptr = base_end_ptr;
rational_base = FN_to_UV(base); ptr = base_end_ptr;
break;
}
ptr++;

2
src/real/input/cl_R_read.cc

@ -95,7 +95,7 @@ const cl_R read_real (const cl_read_flags& flags, const char * string, const cha
fprint(std::cerr, "\n");
cl_abort();
}
rational_base = FN_to_UL(base); ptr = base_end_ptr;
rational_base = FN_to_UV(base); ptr = base_end_ptr;
break;
}
ptr++;

20
src/vector/cl_GV_I.cc

@ -283,9 +283,9 @@ static const cl_I bits1_element (const cl_GV_inner<cl_I>* vec, uintL index)
}
static void bits1_set_element (cl_GV_inner<cl_I>* vec, uintL index, const cl_I& x)
{
var uint32 xval;
var uintV xval;
if (fixnump(x)) {
xval = FN_to_UL(x);
xval = FN_to_UV(x);
if (xval <= 0x1) {
var uintD* ptr = &((cl_heap_GV_I_bits1 *) outcast(vec))->data[index/intDsize];
index = index%intDsize;
@ -305,9 +305,9 @@ static const cl_I bits2_element (const cl_GV_inner<cl_I>* vec, uintL index)
}
static void bits2_set_element (cl_GV_inner<cl_I>* vec, uintL index, const cl_I& x)
{
var uint32 xval;
var uintV xval;
if (fixnump(x)) {
xval = FN_to_UL(x);
xval = FN_to_UV(x);
if (xval <= 0x3) {
var uintD* ptr = &((cl_heap_GV_I_bits2 *) outcast(vec))->data[index/(intDsize/2)];
index = 2*(index%(intDsize/2));
@ -327,9 +327,9 @@ static const cl_I bits4_element (const cl_GV_inner<cl_I>* vec, uintL index)
}
static void bits4_set_element (cl_GV_inner<cl_I>* vec, uintL index, const cl_I& x)
{
var uint32 xval;
var uintV xval;
if (fixnump(x)) {
xval = FN_to_UL(x);
xval = FN_to_UV(x);
if (xval <= 0xF) {
var uintD* ptr = &((cl_heap_GV_I_bits4 *) outcast(vec))->data[index/(intDsize/4)];
index = 4*(index%(intDsize/4));
@ -354,9 +354,9 @@ static const cl_I bits8_element (const cl_GV_inner<cl_I>* vec, uintL index)
}
static void bits8_set_element (cl_GV_inner<cl_I>* vec, uintL index, const cl_I& x)
{
var uint32 xval;
var uintV xval;
if (fixnump(x)) {
xval = FN_to_UL(x);
xval = FN_to_UV(x);
if (xval <= 0xFF) {
#if CL_CPU_BIG_ENDIAN_P
var uintD* ptr = &((cl_heap_GV_I_bits8 *) outcast(vec))->data[index/(intDsize/8)];
@ -386,9 +386,9 @@ static const cl_I bits16_element (const cl_GV_inner<cl_I>* vec, uintL index)
}
static void bits16_set_element (cl_GV_inner<cl_I>* vec, uintL index, const cl_I& x)
{
var uint32 xval;
var uintV xval;
if (fixnump(x)) {
xval = FN_to_UL(x);
xval = FN_to_UV(x);
if (xval <= 0xFFFF) {
#if CL_CPU_BIG_ENDIAN_P
var uintD* ptr = &((cl_heap_GV_I_bits16 *) outcast(vec))->data[index/(intDsize/16)];

2
tests/timefact.cc

@ -20,6 +20,6 @@ int main (int argc, char * argv[])
cl_I m = cl_I(argv[1]);
{ CL_TIMING;
for (int rep = repetitions; rep > 0; rep--)
cl_I f = factorial(FN_to_L(m));
cl_I f = factorial(FN_to_V(m));
}
}
Loading…
Cancel
Save