]> www.ginac.de Git - cln.git/blob - src/base/cl_low.h
Avoid gcc warnings about "statement without effect".
[cln.git] / src / base / cl_low.h
1 // Low-level arithmetic: operations on 16-bit and 32-bit words
2
3 #ifndef _CL_LOW_H
4 #define _CL_LOW_H
5
6 namespace cln {
7
8 // Determines the sign of a 16-bit number.
9 // sign_of(wert)
10 // > wert: eine 16-Bit-Zahl
11 // < sint16 ergebnis: 0 falls wert>=0, -1 falls wert<0.
12 inline sint16 sign_of (sint16 wert)
13 {
14 #if defined(__sparc64__)
15         return (sint64)wert >> 63;
16 #elif defined(__sparc__) || defined(__arm__)
17         return (sint32)wert >> 31;
18 #else
19         return (wert >= 0 ? 0 : -1);
20 #endif
21 }
22
23 // Determines the sign of a 32-bit number.
24 // sign_of(wert)
25 // > wert: eine 32-Bit-Zahl
26 // < sint32 ergebnis: 0 falls wert>=0, -1 falls wert<0.
27 inline sint32 sign_of (sint32 wert)
28 {
29 #if defined(__sparc64__)
30         return (sint64)wert >> 63;
31 #elif defined(__sparc__) || defined(__arm__)
32         return wert >> 31;
33 #else
34         return (wert >= 0 ? 0 : -1);
35 #endif
36 }
37
38 #ifdef HAVE_FAST_LONGLONG
39
40 // Determines the sign of a 64-bit number.
41 // sign_of(wert)
42 // > wert: eine 64-Bit-Zahl
43 // < sint64 ergebnis: 0 falls wert>=0, -1 falls wert<0.
44 inline sint64 sign_of (sint64 wert)
45 {
46         return wert >> 63;
47 }
48
49 #endif /* HAVE_FAST_LONGLONG */
50
51
52 // High-Word einer 32-Bit-Zahl bestimmen
53 // high16(wert)
54 inline uint16 high16 (uint32 wert)
55 {
56         return wert >> 16;
57 }
58
59 // Low-Word einer 32-Bit-Zahl bestimmen
60 // low16(wert)
61 inline uint16 low16 (uint32 wert)
62 {
63         return (uint16)wert;
64 }
65
66 // Eine 32-Bit-Zahl aus ihrem High-Word und ihrem Low-Word bestimmen:
67 // highlow32(uint16 high, uint16 low)
68 inline uint32 highlow32 (uint16 high, uint16 low)
69 {
70         return ((uint32)high << 16) | (uint32)low;
71 }
72
73 // Eine 32-Bit-Zahl aus ihrem High-Word und ihrem Low-Word 0 bestimmen:
74 // highlow32_0(uint16 high)
75 inline uint32 highlow32_0 (uint16 high)
76 {
77         return (uint32)high << 16;
78 }
79
80 #ifdef HAVE_LONGLONG
81
82 // High-Word einer 64-Bit-Zahl bestimmen
83 // high32(wert)
84 inline uint32 high32 (uint64 wert)
85 {
86         return wert >> 32;
87 }
88
89 // Low-Word einer 64-Bit-Zahl bestimmen
90 // low32(wert)
91 inline uint32 low32 (uint64 wert)
92 {
93         return (uint32)wert;
94 }
95
96 // Eine 64-Bit-Zahl aus ihrem High-Word und ihrem Low-Word bestimmen:
97 // highlow64(uint32 high, uint32 low)
98 inline uint64 highlow64 (uint32 high, uint32 low)
99 {
100         return ((uint64)high << 32) | (uint64)low;
101 }
102
103 // Eine 64-Bit-Zahl aus ihrem High-Word und ihrem Low-Word 0 bestimmen:
104 // highlow64_0(uint32 high)
105 inline uint64 highlow64_0 (uint32 high)
106 {
107         return (uint64)high << 32;
108 }
109
110 #endif /* HAVE_LONGLONG */
111
112
113 // Multipliziert zwei 16-Bit-Zahlen miteinander und liefert eine 32-Bit-Zahl:
114 // mulu16(arg1,arg2)
115 // > arg1, arg2 : zwei 16-Bit-Zahlen
116 // < ergebnis: eine 32-Bit-Zahl
117 #if defined(__GNUC__) && defined(__sparc__) && !defined(__sparc64__) && defined(FAST_DOUBLE)
118 // Ist das schneller als mulu16_ ??
119 inline uint32 mulu16 (uint16 arg1, uint16 arg2)
120 {
121         union { double f; uint32 i[2]; } __fi;
122         __fi.f = (double)(sint32)arg1 * (double)(sint32)arg2
123                  + (double)(4503599627370496.0L); // + 2^52, zum Normalisieren
124         return __fi.i[1]; // untere 32 Bit herausholen (benutzt CL_CPU_BIG_ENDIAN_P !)
125 }
126 #elif defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
127 inline uint32 mulu16 (uint16 arg1, uint16 arg2)
128 {
129         register uint64 _prod;
130         __asm__("umul %1,%2,%0"
131                 : "=r" (_prod)
132                 : "r" (arg1), "r" (arg2)
133                );
134         return _prod;
135 }
136 #elif defined(__GNUC__) && (defined(__i386__) || defined(__x86_64__)) && !defined(NO_ASM)
137 inline uint32 mulu16 (uint16 arg1, uint16 arg2)
138 {
139         register uint16 _hi;
140         register uint16 _lo;
141         __asm__("mulw %2"
142                 : "=d" /* %dx */ (_hi), "=a" /* %ax */ (_lo)
143                 : "rm" (arg1), "1" /* %eax */ (arg2)
144                );
145         return highlow32(_hi,_lo);
146 }
147 #elif defined(__sparc__) || defined(__sparc64__)
148   extern "C" uint32 mulu16_ (uint16 arg1, uint16 arg2);
149   #define mulu16  mulu16_  // extern in Assembler
150 #else
151 inline uint32 mulu16 (uint16 arg1, uint16 arg2)
152 {
153         return arg1 * arg2;
154 }
155 #endif
156
157 // Multipliziert zwei 24-Bit-Zahlen zusammen und liefert eine 48-Bit-Zahl.
158 // mulu24(arg1,arg2,hi=,lo=);
159 // > arg1, arg2 : zwei 24-Bit-Zahlen
160 // < 2^32*hi+lo : eine 48-Bit-Zahl
161 #if defined(__sparc__) && !defined(__sparc64__) && defined(FAST_DOUBLE)
162   #define mulu24(x,y,hi_zuweisung,lo_zuweisung)  \
163     { var uint32 _x = (x);                                      \
164       var uint32 _y = (y);                                      \
165       var union { double f; uint32 i[2]; uint16 s[4]; } __fi;   \
166       __fi.f = (double)(sint32)(_x)*(double)(sint32)(_y)        \
167                + (double)(4503599627370496.0L); /* + 2^52, zum Normalisieren */\
168       unused (hi_zuweisung __fi.s[1]); /* mittlere 16 Bit herausholen, (benutzt CL_CPU_BIG_ENDIAN_P !) */\
169       lo_zuweisung __fi.i[1]; /* untere 32 Bit herausholen (benutzt CL_CPU_BIG_ENDIAN_P !)    */\
170     }
171 #else
172   #define mulu24  mulu32
173 #endif
174
175 // Multipliziert zwei 32-Bit-Zahlen miteinander und liefert eine 32-Bit-Zahl:
176 // mulu32_unchecked(arg1,arg2)
177 // > arg1, arg2 : zwei 32-Bit-Zahlen
178 // < ergebnis : eine 32-Bit-Zahl
179 // Es wird vorausgesetzt, daß arg1*arg2 < 2^32.
180 #if defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
181 inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
182 {
183         register uint64 _prod;
184         __asm__("umul %1,%2,%0"
185                 : "=r" (_prod)
186                 : "r" (arg1), "r" (arg2)
187                );
188         return _prod;
189 }
190 #elif defined(__sparc__) && !defined(NO_ASM)
191   extern "C" uint32 mulu32_unchecked (uint32 x, uint32 y); // extern in Assembler
192 #else
193   // Wir können dafür auch die Bibliotheksroutine des C-Compilers nehmen:
194   inline uint32 mulu32_unchecked (uint32 arg1, uint32 arg2)
195   {
196         return arg1 * arg2;
197   }
198 #endif
199
200 // Multipliziert zwei 32-Bit-Zahlen miteinander und liefert eine 64-Bit-Zahl:
201 // mulu32(arg1,arg2,hi=,lo=);
202 // > arg1, arg2 : zwei 32-Bit-Zahlen
203 // < 2^32*hi+lo : eine 64-Bit-Zahl
204   extern "C" uint32 mulu32_ (uint32 arg1, uint32 arg2); // -> Low-Teil
205   extern "C" uint32 mulu32_high;                        // -> High-Teil
206 #if defined(__GNUC__) && defined(__m68k__) && !defined(NO_ASM)
207   #define mulu32(x,y,hi_zuweisung,lo_zuweisung)  \
208     ({ var uint32 _x = (x);       \
209        var uint32 _y = (y);       \
210        var uint32 _hi;            \
211        var uint32 _lo;            \
212        __asm__("mulul %3,%0:%1" : "=d" (_hi), "=d"(_lo) : "1" (_x), "dm" (_y) ); \
213        unused (hi_zuweisung _hi); \
214        lo_zuweisung _lo;          \
215      })
216 #elif defined(__GNUC__) && defined(__m68k__)
217   #define mulu32(x,y,hi_zuweisung,lo_zuweisung)  \
218     ({ var uint32 _x = (x);                                             \
219        var uint32 _y = (y);                                             \
220        var uint16 _x1 = high16(_x);                                     \
221        var uint16 _x0 = low16(_x);                                      \
222        var uint16 _y1 = high16(_y);                                     \
223        var uint16 _y0 = low16(_y);                                      \
224        var uint32 _hi = mulu16(_x1,_y1); /* obere Portion */            \
225        var uint32 _lo = mulu16(_x0,_y0); /* untere Portion */           \
226        {var uint32 _mid = mulu16(_x0,_y1); /* 1. mittlere Portion */    \
227         _hi += high16(_mid); _mid = highlow32_0(low16(_mid));           \
228         _lo += _mid; if (_lo < _mid) { _hi += 1; } /* 64-Bit-Addition */\
229        }                                                                \
230        {var uint32 _mid = mulu16(_x1,_y0); /* 2. mittlere Portion */    \
231         _hi += high16(_mid); _mid = highlow32_0(low16(_mid));           \
232         _lo += _mid; if (_lo < _mid) { _hi += 1; } /* 64-Bit-Addition */\
233        }                                                                \
234        unused (hi_zuweisung _hi);                                       \
235        lo_zuweisung _lo;                                                \
236      })
237 #elif defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
238   #define mulu32(x,y,hi_zuweisung,lo_zuweisung)  \
239     ({ var register uint64 _prod;                               \
240        __asm__("umul %1,%2,%0"                                  \
241                : "=r" (_prod)                                   \
242                : "r" ((uint32)(x)), "r" ((uint32)(y))           \
243               );                                                \
244        unused (hi_zuweisung (uint32)(_prod>>32));               \
245        lo_zuweisung (uint32)(_prod);                            \
246      })
247 #elif defined(__GNUC__) && defined(__sparc__) && !defined(NO_ASM)
248   #define mulu32(x,y,hi_zuweisung,lo_zuweisung)  \
249     ({ lo_zuweisung mulu32_(x,y); /* extern in Assembler */     \
250       {var register uint32 _hi __asm__("%g1");                  \
251        unused (hi_zuweisung _hi);                               \
252      }})
253 #elif defined(__GNUC__) && defined(__arm__) && 0 // see comment cl_asm_arm.cc
254   #define mulu32(x,y,hi_zuweisung,lo_zuweisung)  \
255     ({ lo_zuweisung mulu32_(x,y); /* extern in Assembler */     \
256       {var register uint32 _hi __asm__("%r1"/*"%a2"*/);         \
257        unused (hi_zuweisung _hi);                               \
258      }})
259 #elif defined(__GNUC__) && (defined(__i386__) || defined(__x86_64__)) && !defined(NO_ASM)
260   #define mulu32(x,y,hi_zuweisung,lo_zuweisung)  \
261     ({ var register uint32 _hi;                                  \
262        var register uint32 _lo;                                  \
263        __asm__("mull %2"                                         \
264                : "=d" /* %edx */ (_hi), "=a" /* %eax */ (_lo)    \
265                : "g" ((uint32)(x)), "1" /* %eax */ ((uint32)(y)) \
266               );                                                 \
267        unused (hi_zuweisung _hi); lo_zuweisung _lo;              \
268      })
269 #elif defined(__GNUC__) && defined(__mips__) && !defined(NO_ASM)
270   #define mulu32(x,y,hi_zuweisung,lo_zuweisung)  \
271     ({ var register uint32 _hi;                       \
272        var register uint32 _lo;                       \
273        __asm__("multu %3,%2 ; mfhi %0 ; mflo %1"      \
274                : "=r" (_hi), "=r" (_lo)               \
275                : "r" ((uint32)(x)), "r" ((uint32)(y)) \
276               );                                      \
277        unused (hi_zuweisung _hi); lo_zuweisung _lo;   \
278      })
279 #elif defined(__GNUC__) && defined(HAVE_LONGLONG) && !defined(__arm__)
280   #define mulu32(x,y,hi_zuweisung,lo_zuweisung)  \
281     ({ var register uint64 _prod = (uint64)(uint32)(x) * (uint64)(uint32)(y); \
282        unused (hi_zuweisung (uint32)(_prod>>32));                             \
283        lo_zuweisung (uint32)(_prod);                                          \
284      })
285 #elif defined(WATCOM) && defined(__i386__) && !defined(NO_ASM)
286   #define mulu32(x,y,hi_zuweisung,lo_zuweisung)  \
287     { var register uint32 _hi;                     \
288       var register uint32 _lo;                     \
289       _lo = mulu32_(x,y), _hi = mulu32_high_();    \
290       unused (hi_zuweisung _hi); lo_zuweisung _lo; \
291     }
292   extern "C" uint32 mulu32_high_ (void);
293   #pragma aux mulu32_ = 0xF7 0xE2 /* mull %edx */ parm [eax] [edx] value [eax] modify [eax edx];
294   #pragma aux mulu32_high_ = /* */ value [edx] modify [];
295 #else
296   #define mulu32(x,y,hi_zuweisung,lo_zuweisung)  \
297     { lo_zuweisung mulu32_(x,y); unused (hi_zuweisung mulu32_high); }
298   #if (defined(__m68k__) || defined(__sparc__) || defined(__sparc64__) || defined(__arm__) || (defined(__i386__) && !defined(WATCOM) && !defined(MICROSOFT)) || defined(__x86_64__) || defined(__mips__) || defined(__hppa__)) && !defined(NO_ASM)
299     // mulu32_ extern in Assembler
300     #if defined(__sparc__) || defined(__sparc64__)
301       extern "C" uint32 _get_g1 (void);
302       #define mulu32_high  (_get_g1()) // Rückgabe im Register %g1
303     #elif !defined(__hppa__)
304       #define NEED_VAR_mulu32_high
305     #endif
306   #else
307     #define NEED_FUNCTION_mulu32_
308   #endif
309 #endif
310
311 #ifdef HAVE_FAST_LONGLONG
312
313 // Multipliziert zwei 32-Bit-Zahlen miteinander und liefert eine 64-Bit-Zahl:
314 // mulu32_w(arg1,arg2)
315 // > arg1, arg2 : zwei 32-Bit-Zahlen
316 // < result : eine 64-Bit-Zahl
317 #if defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
318   // Prefer the umul instruction over the mulx instruction (overkill).
319   #define mulu32_w(x,y)  \
320     ({ var register uint64 _prod;                               \
321        __asm__("umul %1,%2,%0"                                  \
322                : "=r" (_prod)                                   \
323                : "r" ((uint32)(x)), "r" ((uint32)(y))           \
324               );                                                \
325        _prod;                                                   \
326      })
327 #elif defined(__GNUC__)
328   #define mulu32_w(x,y)  ((uint64)(uint32)(x) * (uint64)(uint32)(y))
329 #else
330   extern "C" uint64 mulu32_w (uint32 arg1, uint32 arg2);
331   #define NEED_FUNCTION_mulu32_w
332 #endif
333
334 // Multipliziert zwei 64-Bit-Zahlen miteinander und liefert eine 128-Bit-Zahl:
335 // mulu64(arg1,arg2,hi=,lo=);
336 // > arg1, arg2 : zwei 64-Bit-Zahlen
337 // < 2^64*hi+lo : eine 128-Bit-Zahl
338   extern "C" uint64 mulu64_ (uint64 arg1, uint64 arg2); // -> Low-Teil
339   extern "C" uint64 mulu64_high;                        // -> High-Teil
340 #if defined(__GNUC__) && defined(__alpha__) && !defined(NO_ASM)
341   #define mulu64(x,y,hi_zuweisung,lo_zuweisung)  \
342     ({ var register uint64 _x = (x);    \
343        var register uint64 _y = (y);    \
344        var register uint64 _hi;         \
345        var register uint64 _lo;         \
346        __asm__("mulq %1,%2,%0"          \
347                : "=r" (_lo)             \
348                : "r" (_x), "r" (_y)     \
349               );                        \
350        __asm__("umulh %1,%2,%0"         \
351                : "=r" (_hi)             \
352                : "r" (_x), "r" (_y)     \
353               );                        \
354        hi_zuweisung _hi;                \
355        lo_zuweisung _lo;                \
356      })
357 #elif defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
358   #define mulu64(x,y,hi_zuweisung,lo_zuweisung)  \
359     ({ lo_zuweisung mulu64_(x,y); /* extern in Assembler */     \
360       {var register uint64 _hi __asm__("%g2");                  \
361        hi_zuweisung _hi;                                        \
362      }})
363 #elif defined(__GNUC__) && defined(__x86_64__) && !defined(NO_ASM)
364   #define mulu64(x,y,hi_zuweisung,lo_zuweisung)  \
365     ({ var register uint64 _hi;                                  \
366        var register uint64 _lo;                                  \
367        __asm__("mulq %2"                                         \
368                : "=d" /* %rdx */ (_hi), "=a" /* %rax */ (_lo)    \
369                : "rm" ((uint64)(x)), "1" /* %rax */ ((uint64)(y)) \
370               );                                                 \
371        hi_zuweisung _hi; lo_zuweisung _lo;                       \
372      })
373 #elif defined(__GNUC__) && defined(__ia64__) && !defined(NO_ASM)
374   #define mulu64(x,y,hi_zuweisung,lo_zuweisung)  \
375     ({ var register uint64 _x = (x);                              \
376        var register uint64 _y = (y);                              \
377        var register uint64 _hi;                                   \
378        __asm__("xma.hu %0 = %1, %2, f0"                           \
379                : "=f" (_hi)                                       \
380                : "f" ((uint64)(_x)), "f" ((uint64)(_y))           \
381               );                                                  \
382        hi_zuweisung _hi; lo_zuweisung ((uint64)(_x)*(uint64)(_y));\
383      })
384 #else
385   #define mulu64(x,y,hi_zuweisung,lo_zuweisung)  \
386     { lo_zuweisung mulu64_(x,y); hi_zuweisung mulu64_high; }
387   #if defined(__sparc64__)
388     // mulu64_ extern in Assembler
389     #if defined(__sparc64__)
390       extern "C" uint64 _get_g2 (void);
391       #define mulu64_high  (_get_g2()) // Rückgabe im Register %g2
392     #else
393       #define NEED_VAR_mulu64_high
394     #endif
395   #else
396     #define NEED_FUNCTION_mulu64_
397   #endif
398 #endif
399
400 #endif /* HAVE_FAST_LONGLONG */
401
402
403 // Dividiert eine 16-Bit-Zahl durch eine 16-Bit-Zahl und
404 // liefert einen 16-Bit-Quotienten und einen 16-Bit-Rest.
405 // divu_1616_1616(x,y,q=,r=);
406 // > uint16 x: Zähler
407 // > uint16 y: Nenner
408 // < uint16 q: floor(x/y)
409 // < uint16 r: x mod y
410 // < x = q*y+r
411   #define divu_1616_1616(x,y,q_zuweisung,r_zuweisung)  \
412     { var uint16 __x = (x);                                     \
413       var uint16 __y = (y);                                     \
414       q_zuweisung floor(__x,__y);                               \
415       r_zuweisung (__x % __y);                                  \
416     }
417
418 // Dividiert eine 32-Bit-Zahl durch eine 16-Bit-Zahl und
419 // liefert einen 16-Bit-Quotienten und einen 16-Bit-Rest.
420 // divu_3216_1616(x,y,q=,r=);
421 // > uint32 x: Zähler
422 // > uint16 y: Nenner
423 // > Es sei bekannt, daß 0 <= x < 2^16*y .
424 // < uint16 q: floor(x/y)
425 // < uint16 r: x mod y
426 // < x = q*y+r
427 #if defined(__sparc__)
428   extern "C" uint32 divu_3216_1616_ (uint32 x, uint16 y); // -> Quotient q, Rest r
429 #else
430   extern "C" uint16 divu_3216_1616_ (uint32 x, uint16 y); // -> Quotient q
431   extern "C" uint16 divu_16_rest;                         // -> Rest r
432 #endif
433 #if defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
434   #define divu_3216_1616(x,y,q_zuweisung,r_zuweisung)  \
435     ({var uint32 __x = (x);        \
436       var uint16 __y = (y);        \
437       var uint64 __q;              \
438       var uint64 __r;              \
439       __asm__ __volatile__ (       \
440         "wr %%g0,%%g0,%%y\n\t"     \
441         "udiv %2,%3,%0\n\t"        \
442         "umul %0,%3,%1\n\t"        \
443         "sub %2,%1,%1"             \
444         : "=&r" (__q), "=&r" (__r) \
445         : "r" (__x), "r" (__y));   \
446       q_zuweisung (uint16)__q;     \
447       r_zuweisung (uint16)__r;     \
448      })
449 #elif defined(__GNUC__) && (defined(__sparc__) || defined(__sparc64__)) && !defined(NO_ASM)
450   #define divu_3216_1616(x,y,q_zuweisung,r_zuweisung)  \
451     ({ var uint32 __qr = divu_3216_1616_(x,y); /* extern in Assembler */\
452        q_zuweisung low16(__qr);                                         \
453        r_zuweisung high16(__qr);                                        \
454      })
455 #elif defined(__GNUC__) && defined(__m68k__) && !defined(NO_ASM)
456   #define divu_3216_1616(x,y,q_zuweisung,r_zuweisung)  \
457     ({var uint32 __x = (x);                                             \
458       var uint16 __y = (y);                                             \
459       var uint32 __qr;                                                  \
460       __asm__ __volatile__ ("                                           \
461         divu %2,%0                                                      \
462         " : "=d" (__qr) : "0" (__x), "dm" (__y));                       \
463       q_zuweisung low16(__qr);                                          \
464       r_zuweisung high16(__qr);                                         \
465      })
466 #elif defined(__GNUC__) && (defined(__i386__) || defined(__x86_64__)) && !defined(NO_ASM)
467   #define divu_3216_1616(x,y,q_zuweisung,r_zuweisung)  \
468     ({var uint32 __x = (x);                                             \
469       var uint16 __y = (y);                                             \
470       var uint16 __q;                                                   \
471       var uint16 __r;                                                   \
472       __asm__("divw %4"                                                 \
473               : "=a" /* %ax */ (__q), "=d" /* %dx */ (__r)              \
474               : "1" /* %dx */ ((uint16)(high16(__x))), "0" /* %ax */ ((uint16)(low16(__x))), "rm" (__y) \
475              );                                                         \
476       q_zuweisung __q;                                                  \
477       r_zuweisung __r;                                                  \
478      })
479 #elif defined(__GNUC__) && defined(__arm__) && 0 // see comment cl_asm_arm.cc
480   #define divu_3216_1616(x,y,q_zuweisung,r_zuweisung)  \
481     { var uint32 __q = divu_3216_1616_(x,y); /* extern in Assembler */  \
482       var register uint32 __r __asm__("%r1"/*"%a2"*/);                  \
483       q_zuweisung __q; r_zuweisung __r;                                 \
484     }
485 #elif defined(__GNUC__) && !defined(__arm__)
486   #define divu_3216_1616(x,y,q_zuweisung,r_zuweisung)  \
487     ({var uint32 __x = (x);                                             \
488       var uint16 __y = (y);                                             \
489       var uint16 __q = floor(__x,__y);                                  \
490       q_zuweisung __q;                                                  \
491       r_zuweisung (__x - __q * __y);                                    \
492      })
493 #elif (defined(__sparc__) || defined(__sparc64__)) && !defined(NO_ASM)
494   #define divu_3216_1616(x,y,q_zuweisung,r_zuweisung)  \
495     { var uint32 __qr = divu_3216_1616_(x,y); /* extern in Assembler */ \
496       q_zuweisung low16(__qr);                                          \
497       r_zuweisung high16(__qr);                                         \
498     }
499 #elif defined(__arm__) && !defined(NO_ASM)
500   #define divu_3216_1616(x,y,q_zuweisung,r_zuweisung)  \
501     { q_zuweisung divu_3216_1616_(x,y); /* extern in Assembler */       \
502       r_zuweisung divu_16_rest;                                         \
503     }
504   #define NEED_VAR_divu_16_rest
505 #else
506   #define divu_3216_1616(x,y,q_zuweisung,r_zuweisung)  \
507     { q_zuweisung divu_3216_1616_(x,y); r_zuweisung divu_16_rest; }
508   #define NEED_FUNCTION_divu_3216_1616_
509 #endif
510
511 // Dividiert eine 32-Bit-Zahl durch eine 16-Bit-Zahl und
512 // liefert einen 32-Bit-Quotienten und einen 16-Bit-Rest.
513 // divu_3216_3216(x,y,q=,r=);
514 // > uint32 x: Zähler
515 // > uint16 y: Nenner
516 // Es sei bekannt, daß y>0.
517 // < uint32 q: floor(x/y)
518 // < uint16 r: x mod y
519 // < x = q*y+r
520   extern "C" uint32 divu_3216_3216_ (uint32 x, uint16 y); // -> Quotient q
521   extern "C" uint16 divu_16_rest;                         // -> Rest r
522 #if defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
523   #define divu_3216_3216(x,y,q_zuweisung,r_zuweisung)  \
524     ({var uint32 __x = (x);        \
525       var uint16 __y = (y);        \
526       var uint64 __q;              \
527       var uint64 __r;              \
528       __asm__ __volatile__ (       \
529         "wr %%g0,%%g0,%%y\n\t"     \
530         "udiv %2,%3,%0\n\t"        \
531         "umul %0,%3,%1\n\t"        \
532         "sub %2,%1,%1"             \
533         : "=&r" (__q), "=&r" (__r) \
534         : "r" (__x), "r" (__y));   \
535       q_zuweisung (uint32)__q;     \
536       r_zuweisung (uint16)__r;     \
537      })
538 #elif defined(__sparc__) || defined(__sparc64__) || defined(__i386__) || defined(__x86_64__)
539   #define divu_3216_3216  divu_3232_3232
540 #else
541   // Methode: (beta = 2^16)
542   // x = x1*beta+x0 schreiben.
543   // Division mit Rest: x1 = q1*y + r1, wobei 0 <= x1 < beta <= beta*y.
544   // Also 0 <= q1 < beta, 0 <= r1 < y.
545   // Division mit Rest: (r1*beta+x0) = q0*y + r0, wobei 0 <= r1*beta+x0 < beta*y.
546   // Also 0 <= q0 < beta, 0 <= r0 < y
547   // und x = x1*beta+x0 = (q1*beta+q0)*y + r0.
548   // Setze q := q1*beta+q0 und r := r0.
549   #define divu_3216_3216(x,y,q_zuweisung,r_zuweisung)  \
550     { var uint32 _x = (x);                                              \
551       var uint16 _y = (y);                                              \
552       var uint16 _q1;                                                   \
553       var uint16 _q0;                                                   \
554       var uint16 _r1;                                                   \
555       divu_3216_1616(high16(_x),_y, _q1 = , _r1 = );                    \
556       divu_3216_1616(highlow32(_r1,low16(_x)),_y, _q0 = , r_zuweisung); \
557       q_zuweisung highlow32(_q1,_q0);                                   \
558     }
559 #endif
560
561 // Dividiert eine 32-Bit-Zahl durch eine 32-Bit-Zahl und
562 // liefert einen 32-Bit-Quotienten und einen 32-Bit-Rest.
563 // divu_3232_3232(x,y,q=,r=);
564 // > uint32 x: Zähler
565 // > uint32 y: Nenner
566 // Es sei bekannt, daß y>0.
567 // < uint32 q: floor(x/y)
568 // < uint32 r: x mod y
569 // < x = q*y+r
570   extern "C" uint32 divu_3232_3232_ (uint32 x, uint32 y); // -> Quotient q
571   extern "C" uint32 divu_32_rest;                         // -> Rest r
572 #if defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
573   #define divu_3232_3232(x,y,q_zuweisung,r_zuweisung)  \
574     ({var uint32 __x = (x);        \
575       var uint32 __y = (y);        \
576       var uint64 __q;              \
577       var uint64 __r;              \
578       __asm__ __volatile__ (       \
579         "wr %%g0,%%g0,%%y\n\t"     \
580         "udiv %2,%3,%0\n\t"        \
581         "umul %0,%3,%1\n\t"        \
582         "sub %2,%1,%1"             \
583         : "=&r" (__q), "=&r" (__r) \
584         : "r" (__x), "r" (__y));   \
585       q_zuweisung (uint32)__q;     \
586       r_zuweisung (uint32)__r;     \
587      })
588 #elif defined(__sparc__) || defined(__sparc64__) || defined(__i386__) || defined(__x86_64__)
589   #define divu_3232_3232(x,y,q_zuweisung,r_zuweisung)  \
590     divu_6432_3232(0,x,y,q_zuweisung,r_zuweisung)
591   #define divu_3232_3232_(x,y) divu_6432_3232_(0,x,y)
592 #else
593   // Methode: (beta = 2^n = 2^16, n = 16)
594   // Falls y < beta, handelt es sich um eine 32-durch-16-Bit-Division.
595   // Falls y >= beta:
596   // Quotient  q = floor(x/y) < beta  (da 0 <= x < beta^2, y >= beta).
597   // y habe genau n+k Bits (1 <= k <= n), d.h. 2^(n+k-1) <= y < 2^(n+k).
598   // Schreibe  x = 2^k*x1 + x0  mit  x1 := floor(x/2^k)
599   // und       y = 2^k*y1 + y0  mit  y1 := floor(y/2^k)
600   // und bilde den Näherungs-Quotienten floor(x1/y1)
601   // oder (noch besser) floor(x1/(y1+1)).
602   // Wegen 0 <= x1 < 2^(2n) und 0 < 2^(n-1) <= y1 < 2^n
603   // und  x1/(y1+1) <= x/y < x1/(y1+1) + 2
604   // (denn x1/(y1+1) = (x1*2^k)/((y1+1)*2^k) <= (x1*2^k)/y <= x/y
605   // und x/y - x1/(y1+1) = (x+x*y1-x1*y)/(y*(y1+1))
606   // = (x+x0*y1-x1*y0)/(y*(y1+1)) <= (x+x0*y1)/(y*(y1+1))
607   // <= x/(y*(y1+1)) + x0/y
608   // <= 2^(2n)/(2^(n+k-1)*(2^(n-1)+1)) + 2^k/2^(n+k-1)
609   // = 2^(n-k+1)/(2^(n-1)+1) + 2^(1-n) <= 2^n/(2^(n-1)+1) + 2^(1-n) < 2 )
610   // gilt  floor(x1/(y1+1)) <= floor(x/y) <= floor(x1/(y1+1)) + 2  .
611   // Man bildet also  q:=floor(x1/(y1+1))  (ein Shift um n Bit oder
612   // eine (2n)-durch-n-Bit-Division, mit Ergebnis q <= floor(x/y) < beta)
613   // und x-q*y und muß hiervon noch höchstens 2 mal y abziehen und q
614   // incrementieren, um den Quotienten  q = floor(x/y)  und den Rest
615   // x-floor(x/y)*y  der Division zu bekommen.
616   #define divu_3232_3232(x,y,q_zuweisung,r_zuweisung)  \
617     { var uint32 _x = (x);                                              \
618       var uint32 _y = (y);                                              \
619       if (_y <= (uint32)(bit(16)-1))                                    \
620         { var uint16 _q1;                                               \
621           var uint16 _q0;                                               \
622           var uint16 _r1;                                               \
623           divu_3216_1616(high16(_x),_y, _q1 = , _r1 = );                \
624           divu_3216_1616(highlow32(_r1,low16(_x)),_y, _q0 = , r_zuweisung); \
625           q_zuweisung highlow32(_q1,_q0);                               \
626         }                                                               \
627         else                                                            \
628         { var uint32 _x1 = _x; /* x1 := x */                            \
629           var uint32 _y1 = _y; /* y1 := y */                            \
630           var uint16 _q;                                                \
631           do { _x1 = floor(_x1,2); _y1 = floor(_y1,2); } /* k erhöhen */\
632              until (_y1 <= (uint32)(bit(16)-1)); /* bis y1 < beta */    \
633           { var uint16 _y2 = low16(_y1)+1; /* y1+1 bilden */            \
634             if (_y2==0)                                                 \
635               { _q = high16(_x1); } /* y1+1=beta -> ein Shift */        \
636               else                                                      \
637               { divu_3216_1616(_x1,_y2,_q=,); } /* Division von x1 durch y1+1 */\
638           }                                                             \
639           /* _q = q = floor(x1/(y1+1)) */                               \
640           /* x-q*y bilden (eine 16-mal-32-Bit-Multiplikation ohne Überlauf): */\
641           _x -= highlow32_0(mulu16(_q,high16(_y))); /* q * high16(y) * beta */\
642           /* gefahrlos, da q*high16(y) <= q*y/beta <= x/beta < beta */  \
643           _x -= mulu16(_q,low16(_y)); /* q * low16(y) */                \
644           /* gefahrlos, da q*high16(y)*beta + q*low16(y) = q*y <= x */  \
645           /* Noch höchstens 2 mal y abziehen: */                       \
646           if (_x >= _y)                                                 \
647             { _q += 1; _x -= _y;                                        \
648               if (_x >= _y)                                             \
649                 { _q += 1; _x -= _y; }                                  \
650             }                                                           \
651           r_zuweisung _x;                                               \
652           q_zuweisung (uint32)(_q);                                     \
653     }   }
654   #define NEED_FUNCTION_divu_3232_3232_
655 #endif
656
657 // Dividiert eine 64-Bit-Zahl durch eine 32-Bit-Zahl und
658 // liefert einen 32-Bit-Quotienten und einen 32-Bit-Rest.
659 // divu_6432_3232(xhi,xlo,y,q=,r=);
660 // > uint32 xhi,xlo: x = 2^32*xhi+xlo = Zähler
661 // > uint32 y: Nenner
662 // > Es sei bekannt, daß 0 <= x < 2^32*y .
663 // < uint32 q: floor(x/y)
664 // < uint32 r: x mod y
665 // < x = q*y+r
666   extern "C" uint32 divu_6432_3232_ (uint32 xhi, uint32 xlo, uint32 y); // -> Quotient q
667   extern "C" uint32 divu_32_rest;                                       // -> Rest r
668 #if defined(__GNUC__) && defined(__m68k__) && !defined(NO_ASM)
669   #define divu_6432_3232(xhi,xlo,y,q_zuweisung,r_zuweisung)  \
670     ({var uint32 __xhi = (xhi);                                         \
671       var uint32 __xlo = (xlo);                                         \
672       var uint32 __y = (y);                                             \
673       var uint32 __q;                                                   \
674       var uint32 __r;                                                   \
675       __asm__ __volatile__ ("                                           \
676         divul %4,%1:%0                                                  \
677         " : "=d" (__q), "=d" (__r) : "1" (__xhi), "0" (__xlo), "dm" (__y)); \
678       q_zuweisung __q;                                                  \
679       r_zuweisung __r;                                                  \
680      })
681   #define divu_6432_3232_(xhi,xlo,y) \
682     ({var uint32 ___q; divu_6432_3232(xhi,xlo,y,___q=,); ___q; })
683 #elif defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
684   #define divu_6432_3232(xhi,xlo,y,q_zuweisung,r_zuweisung)  \
685     ({var uint32 __xhi = (xhi);    \
686       var uint32 __xlo = (xlo);    \
687       var uint32 __y = (y);        \
688       var uint64 __q;              \
689       var uint64 __r;              \
690       __asm__ __volatile__ (       \
691         "wr %2,%%g0,%%y\n\t"       \
692         "udiv %3,%4,%0\n\t"        \
693         "umul %0,%4,%1\n\t"        \
694         "sub %3,%1,%1"             \
695         : "=&r" (__q), "=&r" (__r) \
696         : "r" (__xhi), "r" (__xlo), "r" (__y)); \
697       q_zuweisung (uint32)__q;     \
698       r_zuweisung (uint32)__r;     \
699      })
700 #elif defined(__GNUC__) && (defined(__sparc__) || defined(__sparc64__)) && !defined(NO_ASM)
701   #define divu_6432_3232(xhi,xlo,y,q_zuweisung,r_zuweisung)  \
702     ({ var uint32 _q = divu_6432_3232_(xhi,xlo,y); /* extern in Assembler */\
703        var register uint32 _r __asm__("%g1");                               \
704        q_zuweisung _q; r_zuweisung _r;                                      \
705      })
706 #elif defined(__GNUC__) && defined(__arm__) && 0 // see comment cl_asm_arm.cc
707   #define divu_6432_3232(xhi,xlo,y,q_zuweisung,r_zuweisung)  \
708     ({ var uint32 _q = divu_6432_3232_(xhi,xlo,y); /* extern in Assembler */\
709        var register uint32 _r __asm__("%r1"/*"%a2"*/);                      \
710        q_zuweisung _q; r_zuweisung _r;                                      \
711      })
712 #elif defined(__GNUC__) && (defined(__i386__) || defined(__x86_64__)) && !defined(NO_ASM)
713   #define divu_6432_3232(xhi,xlo,y,q_zuweisung,r_zuweisung)  \
714     ({var uint32 __xhi = (xhi);                                         \
715       var uint32 __xlo = (xlo);                                         \
716       var uint32 __y = (y);                                             \
717       var uint32 __q;                                                   \
718       var uint32 __r;                                                   \
719       __asm__ __volatile__ (                                            \
720          "divl %4"                                                      \
721          : "=a" /* %eax */ (__q), "=d" /* %edx */ (__r)                 \
722          : "1" /* %edx */ (__xhi), "0" /* %eax */ (__xlo), "rm" (__y)   \
723          );                                                             \
724       q_zuweisung __q;                                                  \
725       r_zuweisung __r;                                                  \
726      })
727   #define divu_6432_3232_(xhi,xlo,y) \
728     ({var uint32 ___q; divu_6432_3232(xhi,xlo,y,___q=,); ___q; })
729 #elif defined(__GNUC__) && defined(HAVE_LONGLONG) && !defined(__arm__)
730   #define divu_6432_3232(xhi,xlo,y,q_zuweisung,r_zuweisung) \
731     ({var uint32 __xhi = (xhi);                                         \
732       var uint32 __xlo = (xlo);                                         \
733       var uint64 __x = ((uint64)__xhi << 32) | (uint64)__xlo;           \
734       var uint32 __y = (y);                                             \
735       var uint32 __q = floor(__x,(uint64)__y);                          \
736       q_zuweisung __q; r_zuweisung __xlo - __q * __y;                   \
737      })
738   #define divu_6432_3232_(xhi,xlo,y) \
739     ({var uint32 ___q; divu_6432_3232(xhi,xlo,y,___q=,); ___q; })
740 #elif defined(WATCOM) && defined(__i386__) && !defined(NO_ASM)
741   #define divu_6432_3232(xhi,xlo,y,q_zuweisung,r_zuweisung)  \
742     { var uint32 __xhi = (xhi);                                         \
743       var uint32 __xlo = (xlo);                                         \
744       var uint32 __y = (y);                                             \
745       var uint32 __q;                                                   \
746       var uint32 __r;                                                   \
747       __q = divu_6432_3232_(__xhi,__xlo,__y); __r = divu_6432_3232_rest(); \
748       q_zuweisung __q;                                                  \
749       r_zuweisung __r;                                                  \
750     }
751   extern "C" uint32 divu_6432_3232_rest (void);
752   #pragma aux divu_6432_3232_ = 0xF7 0xF1 /* divl %ecx */ parm [edx] [eax] [ecx] value [eax] modify [eax edx];
753   #pragma aux divu_6432_3232_rest = /* */ value [edx] modify [];
754 #else
755   #define divu_6432_3232(xhi,xlo,y,q_zuweisung,r_zuweisung)  \
756     { q_zuweisung divu_6432_3232_(xhi,xlo,y); r_zuweisung divu_32_rest; }
757   #if (defined(__m68k__) || defined(__sparc__) || defined(__sparc64__) || defined(__arm__) || (defined(__i386__) && !defined(WATCOM) && !defined(MICROSOFT)) || defined(__x86_64__) || defined(__hppa__)) && !defined(NO_ASM)
758     // divu_6432_3232_ extern in Assembler
759     #if defined(__sparc__) || defined(__sparc64__)
760       extern "C" uint32 _get_g1 (void);
761       #define divu_32_rest  (_get_g1()) // Rückgabe im Register %g1
762     #else
763       #define NEED_VAR_divu_32_rest
764     #endif
765   #else
766     #define NEED_FUNCTION_divu_6432_3232_
767   #endif
768 #endif
769
770 #ifdef HAVE_FAST_LONGLONG
771
772 // Dividiert eine 64-Bit-Zahl durch eine 32-Bit-Zahl und
773 // liefert einen 32-Bit-Quotienten und einen 32-Bit-Rest.
774 // divu_6432_3232_w(x,y,q=,r=);
775 // > uint64 x: Zähler
776 // > uint32 y: Nenner
777 // > Es sei bekannt, daß 0 <= x < 2^32*y .
778 // < uint32 q: floor(x/y)
779 // < uint32 r: x mod y
780 // < x = q*y+r
781 #if defined(__GNUC__) && defined(__sparc64__) && !defined(NO_ASM)
782   // Prefer the udiv and umul instructions over the udivx and mulx instructions
783   // (overkill).
784   #define divu_6432_3232_w(x,y,q_zuweisung,r_zuweisung)  \
785     ({var uint64 __x = (x);           \
786       var uint32 __xhi = high32(__x); \
787       var uint32 __xlo = low32(__x);  \
788       var uint32 __y = (y);           \
789       var uint64 __q;                 \
790       var uint64 __r;                 \
791       __asm__ __volatile__ (          \
792         "wr %2,%%g0,%%y\n\t"          \
793         "udiv %3,%4,%0\n\t"           \
794         "umul %0,%4,%1\n\t"           \
795         "sub %3,%1,%1"                \
796         : "=&r" (__q), "=&r" (__r)    \
797         : "r" (__xhi), "r" (__xlo), "r" (__y)); \
798       q_zuweisung (uint32)__q;        \
799       r_zuweisung (uint32)__r;        \
800      })
801 #elif defined(__GNUC__) && (defined(__alpha__) || defined(__ia64__) || defined(__mips64__) || defined(__sparc64__))
802   // On __alpha__, computing the remainder by multiplication is just two
803   // instructions, compared to the __remqu (libc) function call for the %
804   // operator.
805   // On __ia64__, computing the remainder by multiplication is just four
806   // instructions, compared to the __umoddi3 (libgcc) function call for the %
807   // operator.
808   // On __mips64__, computing the remainder by multiplication is just two
809   // instructions, compared to the __umoddi3 (libgcc) function call for the %
810   // operator.
811   // On __sparc64__, computing the remainder by multiplication uses a 32-bit
812   // multiplication instruction, compared to a 64-bit multiplication when the %
813   // operator is used.
814   #define divu_6432_3232_w(x,y,q_zuweisung,r_zuweisung)  \
815     ({var uint64 __x = (x);                                             \
816       var uint32 __y = (y);                                             \
817       var uint32 __q = floor(__x,(uint64)__y);                          \
818       q_zuweisung __q; r_zuweisung (uint32)__x - __q * __y;             \
819      })
820 #elif defined(__GNUC__) && defined(__x86_64__)
821   // On __x86_64__, gcc 4.0 performs both quotient and remainder computation
822   // in a single instruction.
823   #define divu_6432_3232_w(x,y,q_zuweisung,r_zuweisung)  \
824     ({var uint64 __x = (x);                                             \
825       var uint32 __y = (y);                                             \
826       var uint32 __q = floor(__x,(uint64)__y);                          \
827       q_zuweisung __q; r_zuweisung __x % (uint64)__y;                   \
828      })
829 #else
830   #define divu_6432_3232_w(x,y,q_zuweisung,r_zuweisung)  \
831     { var uint64 __x = (x);                                               \
832       divu_6432_3232(high32(__x),low32(__x),(y),q_zuweisung,r_zuweisung); \
833     }
834 #endif
835
836 // Dividiert eine 64-Bit-Zahl durch eine 32-Bit-Zahl und
837 // liefert einen 64-Bit-Quotienten und einen 32-Bit-Rest.
838 // divu_6432_6432(x,y,q=,r=);
839 // > uint64 x: Zähler
840 // > uint32 y: Nenner
841 // > Es sei bekannt, daß y>0.
842 // < uint64 q: floor(x/y)
843 // < uint32 r: x mod y
844 // < x = q*y+r
845 #if defined(__GNUC__) && (defined(__alpha__) || defined(__ia64__) || defined(__mips64__) || defined(__sparc64__))
846   // On __alpha__, computing the remainder by multiplication is just two
847   // instructions, compared to the __remqu (libc) function call for the %
848   // operator.
849   // On __ia64__, computing the remainder by multiplication is just four
850   // instructions, compared to the __umoddi3 (libgcc) function call for the %
851   // operator.
852   // On __mips64__, computing the remainder by multiplication is just two
853   // instructions, compared to the __umoddi3 (libgcc) function call for the %
854   // operator.
855   // On __sparc64__, computing the remainder by multiplication uses a 32-bit
856   // multiplication instruction, compared to a 64-bit multiplication when the %
857   // operator is used.
858   #define divu_6432_6432(x,y,q_zuweisung,r_zuweisung)  \
859     ({var uint64 _x = (x);                    \
860       var uint32 _y = (y);                    \
861       var uint64 _q;                          \
862       q_zuweisung _q = floor(_x,(uint64)_y);  \
863       r_zuweisung low32(_x) - low32(_q) * _y; \
864      })
865 #elif defined(__GNUC__) && defined(__x86_64__)
866   // On __x86_64__, gcc 4.0 performs both quotient and remainder computation
867   // in a single instruction.
868   #define divu_6432_6432(x,y,q_zuweisung,r_zuweisung)  \
869     ({var uint64 _x = (x);               \
870       var uint32 _y = (y);               \
871       q_zuweisung floor(_x,(uint64)_y);  \
872       r_zuweisung _x % (uint64)_y;       \
873      })
874 #else
875   // Methode: (beta = 2^32)
876   // x = x1*beta+x0 schreiben.
877   // Division mit Rest: x1 = q1*y + r1, wobei 0 <= x1 < beta <= beta*y.
878   // Also 0 <= q1 < beta, 0 <= r1 < y.
879   // Division mit Rest: (r1*beta+x0) = q0*y + r0, wobei 0 <= r1*beta+x0 < beta*y.
880   // Also 0 <= q0 < beta, 0 <= r0 < y
881   // und x = x1*beta+x0 = (q1*beta+q0)*y + r0.
882   // Setze q := q1*beta+q0 und r := r0.
883   #if defined(__GNUC__)
884     #define divu_6432_6432(x,y,q_zuweisung,r_zuweisung)  \
885       ({var uint64 _x = (x);            \
886         var uint32 _y = (y);            \
887         var uint32 _q1;                 \
888         var uint32 _q0;                 \
889         var uint32 _r1;                 \
890         divu_6432_3232(0,high32(_x),_y, _q1 = , _r1 = ); \
891         divu_6432_3232(_r1,low32(_x),_y, _q0 = , r_zuweisung); \
892         q_zuweisung highlow64(_q1,_q0); \
893        })
894   #else
895     #define divu_6432_6432(x,y,q_zuweisung,r_zuweisung)  \
896       {var uint64 _x = (x);            \
897        var uint32 _y = (y);            \
898        var uint32 _q1;                 \
899        var uint32 _q0;                 \
900        var uint32 _r1;                 \
901        divu_6432_3232(0,high32(_x),_y, _q1 = , _r1 = ); \
902        divu_6432_3232(_r1,low32(_x),_y, _q0 = , r_zuweisung); \
903        q_zuweisung highlow64(_q1,_q0); \
904       }
905   #endif
906 #endif
907
908 // Dividiert eine 64-Bit-Zahl durch eine 64-Bit-Zahl und
909 // liefert einen 64-Bit-Quotienten und einen 64-Bit-Rest.
910 // divu_6464_6464(x,y,q=,r=);
911 // > uint64 x: Zähler
912 // > uint64 y: Nenner
913 // > Es sei bekannt, daß y>0.
914 // < uint64 q: floor(x/y)
915 // < uint64 r: x mod y
916 // < x = q*y+r
917 #if defined(__GNUC__) && (defined(__alpha__) || defined(__ia64__) || defined(__mips64__) || defined(__sparc64__))
918   // On __alpha__, computing the remainder by multiplication is just two
919   // instructions, compared to the __remqu (libc) function call for the %
920   // operator.
921   // On __ia64__, computing the remainder by multiplication is just four
922   // instructions, compared to the __umoddi3 (libgcc) function call for the %
923   // operator.
924   // On __mips64__, computing the remainder by multiplication is just two
925   // instructions, compared to the __umoddi3 (libgcc) function call for the %
926   // operator.
927   // On __sparc64__, it doesn't matter.
928   #define divu_6464_6464(x,y,q_zuweisung,r_zuweisung)  \
929     ({var uint64 _x = (x);           \
930       var uint64 _y = (y);           \
931       var uint64 _q;                 \
932       q_zuweisung _q = floor(_x,_y); \
933       r_zuweisung _x - _q * _y;      \
934      })
935 #elif defined(__GNUC__) && (defined(__sparc64__) || defined(__x86_64__))
936   // On __sparc64__, it doesn't matter.
937   // On __x86_64__, gcc 4.0 performs both quotient and remainder computation
938   // in a single instruction.
939   #define divu_6464_6464(x,y,q_zuweisung,r_zuweisung)  \
940     ({var uint64 _x = (x);      \
941       var uint64 _y = (y);      \
942       q_zuweisung floor(_x,_y); \
943       r_zuweisung _x % _y;      \
944      })
945 #else
946   // For unknown CPUs, we don't know whether gcc's __udivdi3 function plus a
947   // multiplication is slower or faster than our own divu_6464_6464_ routine.
948   // Anyway, call our own routine.
949   extern "C" uint64 divu_6464_6464_ (uint64 x, uint64 y); // -> Quotient q
950   extern "C" uint64 divu_64_rest;                         // -> Rest r
951   #define divu_6464_6464(x,y,q_zuweisung,r_zuweisung)  \
952     { q_zuweisung divu_6464_6464_(x,y); r_zuweisung divu_64_rest; }
953   #define NEED_VAR_divu_64_rest
954   #define NEED_FUNCTION_divu_6464_6464_
955 #endif
956
957 // Dividiert eine 128-Bit-Zahl durch eine 64-Bit-Zahl und
958 // liefert einen 64-Bit-Quotienten und einen 64-Bit-Rest.
959 // divu_12864_6464(xhi,xlo,y,q=,r=);
960 // > uint64 xhi,xlo: x = 2^64*xhi+xlo = Zähler
961 // > uint64 y: Nenner
962 // > Es sei bekannt, daß 0 <= x < 2^64*y .
963 // < uint64 q: floor(x/y)
964 // < uint64 r: x mod y
965 // < x = q*y+r
966   extern "C" uint64 divu_12864_6464_ (uint64 xhi, uint64 xlo, uint64 y); // -> Quotient q
967   extern "C" uint64 divu_64_rest;                                        // -> Rest r
968 #if defined(__GNUC__) && defined(__x86_64__) && !defined(NO_ASM)
969   #define divu_12864_6464(xhi,xlo,y,q_zuweisung,r_zuweisung)  \
970     ({var uint64 __xhi = (xhi);                                         \
971       var uint64 __xlo = (xlo);                                         \
972       var uint64 __y = (y);                                             \
973       var uint64 __q;                                                   \
974       var uint64 __r;                                                   \
975       __asm__ __volatile__ (                                            \
976          "divq %4"                                                      \
977          : "=a" /* %rax */ (__q), "=d" /* %rdx */ (__r)                 \
978          : "1" /* %rdx */ (__xhi), "0" /* %rax */ (__xlo), "rm" (__y)   \
979          );                                                             \
980       q_zuweisung __q;                                                  \
981       r_zuweisung __r;                                                  \
982      })
983   #define divu_12864_64364_(xhi,xlo,y) \
984     ({var uint64 ___q; divu_12864_6464(xhi,xlo,y,___q=,); ___q; })
985 #else
986   #define divu_12864_6464(xhi,xlo,y,q_zuweisung,r_zuweisung)  \
987     { q_zuweisung divu_12864_6464_(xhi,xlo,y); r_zuweisung divu_64_rest; }
988   #define NEED_VAR_divu_64_rest
989   #define NEED_FUNCTION_divu_12864_6464_
990 #endif
991
992 #endif /* HAVE_FAST_LONGLONG */
993
994
995 // Zieht die Ganzzahl-Wurzel aus einer 32-Bit-Zahl und
996 // liefert eine 16-Bit-Wurzel und einen Rest.
997 // isqrt_32_16(x,y=,sqrtp=);
998 // > uint32 x: Radikand, >= 2^30, < 2^32
999 // < uint16 y: floor(sqrt(x)), >= 2^15, < 2^16
1000 // < boolean sqrtp: /=0, falls x=y^2
1001   // Methode:
1002   // y := 2^16 als Anfangswert,
1003   // y := floor((y + floor(x/y))/2) als nächster Wert,
1004   // solange z := floor(x/y) < y, setze y := floor((y+z)/2).
1005   // y ist fertig; x=y^2 genau dann, wenn z=y und die letzte Division aufging.
1006   // (Beweis:
1007   //  1. Die Folge der y ist streng monoton fallend.
1008   //  2. Stets gilt y >= floor(sqrt(x)) (denn für alle y>0 ist
1009   //     y + x/y >= 2*sqrt(x) und daher  floor((y + floor(x/y))/2) =
1010   //     floor(y/2 + x/(2*y)) >= floor(sqrt(x)) ).
1011   //  3. Am Schluß gilt x >= y^2.
1012   // )
1013   #define isqrt_32_16(x,y_zuweisung,sqrtp_zuweisung)  \
1014     { var uint32 _x = (x);                                              \
1015       var uint16 _x1 = high16(_x);                                      \
1016       var uint16 _y = floor(_x1,2) | bit(16-1);                         \
1017       loop                                                              \
1018         { var uint16 _z;                                                \
1019           var uint16 _r;                                                \
1020           if (_x1 >= _y) /* Division _x/_y ergäbe Überlauf -> _z > _y */\
1021             { unused (sqrtp_zuweisung FALSE); break; }                  \
1022           divu_3216_1616(_x,_y, _z=,_r=); /* Dividiere _x/_y */         \
1023           if (_z >= _y)                                                 \
1024             { unused (sqrtp_zuweisung (_z == _y) && (_r == 0)); break; } \
1025           _y = floor((uint16)(_z+_y),2) | bit(16-1); /* _y muß >= 2^15 bleiben */\
1026         }                                                               \
1027       y_zuweisung _y;                                                   \
1028     }
1029
1030 // Zieht die Ganzzahl-Wurzel aus einer 64-Bit-Zahl und
1031 // liefert eine 32-Bit-Wurzel und einen Rest.
1032 // isqrt_64_32(xhi,xlo,y=,sqrtp=);
1033 // > uint32 xhi,xlo: Radikand x = 2^32*xhi+xlo, >= 2^62, < 2^64
1034 // < uint32 y: floor(sqrt(x)), >= 2^31, < 2^32
1035 // < boolean sqrtp: /=0, falls x=y^2
1036 #if defined(__sparc__) || defined(__sparc64__) || defined(__m68k__) || defined(__hppa__)
1037   // Methode:
1038   // y := 2^32 als Anfangswert,
1039   // y := floor((y + floor(x/y))/2) als nächster Wert,
1040   // solange z := floor(x/y) < y, setze y := floor((y+z)/2).
1041   // y ist fertig; x=y^2 genau dann, wenn z=y und die letzte Division aufging.
1042   // (Beweis:
1043   //  1. Die Folge der y ist streng monoton fallend.
1044   //  2. Stets gilt y >= floor(sqrt(x)) (denn für alle y>0 ist
1045   //     y + x/y >= 2*sqrt(x) und daher  floor((y + floor(x/y))/2) =
1046   //     floor(y/2 + x/(2*y)) >= floor(sqrt(x)) ).
1047   //  3. Am Schluß gilt x >= y^2.
1048   // )
1049   #define isqrt_64_32(xhi,xlo,y_zuweisung,sqrtp_zuweisung)  \
1050     { var uint32 _xhi = (xhi);                                          \
1051       var uint32 _xlo = (xlo);                                          \
1052       var uint32 _y = floor(_xhi,2) | bit(32-1);                        \
1053       loop                                                              \
1054         { var uint32 _z;                                                \
1055           var uint32 _rest;                                             \
1056           if (_xhi >= _y) /* Division _x/_y ergäbe Überlauf -> _z > _y */\
1057             { sqrtp_zuweisung FALSE; break; }                           \
1058           divu_6432_3232(_xhi,_xlo,_y, _z=,_rest=); /* Dividiere _x/_y */\
1059           if (_z >= _y)                                                 \
1060             { sqrtp_zuweisung (_z == _y) && (_rest == 0); break; }      \
1061           _y = floor(_z+_y,2) | bit(32-1); /* _y muß >= 2^31 bleiben */        \
1062         }                                                               \
1063       y_zuweisung _y;                                                   \
1064     }
1065 #else
1066   // Methode:
1067   // Wie bei UDS_sqrt mit n=2.
1068   // y = 2^16*yhi + ylo ansetzen.
1069   // Dann muß
1070   //   yhi = floor(y/2^16) = floor(floor(sqrt(x))/2^16)
1071   //       = floor(sqrt(x)/2^16) = floor(sqrt(x/2^32)) = isqrt(xhi)
1072   // sein. Es folgt yhi >= 2^15.
1073   // Danach sucht man das größte ylo >=0 mit
1074   // x - 2^32*yhi^2 >= 2*2^16*yhi*ylo + ylo^2.
1075   // Dazu setzen wir  xhi*2^32+xlo := x - 2^32*yhi^2
1076   // (also xhi := xhi - yhi^2, das ist >=0, <=2*yhi).
1077   // Die Schätzung für die zweite Ziffer
1078   //     ylo' := min(2^16-1,floor((xhi*2^32+xlo)/(2*2^16*yhi)))
1079   // erfüllt ylo'-1 <= ylo <= ylo', ist also um höchstens 1 zu groß.
1080   // (Beweis: Rechte Ungleichung klar, da  ylo < 2^16  und
1081   //   xhi*2^32+xlo >= 2*2^16*yhi*ylo + ylo^2 >= 2*2^16*yhi*ylo
1082   //   ==> (xhi*2^32+xlo)/(2*2^16*yhi) >= ylo  gelten muß.
1083   //   Linke Ungleichung: Falls floor(...)>=2^16, ist
1084   //   xhi*2^32+xlo >= 2*2^16*2^16*yhi >= 2*2^16*yhi*(2^16-1) + 2^32
1085   //                >= 2*2^16*yhi*(2^16-1) + (2^16-1)^2
1086   //   und xhi*2^32+xlo < 2*2^16*2^16*yhi + (2^16)^2, also
1087   //   ylo = 2^16-1 = ylo'.
1088   //   Sonst ist ylo' = floor((xhi*2^32+xlo)/(2*2^16*yhi)), also
1089   //   xhi*2^32+xlo >= 2*2^16*yhi*ylo' >= 2*2^16*yhi*(ylo'-1) + 2^32
1090   //                >= 2*2^16*yhi*(ylo'-1) + (ylo'-1)^2,
1091   //   also ylo >= ylo'-1 nach Definition von ylo.)
1092   #define isqrt_64_32(xhi,xlo,y_zuweisung,sqrtp_zuweisung)  \
1093     { var uint32 _xhi = (xhi);                                          \
1094       var uint32 _xlo = (xlo);                                          \
1095       var uint16 _yhi;                                                  \
1096       var uint16 _ylo;                                                  \
1097       /* erste Ziffer berechnen: */                                     \
1098       isqrt_32_16(_xhi,_yhi=,); /* yhi := isqrt(xhi) */                 \
1099       _xhi -= mulu16(_yhi,_yhi); /* jetzt 0 <= xhi <= 2*yhi */          \
1100       /* x = 2^32*yhi^2 + 2^32*xhi + xlo */                             \
1101       /* Schätzung für die zweite Ziffer berechnen: */                        \
1102       /* ylo := min(2^16-1,floor((xhi*2^32+xlo)/(2*2^16*yhi))) bilden: */\
1103      {var uint32 _z = (_xhi << 15) | (_xlo >> 17); /* < 2^15*(2*yhi+1) */\
1104       var uint32 _r = highlow32_0(_yhi);                                \
1105       if (_z >= _r)                                                     \
1106         { _ylo = bit(16)-1; _r = _z - _r + (uint32)_yhi; }              \
1107         else                                                            \
1108         { divu_3216_1616(_z,_yhi, _ylo=,_r=); }                         \
1109       /* x = 2^32*yhi^2 + 2*2^16*yhi*ylo + 2^17*r + (xlo mod 2^17), */  \
1110       /* 0 <= r < yhi + 2^15 */                                         \
1111       _xlo = (_r << 17) | (_xlo & (bit(17)-1));                         \
1112       /* x = 2^32*yhi^2 + 2*2^16*yhi*ylo + 2^32*floor(r/2^15) + xlo */  \
1113       _z = mulu16(_ylo,_ylo); /* z = ylo^2 */                           \
1114       /* Versuche vom Rest 2^32*floor(r/2^15) + xlo  z zu subtrahieren. */\
1115       /* Falls Rest >= z (d.h. r>=2^15 oder xlo>=z), ist ylo fertig, */ \
1116       /* und es gilt x=y^2 genau dann, wenn r<2^15 und xlo=z. */        \
1117       /* Sonst (d.h. r<2^15 und xlo<z), muß man ylo erniedrigen. Dazu */\
1118       /* setzt man  ylo := ylo-1, z := z-(2*ylo+1), */                  \
1119       /* Rest := Rest + 2^17*yhi = xlo + 2^17*yhi >= 2^32 > z, also x>y^2. */\
1120       if (_r < bit(15))                                                 \
1121         { if (_xlo < _z)                                                \
1122             { _ylo -= 1; sqrtp_zuweisung FALSE; }                       \
1123             else                                                        \
1124             { sqrtp_zuweisung (_xlo == _z); }                           \
1125         }                                                               \
1126         else                                                            \
1127         { sqrtp_zuweisung FALSE; }                                      \
1128       y_zuweisung highlow32(_yhi,_ylo);                                 \
1129     }}
1130 #endif
1131
1132 #ifdef HAVE_FAST_LONGLONG
1133
1134 // Zieht die Ganzzahl-Wurzel aus einer 128-Bit-Zahl und
1135 // liefert eine 64-Bit-Wurzel und einen Rest.
1136 // isqrt_128_64(xhi,xlo,y=,sqrtp=);
1137 // > uint64 xhi,xlo: Radikand x = 2^64*xhi+xlo, >= 2^126, < 2^128
1138 // < uint64 y: floor(sqrt(x)), >= 2^63, < 2^64
1139 // < boolean sqrtp: /=0, falls x=y^2
1140   // Methode:
1141   // Wie bei UDS_sqrt mit n=2.
1142   // y = 2^32*yhi + ylo ansetzen.
1143   // Dann muß
1144   //   yhi = floor(y/2^32) = floor(floor(sqrt(x))/2^32)
1145   //       = floor(sqrt(x)/2^32) = floor(sqrt(x/2^64)) = isqrt(xhi)
1146   // sein. Es folgt yhi >= 2^31.
1147   // Danach sucht man das größte ylo >=0 mit
1148   // x - 2^64*yhi^2 >= 2*2^32*yhi*ylo + ylo^2.
1149   // Dazu setzen wir  xhi*2^64+xlo := x - 2^64*yhi^2
1150   // (also xhi := xhi - yhi^2, das ist >=0, <=2*yhi).
1151   // Die Schätzung für die zweite Ziffer
1152   //     ylo' := min(2^32-1,floor((xhi*2^64+xlo)/(2*2^32*yhi)))
1153   // erfüllt ylo'-1 <= ylo <= ylo', ist also um höchstens 1 zu groß.
1154   // (Beweis: Rechte Ungleichung klar, da  ylo < 2^32  und
1155   //   xhi*2^64+xlo >= 2*2^32*yhi*ylo + ylo^2 >= 2*2^32*yhi*ylo
1156   //   ==> (xhi*2^64+xlo)/(2*2^32*yhi) >= ylo  gelten muß.
1157   //   Linke Ungleichung: Falls floor(...)>=2^32, ist
1158   //   xhi*2^64+xlo >= 2*2^32*2^32*yhi >= 2*2^32*yhi*(2^32-1) + 2^64
1159   //                >= 2*2^32*yhi*(2^32-1) + (2^32-1)^2
1160   //   und xhi*2^64+xlo < 2*2^32*2^32*yhi + (2^32)^2, also
1161   //   ylo = 2^32-1 = ylo'.
1162   //   Sonst ist ylo' = floor((xhi*2^64+xlo)/(2*2^32*yhi)), also
1163   //   xhi*2^64+xlo >= 2*2^32*yhi*ylo' >= 2*2^32*yhi*(ylo'-1) + 2^64
1164   //                >= 2*2^32*yhi*(ylo'-1) + (ylo'-1)^2,
1165   //   also ylo >= ylo'-1 nach Definition von ylo.)
1166   #define isqrt_128_64(x_hi,x_lo,y_zuweisung,sqrtp_zuweisung)  \
1167     { var uint64 xhi = (x_hi);                                          \
1168       var uint64 xlo = (x_lo);                                          \
1169       var uint32 yhi;                                                   \
1170       var uint32 ylo;                                                   \
1171       /* erste Ziffer berechnen: */                                     \
1172       isqrt_64_32(high32(xhi),low32(xhi),yhi=,); /* yhi := isqrt(xhi) */\
1173       xhi -= mulu32_w(yhi,yhi); /* jetzt 0 <= xhi <= 2*yhi */           \
1174       /* x = 2^64*yhi^2 + 2^64*xhi + xlo */                             \
1175       /* Schätzung für die zweite Ziffer berechnen: */                        \
1176       /* ylo := min(2^32-1,floor((xhi*2^64+xlo)/(2*2^32*yhi))) bilden: */\
1177      {var uint64 z = (xhi << 31) | (xlo >> 33); /* < 2^31*(2*yhi+1) */  \
1178       var uint64 r = highlow64_0(yhi);                                  \
1179       if (z >= r)                                                       \
1180         { ylo = bit(32)-1; r = z - r + (uint64)yhi; }                   \
1181         else                                                            \
1182         { divu_6432_3232_w(z,yhi, ylo=,r=); }                           \
1183       /* x = 2^64*yhi^2 + 2*2^32*yhi*ylo + 2^33*r + (xlo mod 2^33), */  \
1184       /* 0 <= r < yhi + 2^31 */                                         \
1185       xlo = (r << 33) | (xlo & (bit(33)-1));                            \
1186       /* x = 2^64*yhi^2 + 2*2^32*yhi*ylo + 2^64*floor(r/2^31) + xlo */  \
1187       z = mulu32_w(ylo,ylo); /* z = ylo^2 */                            \
1188       /* Versuche vom Rest 2^64*floor(r/2^31) + xlo  z zu subtrahieren. */\
1189       /* Falls Rest >= z (d.h. r>=2^31 oder xlo>=z), ist ylo fertig, */ \
1190       /* und es gilt x=y^2 genau dann, wenn r<2^31 und xlo=z. */        \
1191       /* Sonst (d.h. r<2^31 und xlo<z), muß man ylo erniedrigen. Dazu */\
1192       /* setzt man  ylo := ylo-1, z := z-(2*ylo+1), */                  \
1193       /* Rest := Rest + 2^33*yhi = xlo + 2^33*yhi >= 2^64 > z, also x>y^2. */\
1194       if (r < bit(31))                                                  \
1195         { if (xlo < z)                                                  \
1196             { ylo -= 1; sqrtp_zuweisung FALSE; }                        \
1197             else                                                        \
1198             { sqrtp_zuweisung (xlo == z); }                             \
1199         }                                                               \
1200         else                                                            \
1201         { sqrtp_zuweisung FALSE; }                                      \
1202       y_zuweisung highlow64(yhi,ylo);                                   \
1203     }}
1204
1205 #endif /* HAVE_FAST_LONGLONG */
1206
1207 // Zieht die Ganzzahl-Wurzel aus einer 32-Bit-Zahl und
1208 // liefert eine 16-Bit-Wurzel.
1209 // isqrt(x)
1210 // > uintL x : Radikand, >=0, <2^32
1211 // < uintL ergebnis : Wurzel, >=0, <2^16
1212   extern uintL isqrt (uintL x);
1213
1214 #ifdef HAVE_LONGLONG
1215 // Extracts integer root of a 64-bit number and returns a 32-bit number.
1216 // isqrt(x)
1217 // > uintQ x : radicand, >=0, <2^64
1218 // < uintL result : square root, >=0, <2^32
1219   extern uintL isqrt (uintQ x);
1220 #endif
1221
1222 // Sorry for this. We need an isqrt function taking uintC arguments but we
1223 // cannot use overloading since this would lead to ambiguities with any of the
1224 // two signatures above.
1225   inline uintL isqrtC (uintC x)
1226   {
1227 #if (intCsize==32)
1228       return isqrt((uintL)x);
1229 #else
1230       return isqrt((uintQ)x);
1231 #endif
1232   }
1233
1234
1235 // Zieht die Ganzzahl-Wurzel aus einer 64-Bit-Zahl und
1236 // liefert eine 32-Bit-Wurzel.
1237 // isqrt(x1,x0)
1238 // > uintL2 x = x1*2^32+x0 : Radikand, >=0, <2^64
1239 // < uintL ergebnis : Wurzel, >=0, <2^32
1240   extern uintL isqrt (uintL x1, uintL x0);
1241
1242
1243 // Bits einer 8-Bit-Zahl zählen:
1244 // integerlength8(digit,size=);
1245 // setzt size auf die höchste in digit vorkommende Bitnummer.
1246 // > digit: ein uint8 >0
1247 // < size: >0, <=8, mit 2^(size-1) <= digit < 2^size
1248 #if defined(__GNUC__) && defined(__m68k__) && !defined(NO_ASM)
1249   #define integerlength8(digit,size_zuweisung)  \
1250     { var uintL _zero_counter; /* zählt die führenden Nullbits in digit            */\
1251       __asm__("bfffo %1{#0:#8},%0" : "=d" (_zero_counter) : "dm" ((uint8)(digit)) ); \
1252       size_zuweisung (8-_zero_counter);                                              \
1253     }
1254 #elif defined(__sparc__) && !defined(__sparc64__)
1255   #define integerlength8(digit,size_zuweisung)  \
1256     integerlength32((uint32)(digit),size_zuweisung) // siehe unten
1257 #elif defined(__GNUC__) && defined(__i386__) && !defined(NO_ASM)
1258   #define integerlength8(digit,size_zuweisung)  \
1259     integerlength16((uint16)(digit),size_zuweisung)
1260 #else
1261   #define integerlength8(digit,size_zuweisung)  \
1262     { var uintC _bitsize = 1;                                   \
1263       var uintL _x8 = (uint8)(digit);                           \
1264       /* _x8 hat höchstens 8 Bits.                             */\
1265       if (_x8 >= bit(4)) { _x8 = _x8>>4; _bitsize += 4; }               \
1266       /* _x8 hat höchstens 4 Bits.                             */\
1267       if (_x8 >= bit(2)) { _x8 = _x8>>2; _bitsize += 2; }               \
1268       /* _x8 hat höchstens 2 Bits.                             */\
1269       if (_x8 >= bit(1)) { /* _x8 = _x8>>1; */ _bitsize += 1; } \
1270       /* _x8 hat höchstens 1 Bit. Dieses Bit muß gesetzt sein. */\
1271       size_zuweisung _bitsize;                                  \
1272     }
1273 #endif
1274
1275 // Bits einer 16-Bit-Zahl zählen:
1276 // integerlength16(digit,size=);
1277 // setzt size auf die höchste in digit vorkommende Bitnummer.
1278 // > digit: ein uint16 >0
1279 // < size: >0, <=16, mit 2^(size-1) <= digit < 2^size
1280 #if defined(__GNUC__) && defined(__m68k__) && !defined(NO_ASM)
1281   #define integerlength16(digit,size_zuweisung)  \
1282     { var uintL _zero_counter; /* zählt die führenden Nullbits in digit              */\
1283       __asm__("bfffo %1{#0:#16},%0" : "=d" (_zero_counter) : "dm" ((uint16)(digit)) ); \
1284       size_zuweisung (16-_zero_counter);                                               \
1285     }
1286 #elif defined(__sparc__) && !defined(__sparc64__)
1287   #define integerlength16(digit,size_zuweisung)  \
1288     integerlength32((uint32)(digit),size_zuweisung) // siehe unten
1289 #elif defined(__GNUC__) && defined(__i386__) && !defined(NO_ASM)
1290   #define integerlength16(digit,size_zuweisung)  \
1291     { var uintW _one_position; /* Position der führenden 1                 */\
1292       __asm__("bsrw %1,%0" : "=r" (_one_position) : "r" ((uint16)(digit)) ); \
1293       size_zuweisung (1+_one_position);                                      \
1294     }
1295 // Die weiteren kommen von gcc/longlong.h :
1296 #elif defined(__GNUC__) && defined(__ibm032__) && !defined(NO_ASM) // RT/ROMP
1297   #define integerlength16(digit,size_zuweisung)  \
1298     { var uintL _zero_counter; /* zählt die führenden Nullbits in digit   */\
1299       __asm__("clz %0,%1" : "=r" (_zero_counter) : "r" ((uint32)(digit)) ); \
1300       size_zuweisung (16-_zero_counter);                                    \
1301     }
1302 #else
1303   #define integerlength16(digit,size_zuweisung)  \
1304     { var uintC _bitsize = 1;                                           \
1305       var uintWL _x16 = (uint16)(digit);                                        \
1306       /* _x16 hat höchstens 16 Bits.                                   */\
1307       if (_x16 >= bit(8)) { _x16 = _x16>>8; _bitsize += 8; }            \
1308       /* _x16 hat höchstens 8 Bits.                                    */\
1309       if (_x16 >= bit(4)) { _x16 = _x16>>4; _bitsize += 4; }            \
1310       /* _x16 hat höchstens 4 Bits.                                    */\
1311       if (_x16 >= bit(2)) { _x16 = _x16>>2; _bitsize += 2; }            \
1312       /* _x16 hat höchstens 2 Bits.                                    */\
1313       if (_x16 >= bit(1)) { /* _x16 = _x16>>1; */ _bitsize += 1; }              \
1314       /* _x16 hat höchstens 1 Bit. Dieses Bit muß gesetzt sein.        */\
1315       size_zuweisung _bitsize;                                          \
1316     }
1317 #endif
1318
1319 // Bits einer 32-Bit-Zahl zählen:
1320 // integerlength32(digit,size=);
1321 // setzt size auf die höchste in digit vorkommende Bitnummer.
1322 // > digit: ein uint32 >0
1323 // < size: >0, <=32, mit 2^(size-1) <= digit < 2^size
1324 #if defined(__GNUC__) && defined(__m68k__) && !defined(NO_ASM)
1325   #define integerlength32(digit,size_zuweisung)  \
1326     { var uintL _zero_counter; /* zählt die führenden Nullbits in digit              */\
1327       __asm__("bfffo %1{#0:#32},%0" : "=d" (_zero_counter) : "dm" ((uint32)(digit)) ); \
1328       size_zuweisung (32-_zero_counter);                                               \
1329     }
1330 #elif defined(__sparc__) && !defined(__sparc64__) && defined(FAST_DOUBLE)
1331   #define integerlength32(digit,size_zuweisung)  \
1332     {var union { double f; uint32 i[2]; } __fi;                         \
1333      const int df_mant_len = 52;  /* mantissa bits (excl. hidden bit) */\
1334      const int df_exp_mid = 1022; /* exponent bias */                   \
1335      /* Bilde 2^52 + digit:                                           */\
1336      __fi.i[0] = (uint32)(df_mant_len+1+df_exp_mid) << (df_mant_len-32); /* Vorzeichen 0, Exponent 53 */\
1337      __fi.i[1] = (digit); /* untere 32 Bits setzen (benutzt CL_CPU_BIG_ENDIAN_P !) */\
1338      /* subtrahiere 2^52:                                             */\
1339      __fi.f = __fi.f - (double)(4503599627370496.0L);                   \
1340      /* Hole davon den Exponenten:                                    */\
1341      size_zuweisung ((__fi.i[0] >> (df_mant_len-32)) - df_exp_mid);     \
1342     }
1343 #elif defined(__GNUC__) && defined(__i386__) && !defined(NO_ASM)
1344   #define integerlength32(digit,size_zuweisung)  \
1345     { var uintL _one_position; /* Position der führenden 1                  */\
1346       __asm__("bsrl %1,%0" : "=r" (_one_position) : "rm" ((uint32)(digit)) ); \
1347       size_zuweisung (1+_one_position);                                       \
1348     }
1349 #elif defined(__hppa__) && !defined(NO_ASM)
1350   #define integerlength32(digit,size_zuweisung)  \
1351     size_zuweisung length32(digit);
1352   extern "C" uintL length32 (uintL digit); // extern in Assembler
1353 // Die weiteren kommen von gcc/longlong.h :
1354 #elif defined(__GNUC__) && (defined(__a29k__) || defined(___AM29K__)) && !defined(NO_ASM)
1355   #define integerlength32(digit,size_zuweisung)  \
1356     { var uintL _zero_counter; /* zählt die führenden Nullbits in digit   */\
1357       __asm__("clz %0,%1" : "=r" (_zero_counter) : "r" ((uint32)(digit)) ); \
1358       size_zuweisung (32-_zero_counter);                                    \
1359     }
1360 #elif defined(__GNUC__) && defined(__gmicro__) && !defined(NO_ASM)
1361   #define integerlength32(digit,size_zuweisung)  \
1362     { var uintL _zero_counter; /* zählt die führenden Nullbits in digit      */\
1363       __asm__("bsch/1 %1,%0" : "=g" (_zero_counter) : "g" ((uint32)(digit)) ); \
1364       size_zuweisung (32-_zero_counter);                                       \
1365     }
1366 #elif defined(__GNUC__) && defined(__rs6000__) && !defined(NO_ASM)
1367  #ifdef _AIX
1368   // old assembler syntax
1369   #define integerlength32(digit,size_zuweisung)  \
1370     { var uintL _zero_counter; /* zählt die führenden Nullbits in digit     */\
1371       __asm__("cntlz %0,%1" : "=r" (_zero_counter) : "r" ((uint32)(digit)) ); \
1372       size_zuweisung (32-_zero_counter);                                      \
1373     }
1374  #else
1375   // new assembler syntax
1376   #define integerlength32(digit,size_zuweisung)  \
1377     { var uintL _zero_counter; /* zählt die führenden Nullbits in digit      */\
1378       __asm__("cntlzw %0,%1" : "=r" (_zero_counter) : "r" ((uint32)(digit)) ); \
1379       size_zuweisung (32-_zero_counter);                                       \
1380     }
1381  #endif
1382 #elif defined(__GNUC__) && defined(__m88k__) && !defined(NO_ASM)
1383   #define integerlength32(digit,size_zuweisung)  \
1384     { var uintL _one_position; /* Position der führenden 1                */\
1385       __asm__("ff1 %0,%1" : "=r" (_one_position) : "r" ((uint32)(digit)) ); \
1386       size_zuweisung (1+_one_position);                                     \
1387     }
1388 #elif defined(__GNUC__) && defined(__ibm032__) && !defined(NO_ASM) // RT/ROMP
1389   #define integerlength32(digit,size_zuweisung)  \
1390     { var uintL _x32 = (uint32)(digit);                         \
1391       if (_x32 >= bit(16))                                      \
1392         { integerlength16(_x32>>16,size_zuweisung 16 + ); }     \
1393         else                                                    \
1394         { integerlength16(_x32,size_zuweisung); }               \
1395     }
1396 #else
1397   #define integerlength32(digit,size_zuweisung)  \
1398     { var uintC _bitsize = 1;                                           \
1399       var uintL _x32 = (uint32)(digit);                                 \
1400       /* _x32 hat höchstens 32 Bits.                                   */\
1401       if (_x32 >= bit(16)) { _x32 = _x32>>16; _bitsize += 16; }         \
1402       /* _x32 hat höchstens 16 Bits.                                   */\
1403       if (_x32 >= bit(8)) { _x32 = _x32>>8; _bitsize += 8; }            \
1404       /* _x32 hat höchstens 8 Bits.                                    */\
1405       if (_x32 >= bit(4)) { _x32 = _x32>>4; _bitsize += 4; }            \
1406       /* _x32 hat höchstens 4 Bits.                                    */\
1407       if (_x32 >= bit(2)) { _x32 = _x32>>2; _bitsize += 2; }            \
1408       /* _x32 hat höchstens 2 Bits.                                    */\
1409       if (_x32 >= bit(1)) { /* _x32 = _x32>>1; */ _bitsize += 1; }      \
1410       /* _x32 hat höchstens 1 Bit. Dieses Bit muß gesetzt sein.        */\
1411       size_zuweisung _bitsize;                                          \
1412     }
1413   #define GENERIC_INTEGERLENGTH32
1414 #endif
1415
1416 // Bits einer 64-Bit-Zahl zählen:
1417 // integerlength64(digit,size=);
1418 // setzt size auf die höchste in digit vorkommende Bitnummer.
1419 // > digit: ein uint64 >0
1420 // < size: >0, <=64, mit 2^(size-1) <= digit < 2^size
1421 #ifdef GENERIC_INTEGERLENGTH32
1422   #define integerlength64(digit,size_zuweisung)  \
1423     { var uintC _bitsize = 1;                                           \
1424       var uint64 _x64 = (uint64)(digit);                                \
1425       /* _x64 hat höchstens 64 Bits.                                   */\
1426       if (_x64 >= bit(32)) { _x64 = _x64>>32; _bitsize += 32; }         \
1427       /* _x64 hat höchstens 32 Bits.                                   */\
1428       if (_x64 >= bit(16)) { _x64 = _x64>>16; _bitsize += 16; }         \
1429       /* _x64 hat höchstens 16 Bits.                                   */\
1430       if (_x64 >= bit(8)) { _x64 = _x64>>8; _bitsize += 8; }            \
1431       /* _x64 hat höchstens 8 Bits.                                    */\
1432       if (_x64 >= bit(4)) { _x64 = _x64>>4; _bitsize += 4; }            \
1433       /* _x64 hat höchstens 4 Bits.                                    */\
1434       if (_x64 >= bit(2)) { _x64 = _x64>>2; _bitsize += 2; }            \
1435       /* _x64 hat höchstens 2 Bits.                                    */\
1436       if (_x64 >= bit(1)) { /* _x64 = _x64>>1; */ _bitsize += 1; }      \
1437       /* _x64 hat höchstens 1 Bit. Dieses Bit muß gesetzt sein.        */\
1438       size_zuweisung _bitsize;                                          \
1439     }
1440 #else
1441   #define integerlength64(digit,size_zuweisung)  \
1442     { var uint64 _x64 = (digit);                                        \
1443       var uintC _bitsize64 = 0;                                         \
1444       var uint32 _x32_from_integerlength64;                             \
1445       if (_x64 >= (1ULL << 32)) {                                       \
1446         _x32_from_integerlength64 = _x64>>32; _bitsize64 += 32;         \
1447       } else {                                                          \
1448         _x32_from_integerlength64 = _x64;                               \
1449       }                                                                 \
1450       integerlength32(_x32_from_integerlength64, size_zuweisung _bitsize64 + ); \
1451     }
1452 #endif
1453
1454 // Bits einer uintC-Zahl zählen:
1455 // integerlengthC(digit,size=);
1456 // setzt size auf die höchste in digit vorkommende Bitnummer.
1457 // > digit: ein uintC >0
1458 // < size: >0, <=intCsize, mit 2^(size-1) <= digit < 2^size
1459   #if (intCsize==32)
1460     #define integerlengthC  integerlength32
1461   #endif
1462   #if (intCsize==64)
1463     #define integerlengthC  integerlength64
1464   #endif
1465
1466 // Hintere Nullbits eines 32-Bit-Wortes zählen:
1467 // ord2_32(digit,count=);
1468 // setzt size auf die kleinste in digit vorkommende Bitnummer.
1469 // > digit: ein uint32 >0
1470 // < count: >=0, <32, mit 2^count | digit, digit/2^count ungerade
1471   #if defined(__GNUC__) && defined(__i386__) && !defined(NO_ASM)
1472     #define ord2_32(digit,count_zuweisung)  \
1473       { var uintL _one_position; /* Position der letzten 1                    */\
1474         __asm__("bsfl %1,%0" : "=r" (_one_position) : "rm" ((uint32)(digit)) ); \
1475         count_zuweisung _one_position;                                          \
1476       }
1477     #define FAST_ORD2
1478   #elif defined(__sparc__) && !defined(__sparc64__)
1479     #define ord2_32(digit,count_zuweisung)  \
1480     { var uint32 n = (digit);                                             \
1481       n = n | -n;                                                         \
1482       n = (n<<4) + n;                                                     \
1483       n = (n<<6) + n;                                                     \
1484       n = n - (n<<16); /* or  n = n ^ (n<<16);  or  n = n &~ (n<<16);  */ \
1485       /* static const char ord2_tab [64] = {-1,0,1,12,2,6,-1,13,3,-1,7,-1,-1,-1,-1,14,10,4,-1,-1,8,-1,-1,25,-1,-1,-1,-1,-1,21,27,15,31,11,5,-1,-1,-1,-1,-1,9,-1,-1,24,-1,-1,20,26,30,-1,-1,-1,-1,23,-1,19,29,-1,22,18,28,17,16,-1}; */ \
1486       /* count_zuweisung ord2_tab[n>>26];                              */ \
1487       count_zuweisung "\377\000\001\014\002\006\377\015\003\377\007\377\377\377\377\016\012\004\377\377\010\377\377\031\377\377\377\377\377\025\033\017\037\013\005\377\377\377\377\377\011\377\377\030\377\377\024\032\036\377\377\377\377\027\377\023\035\377\026\022\034\021\020"[n>>26]; \
1488     }
1489     #define FAST_ORD2
1490   #else
1491     // Sei n = ord2(x). Dann ist logxor(x,x-1) = 2^n + (2^n-1) = 2^(n+1)-1.
1492     // Also  (ord2 x) = (1- (integer-length (logxor x (1- x)))) .
1493     #define ord2_32(digit,count_zuweisung)  \
1494       { var uint32 _digit = (digit) ^ ((digit) - 1);    \
1495         integerlength32(_digit,count_zuweisung -1 + )   \
1496       }
1497   #endif
1498
1499 // Hintere Nullbits eines 64-Bit-Wortes zählen:
1500 // ord2_64(digit,count=);
1501 // setzt size auf die kleinste in digit vorkommende Bitnummer.
1502 // > digit: ein uint64 >0
1503 // < count: >=0, <64, mit 2^count | digit, digit/2^count ungerade
1504   // Sei n = ord2(x). Dann ist logxor(x,x-1) = 2^n + (2^n-1) = 2^(n+1)-1.
1505   // Also  (ord2 x) = (1- (integer-length (logxor x (1- x)))) .
1506   #define ord2_64(digit,count_zuweisung)  \
1507     { var uint64 _digit = (digit) ^ ((digit) - 1);      \
1508       integerlength64(_digit,count_zuweisung -1 + )     \
1509     }
1510
1511
1512 // Bits eines Wortes zählen.
1513 // logcount_NN();
1514 // > xNN: ein uintNN
1515 // < xNN: Anzahl der darin gesetzten Bits
1516   // Bits von x8 zählen: (Input x8, Output x8)
1517   #define logcount_8()  \
1518     ( /* x8 besteht aus 8 1-Bit-Zählern (0,1).        */\
1519       x8 = (x8 & 0x55U) + ((x8 & 0xAAU) >> 1),          \
1520       /* x8 besteht aus 4 2-Bit-Zählern (0,1,2).      */\
1521       x8 = (x8 & 0x33U) + ((x8 & 0xCCU) >> 2),          \
1522       /* x8 besteht aus 2 4-Bit-Zählern (0,1,2,3,4).  */\
1523       x8 = (x8 & 0x0FU) + (x8 >> 4)                     \
1524       /* x8 besteht aus 1 8-Bit-Zähler (0,...,8).     */\
1525     )
1526   // Bits von x16 zählen: (Input x16, Output x16)
1527   #define logcount_16()  \
1528     ( /* x16 besteht aus 16 1-Bit-Zählern (0,1).      */\
1529       x16 = (x16 & 0x5555U) + ((x16 & 0xAAAAU) >> 1),   \
1530       /* x16 besteht aus 8 2-Bit-Zählern (0,1,2).     */\
1531       x16 = (x16 & 0x3333U) + ((x16 & 0xCCCCU) >> 2),   \
1532       /* x16 besteht aus 4 4-Bit-Zählern (0,1,2,3,4). */\
1533       x16 = (x16 & 0x0F0FU) + ((x16 & 0xF0F0U) >> 4),   \
1534       /* x16 besteht aus 2 8-Bit-Zählern (0,...,8).   */\
1535       x16 = (x16 & 0x00FFU) + (x16 >> 8)                \
1536       /* x16 besteht aus 1 16-Bit-Zähler (0,...,16).  */\
1537     )
1538   // Bits von x32 zählen: (Input x32, Output x32)
1539   #define logcount_32()  \
1540     ( /* x32 besteht aus 32 1-Bit-Zählern (0,1).              */\
1541       x32 = (x32 & 0x55555555UL) + ((x32 & 0xAAAAAAAAUL) >> 1), \
1542       /* x32 besteht aus 16 2-Bit-Zählern (0,1,2).            */\
1543       x32 = (x32 & 0x33333333UL) + ((x32 & 0xCCCCCCCCUL) >> 2), \
1544       /* x32 besteht aus 8 4-Bit-Zählern (0,1,2,3,4).         */\
1545       x32 = high16(x32)+low16(x32),                             \
1546       /* x32 besteht aus 4 4-Bit-Zählern (0,...,8).           */\
1547       x32 = (x32 & 0x0F0FU) + ((x32 & 0xF0F0U) >> 4),           \
1548       /* x32 besteht aus 2 8-Bit-Zählern (0,...,16).          */\
1549       x32 = (x32 & 0x00FFU) + (x32 >> 8)                        \
1550       /* x32 besteht aus 1 16-Bit-Zähler (0,...,32).          */\
1551     )
1552   // Bits von x64 zählen: (Input x64, Output x64)
1553   #define logcount_64()  \
1554     ( /* x64 besteht aus 64 1-Bit-Zählern (0,1).                             */\
1555       x64 = (x64 & 0x5555555555555555ULL) + ((x64 & 0xAAAAAAAAAAAAAAAAULL) >> 1),\
1556       /* x64 besteht aus 32 2-Bit-Zählern (0,1,2).                           */\
1557       x64 = (x64 & 0x3333333333333333ULL) + ((x64 & 0xCCCCCCCCCCCCCCCCULL) >> 2),\
1558       /* x64 besteht aus 16 4-Bit-Zählern (0,1,2,3,4).                       */\
1559       x64 = (uint32)(x64 + (x64 >> 32)),                                       \
1560       /* x64 besteht aus 8 4-Bit-Zählern (0,...,8).                          */\
1561       x64 = (x64 & 0x0F0F0F0FUL) + ((x64 & 0xF0F0F0F0UL) >> 4),                \
1562       /* x64 besteht aus 4 8-Bit-Zählern (0,...,16).                         */\
1563       x64 = (x64 & 0x00FF00FFU) + ((x64 & 0xFF00FF00U) >> 8),                  \
1564       /* x64 besteht aus 2 16-Bit-Zählern (0,...,32).                        */\
1565       x64 = (x64 & 0x0000FFFFU) + (x64 >> 16)                                  \
1566       /* x64 besteht aus 1 16-Bit-Zähler (0,...,64).                         */\
1567     )
1568
1569 }  // namespace cln
1570
1571 #endif /* _CL_LOW_H */