not-implemented.tex - libzahl - big integer library | |
git clone git://git.suckless.org/libzahl | |
Log | |
Files | |
Refs | |
README | |
LICENSE | |
--- | |
not-implemented.tex (19552B) | |
--- | |
1 \chapter{Not implemented} | |
2 \label{chap:Not implemented} | |
3 | |
4 In this chapter we maintain a list of | |
5 features we have chosen not to implement at | |
6 the moment, but may add when libzahl matures, | |
7 but to a separate library that could be made | |
8 to support multiple bignum libraries. Functions | |
9 listed herein will only be implemented if it | |
10 is shown that it would be overwhelmingly | |
11 advantageous. For each feature, a sample | |
12 implementation or a mathematical expression | |
13 on which you can base your implementation is | |
14 included. The sample implementations create | |
15 temporary integer references to simplify the | |
16 examples. You should try to use dedicated | |
17 variables; in case of recursion, a robust | |
18 program should store temporary variables on | |
19 a stack, so they can be cleaned up if | |
20 something happens. | |
21 | |
22 Research problems, like prime factorisation | |
23 and discrete logarithms, do not fit in the | |
24 scope of bignum libraries % Unless they are extraordinarily bloated with… | |
25 and therefore do not fit into libzahl, | |
26 and will not be included in this chapter. | |
27 Operators and functions that grow so | |
28 ridiculously fast that a tiny lookup table | |
29 can be constructed to cover all practical | |
30 input will also not be included in this | |
31 chapter, nor in libzahl. | |
32 | |
33 \vspace{1cm} | |
34 \minitoc | |
35 | |
36 | |
37 \newpage | |
38 \section{Extended greatest common divisor} | |
39 \label{sec:Extended greatest common divisor} | |
40 | |
41 \begin{alltt} | |
42 void | |
43 extgcd(z_t bézout_coeff_1, z_t bézout_coeff_2, z_t gcd | |
44 z_t quotient_1, z_t quotient_2, z_t a, z_t b) | |
45 \{ | |
46 #define old_r gcd | |
47 #define old_s bézout_coeff_1 | |
48 #define old_t bézout_coeff_2 | |
49 #define s quotient_2 | |
50 #define t quotient_1 | |
51 z_t r, q, qs, qt; | |
52 int odd = 0; | |
53 zinit(r), zinit(q), zinit(qs), zinit(qt); | |
54 zset(r, b), zset(old_r, a); | |
55 zseti(s, 0), zseti(old_s, 1); | |
56 zseti(t, 1), zseti(old_t, 0); | |
57 while (!zzero(r)) \{ | |
58 odd ^= 1; | |
59 zdivmod(q, old_r, old_r, r), zswap(old_r, r); | |
60 zmul(qs, q, s), zsub(old_s, old_s, qs); | |
61 zmul(qt, q, t), zsub(old_t, old_t, qt); | |
62 zswap(old_s, s), zswap(old_t, t); | |
63 \} | |
64 odd ? abs(s, s) : abs(t, t); | |
65 zfree(r), zfree(q), zfree(qs), zfree(qt); | |
66 \} | |
67 \end{alltt} | |
68 | |
69 Perhaps you are asking yourself ``wait a minute, | |
70 doesn't the extended Euclidean algorithm only | |
71 have three outputs if you include the greatest | |
72 common divisor, what is this shenanigans?'' | |
73 No\footnote{Well, technically yes, but it calculates | |
74 two values for free in the same ways as division | |
75 calculates the remainder for free.}, it has five | |
76 outputs, most implementations just ignore two of | |
77 them. If this confuses you, or you want to know | |
78 more about this, I refer you to Wikipeida. | |
79 | |
80 | |
81 \newpage | |
82 \section{Least common multiple} | |
83 \label{sec:Least common multiple} | |
84 | |
85 \( \displaystyle{ | |
86 \mbox{lcm}(a, b) = \frac{\lvert a \cdot b \rvert}{\mbox{gcd}(a, b)} | |
87 }\) | |
88 \vspace{1em} | |
89 | |
90 $\mbox{lcm}(a, b)$ is undefined when $a$ or | |
91 $b$ is zero, because division by zero is | |
92 undefined. Note however that $\mbox{gcd}(a, b)$ | |
93 is only zero when both $a$ and $b$ is zero. | |
94 | |
95 \newpage | |
96 \section{Modular multiplicative inverse} | |
97 \label{sec:Modular multiplicative inverse} | |
98 | |
99 \begin{alltt} | |
100 int | |
101 modinv(z_t inv, z_t a, z_t m) | |
102 \{ | |
103 z_t x, _1, _2, _3, gcd, mabs, apos; | |
104 int invertible, aneg = zsignum(a) < 0; | |
105 zinit(x), zinit(_1), zinit(_2), zinit(_3), zinit(gcd); | |
106 *mabs = *m; | |
107 zabs(mabs, mabs); | |
108 if (aneg) \{ | |
109 zinit(apos); | |
110 zset(apos, a); | |
111 if (zcmpmag(apos, mabs)) | |
112 zmod(apos, apos, mabs); | |
113 zadd(apos, apos, mabs); | |
114 \} | |
115 extgcd(inv, _1, _2, _3, gcd, apos, mabs); | |
116 if ((invertible = !zcmpi(gcd, 1))) \{ | |
117 if (zsignum(inv) < 0) | |
118 (zsignum(m) < 0 ? zsub : zadd)(x, x, m); | |
119 zswap(x, inv); | |
120 \} | |
121 if (aneg) | |
122 zfree(apos); | |
123 zfree(x), zfree(_1), zfree(_2), zfree(_3), zfree(gcd); | |
124 return invertible; | |
125 \} | |
126 \end{alltt} | |
127 | |
128 | |
129 \newpage | |
130 \section{Random prime number generation} | |
131 \label{sec:Random prime number generation} | |
132 | |
133 TODO | |
134 | |
135 | |
136 \newpage | |
137 \section{Symbols} | |
138 \label{sec:Symbols} | |
139 | |
140 \subsection{Legendre symbol} | |
141 \label{sec:Legendre symbol} | |
142 | |
143 \( \displaystyle{ | |
144 \left ( \frac{a}{p} \right ) \equiv a^{\frac{p - 1}{2}} ~(\text{Mod}~p… | |
145 \left ( \frac{a}{p} \right ) \in \{-1,~0,~1\},~ | |
146 p \in \textbf{P},~ p > 2 | |
147 }\) | |
148 \vspace{1em} | |
149 | |
150 \noindent | |
151 That is, unless $\displaystyle{a^{\frac{p - 1}{2}} \mod p \le 1}$, | |
152 $\displaystyle{a^{\frac{p - 1}{2}} \mod p = p - 1}$, so | |
153 $\displaystyle{\left ( \frac{a}{p} \right ) = -1}$. | |
154 | |
155 It should be noted that | |
156 \( \displaystyle{ | |
157 \left ( \frac{a}{p} \right ) = | |
158 \left ( \frac{a ~\text{Mod}~ p}{p} \right ), | |
159 }\) | |
160 so a compressed lookup table can be used for small $p$. | |
161 | |
162 | |
163 \subsection{Jacobi symbol} | |
164 \label{sec:Jacobi symbol} | |
165 | |
166 \( \displaystyle{ | |
167 \left ( \frac{a}{n} \right ) = | |
168 \prod_k \left ( \frac{a}{p_k} \right )^{n_k}, | |
169 }\) | |
170 where $\displaystyle{n = \prod_k p_k^{n_k} > 0}$, | |
171 and $p_k \in \textbf{P}$. | |
172 \vspace{1em} | |
173 | |
174 Like the Legendre symbol, the Jacobi symbol is $n$-periodic over $a$. | |
175 If $n$, is prime, the Jacobi symbol is the Legendre symbol, but | |
176 the Jacobi symbol is defined for all odd numbers $n$, where the | |
177 Legendre symbol is defined only for odd primes $n$. | |
178 | |
179 Use the following algorithm to calculate the Jacobi symbol: | |
180 | |
181 \vspace{1em} | |
182 \hspace{-2.8ex} | |
183 \begin{minipage}{\linewidth} | |
184 \begin{algorithmic} | |
185 \STATE $a \gets a \mod n$ | |
186 \STATE $r \gets \mbox{lsb}~ a$ | |
187 \STATE $a \gets a \cdot 2^{-r}$ | |
188 \STATE \(\displaystyle{ | |
189 r \gets \left \lbrace \begin{array}{rl} | |
190 1 & | |
191 \text{if}~ n \equiv 1, 7 ~(\text{Mod}~ 8) | |
192 ~\text{or}~ r \equiv 0 ~(\text{Mod}~ 2) \\ | |
193 -1 & \text{otherwise} \\ | |
194 \end{array} \right . | |
195 }\) | |
196 \IF{$a = 1$} | |
197 \RETURN $r$ | |
198 \ELSIF{$\gcd(a, n) \neq 1$} | |
199 \RETURN $0$ | |
200 \ENDIF | |
201 \STATE $(a, n) = (n, a)$ | |
202 \STATE \textbf{start over} | |
203 \end{algorithmic} | |
204 \end{minipage} | |
205 | |
206 | |
207 \subsection{Kronecker symbol} | |
208 \label{sec:Kronecker symbol} | |
209 | |
210 The Kronecker symbol | |
211 $\displaystyle{\left ( \frac{a}{n} \right )}$ | |
212 is a generalisation of the Jacobi symbol, | |
213 where $n$ can be any integer. For positive | |
214 odd $n$, the Kronecker symbol is equal to | |
215 the Jacobi symbol. For even $n$, the | |
216 Kronecker symbol is $2n$-periodic over $a$, | |
217 the Kronecker symbol is zero for all | |
218 $(a, n)$ with both $a$ and $n$ are even. | |
219 | |
220 \vspace{1em} | |
221 \noindent | |
222 \( \displaystyle{ | |
223 \left ( \frac{a}{2^k \cdot n} \right ) = | |
224 \left ( \frac{a}{n} \right ) \cdot \left ( \frac{a}{2} \right )^k, | |
225 }\) | |
226 where | |
227 \( \displaystyle{ | |
228 \left ( \frac{a}{2} \right ) = | |
229 \left \lbrace \begin{array}{rl} | |
230 1 & \text{if}~ a \equiv 1, 7 ~(\text{Mod}~ 8) \\ | |
231 -1 & \text{if}~ a \equiv 3, 5 ~(\text{Mod}~ 8) \\ | |
232 0 & \text{otherwise} | |
233 \end{array} \right . | |
234 }\) | |
235 | |
236 \vspace{1em} | |
237 \noindent | |
238 \( \displaystyle{ | |
239 \left ( \frac{-a}{n} \right ) = | |
240 \left ( \frac{a}{n} \right ) \cdot \left ( \frac{a}{-1} \right ), | |
241 }\) | |
242 where | |
243 \( \displaystyle{ | |
244 \left ( \frac{a}{-1} \right ) = | |
245 \left \lbrace \begin{array}{rl} | |
246 1 & \text{if}~ a \ge 0 \\ | |
247 -1 & \text{if}~ a < 0 | |
248 \end{array} \right . | |
249 }\) | |
250 \vspace{1em} | |
251 | |
252 \noindent | |
253 However, for $n = 0$, the symbol is defined as | |
254 | |
255 \vspace{1em} | |
256 \noindent | |
257 \( \displaystyle{ | |
258 \left ( \frac{a}{0} \right ) = | |
259 \left \lbrace \begin{array}{rl} | |
260 1 & \text{if}~ a = \pm 1 \\ | |
261 0 & \text{otherwise.} | |
262 \end{array} \right . | |
263 }\) | |
264 | |
265 | |
266 \subsection{Power residue symbol} | |
267 \label{sec:Power residue symbol} | |
268 | |
269 TODO | |
270 | |
271 | |
272 \subsection{Pochhammer \emph{k}-symbol} | |
273 \label{sec:Pochhammer k-symbol} | |
274 | |
275 \( \displaystyle{ | |
276 (x)_{n,k} = \prod_{i = 1}^n (x + (i - 1)k) | |
277 }\) | |
278 | |
279 | |
280 \newpage | |
281 \section{Logarithm} | |
282 \label{sec:Logarithm} | |
283 | |
284 TODO | |
285 | |
286 | |
287 \newpage | |
288 \section{Roots} | |
289 \label{sec:Roots} | |
290 | |
291 TODO | |
292 | |
293 | |
294 \newpage | |
295 \section{Modular roots} | |
296 \label{sec:Modular roots} | |
297 | |
298 TODO % Square: Cipolla's algorithm, Pocklington's algorithm, Tonelli–S… | |
299 | |
300 | |
301 \newpage | |
302 \section{Combinatorial} | |
303 \label{sec:Combinatorial} | |
304 | |
305 \subsection{Factorial} | |
306 \label{sec:Factorial} | |
307 | |
308 \( \displaystyle{ | |
309 n! = \left \lbrace \begin{array}{ll} | |
310 \displaystyle{\prod_{i = 1}^n i} & \textrm{if}~ n \ge 0 \\ | |
311 \textrm{undefined} & \textrm{otherwise} | |
312 \end{array} \right . | |
313 }\) | |
314 \vspace{1em} | |
315 | |
316 This can be implemented much more efficiently | |
317 than using the naïve method, and is a very | |
318 important function for many combinatorial | |
319 applications, therefore it may be implemented | |
320 in the future if the demand is high enough. | |
321 | |
322 An efficient, yet not optimal, implementation | |
323 of factorials that about halves the number of | |
324 required multiplications compared to the naïve | |
325 method can be derived from the observation | |
326 | |
327 \vspace{1em} | |
328 \( \displaystyle{ | |
329 n! = n!! ~ \lfloor n / 2 \rfloor! ~ 2^{\lfloor n / 2 \rfloor} | |
330 }\), $n$ odd. | |
331 \vspace{1em} | |
332 | |
333 \noindent | |
334 The resulting algorithm can be expressed as | |
335 | |
336 \begin{alltt} | |
337 void | |
338 fact(z_t r, uint64_t n) | |
339 \{ | |
340 z_t p, f, two; | |
341 uint64_t *ns, s = 1, i = 1; | |
342 zinit(p), zinit(f), zinit(two); | |
343 zseti(r, 1), zseti(p, 1), zseti(f, n), zseti(two, 2); | |
344 ns = alloca(zbits(f) * sizeof(*ns)); | |
345 while (n > 1) \{ | |
346 if (n & 1) \{ | |
347 ns[i++] = n; | |
348 s += n >>= 1; | |
349 \} else \{ | |
350 zmul(r, r, (zsetu(f, n), f)); | |
351 n -= 1; | |
352 \} | |
353 \} | |
354 for (zseti(f, 1); i-- > 0; zmul(r, r, p);) | |
355 for (n = ns[i]; zcmpu(f, n); zadd(f, f, two)) | |
356 zmul(p, p, f); | |
357 zlsh(r, r, s); | |
358 zfree(two), zfree(f), zfree(p); | |
359 \} | |
360 \end{alltt} | |
361 | |
362 | |
363 \subsection{Subfactorial} | |
364 \label{sec:Subfactorial} | |
365 | |
366 \( \displaystyle{ | |
367 !n = \left \lbrace \begin{array}{ll} | |
368 n(!(n - 1)) + (-1)^n & \textrm{if}~ n > 0 \\ | |
369 1 & \textrm{if}~ n = 0 \\ | |
370 \textrm{undefined} & \textrm{otherwise} | |
371 \end{array} \right . = | |
372 n! \sum_{i = 0}^n \frac{(-1)^i}{i!} | |
373 }\) | |
374 | |
375 | |
376 \subsection{Alternating factorial} | |
377 \label{sec:Alternating factorial} | |
378 | |
379 \( \displaystyle{ | |
380 \mbox{af}(n) = \sum_{i = 1}^n {(-1)^{n - i} i!} | |
381 }\) | |
382 | |
383 | |
384 \subsection{Multifactorial} | |
385 \label{sec:Multifactorial} | |
386 | |
387 \( \displaystyle{ | |
388 n!^{(k)} = \left \lbrace \begin{array}{ll} | |
389 1 & \textrm{if}~ n = 0 \\ | |
390 n & \textrm{if}~ 0 < n \le k \\ | |
391 n((n - k)!^{(k)}) & \textrm{if}~ n > k \\ | |
392 \textrm{undefined} & \textrm{otherwise} | |
393 \end{array} \right . | |
394 }\) | |
395 | |
396 | |
397 \subsection{Quadruple factorial} | |
398 \label{sec:Quadruple factorial} | |
399 | |
400 \( \displaystyle{ | |
401 (4n - 2)!^{(4)} | |
402 }\) | |
403 | |
404 | |
405 \subsection{Superfactorial} | |
406 \label{sec:Superfactorial} | |
407 | |
408 \( \displaystyle{ | |
409 \mbox{sf}(n) = \prod_{k = 1}^n k^{1 + n - k} | |
410 }\), undefined for $n < 0$. | |
411 | |
412 | |
413 \subsection{Hyperfactorial} | |
414 \label{sec:Hyperfactorial} | |
415 | |
416 \( \displaystyle{ | |
417 H(n) = \prod_{k = 1}^n k^k | |
418 }\), undefined for $n < 0$. | |
419 | |
420 | |
421 \subsection{Raising factorial} | |
422 \label{sec:Raising factorial} | |
423 | |
424 \( \displaystyle{ | |
425 x^{(n)} = \frac{(x + n - 1)!}{(x - 1)!} | |
426 }\), undefined for $n < 0$. | |
427 | |
428 | |
429 \subsection{Falling factorial} | |
430 \label{sec:Falling factorial} | |
431 | |
432 \( \displaystyle{ | |
433 (x)_n = \frac{x!}{(x - n)!} | |
434 }\), undefined for $n < 0$. | |
435 | |
436 | |
437 \subsection{Primorial} | |
438 \label{sec:Primorial} | |
439 | |
440 \( \displaystyle{ | |
441 n\# = \prod_{\lbrace i \in \textbf{P} ~:~ i \le n \rbrace} i | |
442 }\) | |
443 \vspace{1em} | |
444 | |
445 \noindent | |
446 \( \displaystyle{ | |
447 p_n\# = \prod_{i \in \textbf{P}_{\pi(n)}} i | |
448 }\) | |
449 | |
450 | |
451 \subsection{Gamma function} | |
452 \label{sec:Gamma function} | |
453 | |
454 $\Gamma(n) = (n - 1)!$, undefined for $n \le 0$. | |
455 | |
456 | |
457 \subsection{K-function} | |
458 \label{sec:K-function} | |
459 | |
460 \( \displaystyle{ | |
461 K(n) = \left \lbrace \begin{array}{ll} | |
462 \displaystyle{\prod_{i = 1}^{n - 1} i^i} & \textrm{if}~ n \ge 0 \\ | |
463 1 & \textrm{if}~ n = -1 \\ | |
464 0 & \textrm{otherwise (result is truncated)} | |
465 \end{array} \right . | |
466 }\) | |
467 | |
468 | |
469 \subsection{Binomial coefficient} | |
470 \label{sec:Binomial coefficient} | |
471 | |
472 \( \displaystyle{ | |
473 \binom{n}{k} = \frac{n!}{k!(n - k)!} | |
474 = \frac{1}{(n - k)!} \prod_{i = k + 1}^n i | |
475 = \frac{1}{k!} \prod_{i = n - k + 1}^n i | |
476 }\) | |
477 | |
478 | |
479 \subsection{Catalan number} | |
480 \label{sec:Catalan number} | |
481 | |
482 \( \displaystyle{ | |
483 C_n = \left . \binom{2n}{n} \middle / (n + 1) \right . | |
484 }\) | |
485 | |
486 | |
487 \subsection{Fuss–Catalan number} | |
488 \label{sec:Fuss-Catalan number} % not en dash | |
489 | |
490 \( \displaystyle{ | |
491 A_m(p, r) = \frac{r}{mp + r} \binom{mp + r}{m} | |
492 }\) | |
493 | |
494 | |
495 \newpage | |
496 \section{Fibonacci numbers} | |
497 \label{sec:Fibonacci numbers} | |
498 | |
499 Fibonacci numbers can be computed efficiently | |
500 using the following algorithm: | |
501 | |
502 \begin{alltt} | |
503 static void | |
504 fib_ll(z_t f, z_t g, z_t n) | |
505 \{ | |
506 z_t a, k; | |
507 int odd; | |
508 if (zcmpi(n, 1) <= 0) \{ | |
509 zseti(f, !zzero(n)); | |
510 zseti(f, zzero(n)); | |
511 return; | |
512 \} | |
513 zinit(a), zinit(k); | |
514 zrsh(k, n, 1); | |
515 if (zodd(n)) \{ | |
516 odd = zodd(k); | |
517 fib_ll(a, g, k); | |
518 zadd(f, a, a); | |
519 zadd(k, f, g); | |
520 zsub(f, f, g); | |
521 zmul(f, f, k); | |
522 zseti(k, odd ? -2 : +2); | |
523 zadd(f, f, k); | |
524 zadd(g, g, g); | |
525 zadd(g, g, a); | |
526 zmul(g, g, a); | |
527 \} else \{ | |
528 fib_ll(g, a, k); | |
529 zadd(f, a, a); | |
530 zadd(f, f, g); | |
531 zmul(f, f, g); | |
532 zsqr(a, a); | |
533 zsqr(g, g); | |
534 zadd(g, a); | |
535 \} | |
536 zfree(k), zfree(a); | |
537 \} | |
538 \end{alltt} | |
539 | |
540 \newpage | |
541 \begin{alltt} | |
542 void | |
543 fib(z_t f, z_t n) | |
544 \{ | |
545 z_t tmp, k; | |
546 zinit(tmp), zinit(k); | |
547 zset(k, n); | |
548 fib_ll(f, tmp, k); | |
549 zfree(k), zfree(tmp); | |
550 \} | |
551 \end{alltt} | |
552 | |
553 \noindent | |
554 This algorithm is based on the rules | |
555 | |
556 \vspace{1em} | |
557 \( \displaystyle{ | |
558 F_{2k + 1} = 4F_k^2 - F_{k - 1}^2 + 2(-1)^k = (2F_k + F_{k-1})(2F_k … | |
559 }\) | |
560 \vspace{1em} | |
561 | |
562 \( \displaystyle{ | |
563 F_{2k} = F_k \cdot (F_k + 2F_{k - 1}) | |
564 }\) | |
565 \vspace{1em} | |
566 | |
567 \( \displaystyle{ | |
568 F_{2k - 1} = F_k^2 + F_{k - 1}^2 | |
569 }\) | |
570 \vspace{1em} | |
571 | |
572 \noindent | |
573 Each call to {\tt fib\_ll} returns $F_n$ and $F_{n - 1}$ | |
574 for any input $n$. $F_{k}$ is only correctly returned | |
575 for $k \ge 0$. $F_n$ and $F_{n - 1}$ is used for | |
576 calculating $F_{2n}$ or $F_{2n + 1}$. The algorithm | |
577 can be sped up with a larger lookup table than one | |
578 covering just the base cases. Alternatively, a naïve | |
579 calculation could be used for sufficiently small input. | |
580 | |
581 | |
582 \newpage | |
583 \section{Lucas numbers} | |
584 \label{sec:Lucas numbers} | |
585 | |
586 Lucas [lyk\textscripta] numbers can be calculated by utilising | |
587 {\tt fib\_ll} from \secref{sec:Fibonacci numbers}: | |
588 | |
589 \begin{alltt} | |
590 void | |
591 lucas(z_t l, z_t n) | |
592 \{ | |
593 z_t k; | |
594 int odd; | |
595 if (zcmp(n, 1) <= 0) \{ | |
596 zset(l, 1 + zzero(n)); | |
597 return; | |
598 \} | |
599 zinit(k); | |
600 zrsh(k, n, 1); | |
601 if (zeven(n)) \{ | |
602 lucas(l, k); | |
603 zsqr(l, l); | |
604 zseti(k, zodd(k) ? +2 : -2); | |
605 zadd(l, k); | |
606 \} else \{ | |
607 odd = zodd(k); | |
608 fib_ll(l, k, k); | |
609 zadd(l, l, l); | |
610 zadd(l, l, k); | |
611 zmul(l, l, k); | |
612 zseti(k, 5); | |
613 zmul(l, l, k); | |
614 zseti(k, odd ? +4 : -4); | |
615 zadd(l, l, k); | |
616 \} | |
617 zfree(k); | |
618 \} | |
619 \end{alltt} | |
620 | |
621 \noindent | |
622 This algorithm is based on the rules | |
623 | |
624 \vspace{1em} | |
625 \( \displaystyle{ | |
626 L_{2k} = L_k^2 - 2(-1)^k | |
627 }\) | |
628 \vspace{1ex} | |
629 | |
630 \( \displaystyle{ | |
631 L_{2k + 1} = 5F_{k - 1} \cdot (2F_k + F_{k - 1}) - 4(-1)^k | |
632 }\) | |
633 \vspace{1em} | |
634 | |
635 \noindent | |
636 Alternatively, the function can be implemented | |
637 trivially using the rule | |
638 | |
639 \vspace{1em} | |
640 \( \displaystyle{ | |
641 L_k = F_k + 2F_{k - 1} | |
642 }\) | |
643 | |
644 | |
645 \newpage | |
646 \section{Bit operation} | |
647 \label{sec:Bit operation unimplemented} | |
648 | |
649 \subsection{Bit scanning} | |
650 \label{sec:Bit scanning} | |
651 | |
652 Scanning for the next set or unset bit can be | |
653 trivially implemented using {\tt zbtest}. A | |
654 more efficient, although not optimally efficient, | |
655 implementation would be | |
656 | |
657 \begin{alltt} | |
658 size_t | |
659 bscan(z_t a, size_t whence, int direction, int value) | |
660 \{ | |
661 size_t ret; | |
662 z_t t; | |
663 zinit(t); | |
664 value ? zset(t, a) : znot(t, a); | |
665 ret = direction < 0 | |
666 ? (ztrunc(t, t, whence + 1), zbits(t) - 1) | |
667 : (zrsh(t, t, whence), zlsb(t) + whence); | |
668 zfree(t); | |
669 return ret; | |
670 \} | |
671 \end{alltt} | |
672 | |
673 | |
674 \subsection{Population count} | |
675 \label{sec:Population count} | |
676 | |
677 The following function can be used to compute | |
678 the population count, the number of set bits, | |
679 in an integer, counting the sign bit: | |
680 | |
681 \begin{alltt} | |
682 size_t | |
683 popcount(z_t a) | |
684 \{ | |
685 size_t i, ret = zsignum(a) < 0; | |
686 for (i = 0; i < a->used; i++) \{ | |
687 ret += __builtin_popcountll(a->chars[i]); | |
688 \} | |
689 return ret; | |
690 \} | |
691 \end{alltt} | |
692 | |
693 \noindent | |
694 It requires a compiler extension; if it's not | |
695 available, there are other ways to computer the | |
696 population count for a word: manually bit-by-bit, | |
697 or with a fully unrolled | |
698 | |
699 \begin{alltt} | |
700 int s; | |
701 for (s = 1; s < 64; s <<= 1) | |
702 w = (w >> s) + w; | |
703 \end{alltt} | |
704 | |
705 | |
706 \subsection{Hamming distance} | |
707 \label{sec:Hamming distance} | |
708 | |
709 A simple way to compute the Hamming distance, | |
710 the number of differing bits between two | |
711 numbers is with the function | |
712 | |
713 \begin{alltt} | |
714 size_t | |
715 hammdist(z_t a, z_t b) | |
716 \{ | |
717 size_t ret; | |
718 z_t t; | |
719 zinit(t); | |
720 zxor(t, a, b); | |
721 ret = popcount(t); | |
722 zfree(t); | |
723 return ret; | |
724 \} | |
725 \end{alltt} | |
726 | |
727 \noindent | |
728 The performance of this function could | |
729 be improved by comparing character by | |
730 character manually using {\tt zxor}. | |
731 | |
732 | |
733 \newpage | |
734 \section{Miscellaneous} | |
735 \label{sec:Miscellaneous} | |
736 | |
737 | |
738 \subsection{Character retrieval} | |
739 \label{sec:Character retrieval} | |
740 | |
741 \begin{alltt} | |
742 uint64_t | |
743 getu(z_t a) | |
744 \{ | |
745 return zzero(a) ? 0 : a->chars[0]; | |
746 \} | |
747 \end{alltt} | |
748 | |
749 \subsection{Fit test} | |
750 \label{sec:Fit test} | |
751 | |
752 Some libraries have functions for testing | |
753 whether a big integer is small enough to | |
754 fit into an intrinsic type. Since libzahl | |
755 does not provide conversion to intrinsic | |
756 types this is irrelevant. But additionally, | |
757 it can be implemented with a single | |
758 one-line macro that does not have any | |
759 side-effects. | |
760 | |
761 \begin{alltt} | |
762 #define fits_in(a, type) (zbits(a) <= 8 * sizeof(type)) | |
763 \textcolor{c}{/* \textrm{Just be sure the type is integral.} */} | |
764 \end{alltt} | |
765 | |
766 | |
767 \subsection{Reference duplication} | |
768 \label{sec:Reference duplication} | |
769 | |
770 This could be useful for creating duplicates | |
771 with modified sign, but only if neither | |
772 {\tt r} nor {\tt a} will be modified whilst | |
773 both are in use. Because it is unsafe, | |
774 fairly simple to create an implementation | |
775 with acceptable performance — {\tt *r = *a}, | |
776 — and probably seldom useful, this has not | |
777 been implemented. | |
778 | |
779 \begin{alltt} | |
780 void | |
781 refdup(z_t r, z_t a) | |
782 \{ | |
783 \textcolor{c}{/* \textrm{Almost fully optimised, but perfectly po… | |
784 r->sign = a->sign; | |
785 r->used = a->used; | |
786 r->alloced = a->alloced; | |
787 r->chars = a->chars; | |
788 \} | |
789 \end{alltt} | |
790 | |
791 | |
792 \subsection{Variadic initialisation} | |
793 \label{sec:Variadic initialisation} | |
794 | |
795 Most bignum libraries have variadic functions | |
796 for initialisation and uninitialisation. This | |
797 is not available in libzahl, because it is | |
798 not useful enough and has performance overhead. | |
799 And what's next, support {\tt va\_list}, | |
800 variadic addition, variadic multiplication, | |
801 power towers, set manipulation? Anyone can | |
802 implement variadic wrapper for {\tt zinit} and | |
803 {\tt zfree} if they really need it. But if | |
804 you want to avoid the overhead, you can use | |
805 something like this: | |
806 | |
807 \begin{alltt} | |
808 /* \textrm{Call like this:} MANY(zinit, (a), (b), (c)) */ | |
809 #define MANY(f, ...) (_MANY1(f, __VA_ARGS__,,,,,,,,,)) | |
810 | |
811 #define _MANY1(f, a, ...) (void)f a, _MANY2(f, __VA_ARGS__) | |
812 #define _MANY2(f, a, ...) (void)f a, _MANY3(f, __VA_ARGS__) | |
813 #define _MANY3(f, a, ...) (void)f a, _MANY4(f, __VA_ARGS__) | |
814 #define _MANY4(f, a, ...) (void)f a, _MANY5(f, __VA_ARGS__) | |
815 #define _MANY5(f, a, ...) (void)f a, _MANY6(f, __VA_ARGS__) | |
816 #define _MANY6(f, a, ...) (void)f a, _MANY7(f, __VA_ARGS__) | |
817 #define _MANY7(f, a, ...) (void)f a, _MANY8(f, __VA_ARGS__) | |
818 #define _MANY8(f, a, ...) (void)f a, _MANY9(f, __VA_ARGS__) | |
819 #define _MANY9(f, a, ...) (void)f a | |
820 \end{alltt} |