4 #include "base/cl_sysdep.h"
7 #include "base/cl_low.h"
14 // Zieht die Ganzzahl-Wurzel aus einer 32-Bit-Zahl und
15 // liefert eine 16-Bit-Wurzel.
17 // > uintL x : Radikand, >=0, <2^32
18 // < uintL ergebnis : Wurzel, >=0, <2^16
22 // x=0 -> y=0, fertig.
23 // y := 2^k als Anfangswert, wobei k>0, k<=16 mit 2^(2k-2) <= x < 2^(2k) sei.
24 // y := floor((y + floor(x/y))/2) als nächster Wert,
25 // solange z := floor(x/y) < y, setze y := floor((y+z)/2).
28 // 1. Die Folge der y ist streng monoton fallend.
29 // 2. Stets gilt y >= floor(sqrt(x)) (denn für alle y>0 ist
30 // y + x/y >= 2*sqrt(x) und daher floor((y + floor(x/y))/2) =
31 // floor(y/2 + x/(2*y)) >= floor(sqrt(x)) ).
32 // 3. Am Schluß gilt x >= y^2.
34 if (x==0) return 0; // x=0 -> y=0
35 { var uintC k2; integerlength32(x,k2=); // 2^(k2-1) <= x < 2^k2
36 {var uintC k1 = floor(k2-1,2); // k1 = k-1, k wie oben
39 { var uintL y = (x >> (k1+2)) | bit(k1); // stets 2^(k-1) <= y < 2^k
42 divu_3216_1616(x,y, z=,); // Dividiere x/y (geht, da x/y < 2^(2k)/2^(k-1) = 2^(k+1) <= 2^16)
44 y = floor(z+y,2); // geht, da z+y < 2*y < 2^(k+1) <= 2^16
50 { var uintL x1 = high16(x);
51 var uintL y = (x >> (16+1)) | bit(16-1); // stets 2^(k-1) <= y < 2^k
54 if (x1 >= y) break; // Division x/y ergäbe Überlauf -> z > y
55 divu_3216_1616(x,y, z=,); // Dividiere x/y
66 // As isqrt (uintL) above, but with 64-bit numbers.
67 if (x==0) return 0; // x=0 -> y=0
68 { var uintC k2; integerlength64(x,k2=); // 2^(k2-1) <= x < 2^k2
69 {var uintC k1 = floor(k2-1,2); // k1 = k-1, k as above
72 { var uintL y = (x >> (k1+2)) | bit(k1); // always 2^(k-1) <= y < 2^k
75 divu_6432_3232(high32(x),low32(x),y, z=,); // z := x/y (works, since x/y < 2^(2k)/2^(k-1) = 2^(k+1) <= 2^32)
77 y = floor(z+y,2); // geht, da z+y < 2*y < 2^(k+1) <= 2^32
83 { var uintL x1 = high32(x);
84 var uintL y = (x >> (32+1)) | bit(32-1); // stets 2^(k-1) <= y < 2^k
87 if (x1 >= y) break; // division x/y would overflow -> z > y
88 divu_6432_3232(high32(x),low32(x),y, z=,); // divide x/y