path: root/doc/not-implemented.tex

\chapter{Not implemented}
\label{chap:Not implemented}

In this chapter we maintain a list of
features we have choosen not to implement,
but would fit into libzahl had we not have
our priorities straight. Functions listed
herein will only be implemented if there
is shown that it would be overwhelmingly
advantageous. For each feature, a sample
implementation or a mathematical expression
on which you can base your implemention.
The sample implementation create temporary
integere references, this is for simplicity
in the example, your should try to use
dedicated variables, in case of recursion,
a robust program should store temporary
variables on a stack, so they can be
clean up of something happens.

\vspace{1cm}
\minitoc


\newpage
\section{Extended greatest common divisor}
\label{sec:Extended greatest common divisor}

\begin{alltt}
void
extgcd(z_t bézout_coeff_1, z_t bézout_coeff_2, z_t gcd
       z_t quotient_1, z_t quotient_2, z_t a, z_t b)
\{
#define old_r gcd
#define old_s bézout_coeff_1
#define old_t bézout_coeff_2
#define s quotient_2
#define t quotient_1
    z_t r, q, qs, qt;
    int odd = 0;
    zinit(r), zinit(q), zinit(qs), zinit(qt);
    zset(r, b), zset(old_r, a);
    zseti(s, 0), zseti(old_s, 1);
    zseti(t, 1), zseti(old_t, 0);
    while (!zzero(r)) \{
        odd ^= 1;
        zdivmod(q, old_r, old_r, r), zswap(old_r, r);
        zmul(qs, q, s), zsub(old_s, old_s, qs);
        zmul(qt, q, t), zsub(old_t, old_t, qt);
        zswap(old_s, s), zswap(old_t, t);
    \}
    odd ? abs(s, s) : abs(t, t);
    zfree(r), zfree(q), zfree(qs), zfree(qt);
\}
\end{alltt}


\newpage
\section{Least common multiple}
\label{sec:Least common multiple}

\( \displaystyle{
    \mbox{lcm}(a, b) = \frac{\lvert a \cdot b \rvert}{\mbox{gcd}(a, b)}
}\)


\newpage
\section{Modular multiplicative inverse}
\label{sec:Modular multiplicative inverse}

\begin{alltt}
int
modinv(z_t inv, z_t a, z_t m)
\{
    z_t x, _1, _2, _3, gcd, mabs, apos;
    int invertible, aneg = zsignum(a) < 0;
    zinit(x), zinit(_1), zinit(_2), zinit(_3), zinit(gcd);
    *mabs = *m;
    zabs(mabs, mabs);
    if (aneg) \{
        zinit(apos);
        zset(apos, a);
        if (zcmpmag(apos, mabs))
            zmod(apos, apos, mabs);
        zadd(apos, mabs, apos);
    \}
    extgcd(inv, _1, _2, _3, gcd, apos, mabs);
    if ((invertible = !zcmpi(gcd, 1))) \{
        if (zsignum(inv) < 0)
            (zsignum(m) < 0 ? zsub : zadd)(x, x, m);
        zswap(x, inv);
    \}
    if (aneg)
        zfree(apos);
    zfree(x), zfree(_1), zfree(_2), zfree(_3), zfree(gcd);
    return invertible;
\}
\end{alltt}


\newpage
\section{Random prime number generation}
\label{sec:Random prime number generation}

TODO


\newpage
\section{Symbols}
\label{sec:Symbols}

\subsection{Legendre symbol}
\label{sec:Legendre symbol}

TODO


\subsection{Jacobi symbol}
\label{sec:Jacobi symbol}

TODO


\subsection{Kronecker symbol}
\label{sec:Kronecker symbol}

TODO


\subsection{Power residue symbol}
\label{sec:Power residue symbol}

TODO


\subsection{Pochhammer \emph{k}-symbol}
\label{sec:Pochhammer k-symbol}

\( \displaystyle{
    (x)_{n,k} = \prod_{i = 1}^n (x + (i - 1)k)
}\)


\newpage
\section{Logarithm}
\label{sec:Logarithm}

TODO


\newpage
\section{Roots}
\label{sec:Roots}

TODO


\newpage
\section{Modular roots}
\label{sec:Modular roots}

TODO % Square: Cipolla's algorithm, Pocklington's algorithm, Tonelli–Shanks algorithm


\newpage
\section{Combinatorial}
\label{sec:Combinatorial}

\subsection{Factorial}
\label{sec:Factorial}

\( \displaystyle{
    n! = \left \lbrace \begin{array}{ll}
        \displaystyle{\prod_{i = 0}^n i} & \textrm{if}~ n \ge 0 \\
        \textrm{undefined} & \textrm{otherwise}
    \end{array} \right .
}\)
\vspace{1em}

This can be implemented much more efficently
than using the naïve method, and is a very
important function for many combinatorial
applications, therefore it may be implemented
in the future if the demand is high enough.

An efficient, yet not optimal, implementation
of factorials that about halves the number of
required multiplications compared to the naïve
method can be derived from the observation

\vspace{1em}
\( \displaystyle{
    n! = n!! ~ \lfloor n / 2 \rfloor! ~ 2^{\lfloor n / 2 \rfloor}
}\), $n$ odd.
\vspace{1em}

\noindent
The resulting algorithm can be expressed

\begin{alltt}
   void
   fact(z_t r, uint64_t n)
   \{
       z_t p, f, two;
       uint64_t *ns, s = 1, i = 1;
       zinit(p), zinit(f), zinit(two);
       zseti(r, 1), zseti(p, 1), zseti(f, n), zseti(two, 2);
       ns = alloca(zbits(f) * sizeof(*ns));
       while (n > 1) \{
           if (n & 1) \{
               ns[i++] = n;
               s += n >>= 1;
           \} else \{
               zmul(r, r, (zsetu(f, n), f));
               n -= 1;
           \}
       \}
       for (zseti(f, 1); i-- > 0; zmul(r, r, p);)
           for (n = ns[i]; zcmpu(f, n); zadd(f, f, two))
               zmul(p, p, f);
       zlsh(r, r, s);
       zfree(two), zfree(f), zfree(p);
   \}
\end{alltt}


\subsection{Subfactorial}
\label{sec:Subfactorial}

\( \displaystyle{
    !n = \left \lbrace \begin{array}{ll}
      n(!(n - 1)) + (-1)^n & \textrm{if}~ n > 0 \\
      1 & \textrm{if}~ n = 0 \\
      \textrm{undefined} & \textrm{otherwise}
    \end{array} \right . =
    n! \sum_{i = 0}^n \frac{(-1)^i}{i!}
}\)


\subsection{Alternating factorial}
\label{sec:Alternating factorial}

\( \displaystyle{
    \mbox{af}(n) = \sum_{i = 1}^n {(-1)^{n - i} i!}
}\)


\subsection{Multifactorial}
\label{sec:Multifactorial}

\( \displaystyle{
    n!^{(k)} = \left \lbrace \begin{array}{ll}
      1 & \textrm{if}~ n = 0 \\
      n & \textrm{if}~ 0 < n \le k \\
      n((n - k)!^{(k)}) & \textrm{if}~ n > k \\
      \textrm{undefined} & \textrm{otherwise}
    \end{array} \right .
}\)


\subsection{Quadruple factorial}
\label{sec:Quadruple factorial}

\( \displaystyle{
    (4n - 2)!^{(4)}
}\)


\subsection{Superfactorial}
\label{sec:Superfactorial}

\( \displaystyle{
    \mbox{sf}(n) = \prod_{k = 1}^n k^{1 + n - k}
}\), undefined for $n < 0$.


\subsection{Hyperfactorial}
\label{sec:Hyperfactorial}

\( \displaystyle{
    H(n) = \prod_{k = 1}^n k^k
}\), undefined for $n < 0$.


\subsection{Raising factorial}
\label{sec:Raising factorial}

\( \displaystyle{
    x^{(n)} = \frac{(x + n - 1)!}{(x - 1)!}
}\), undefined for $n < 0$.


\subsection{Falling factorial}
\label{sec:Falling factorial}

\( \displaystyle{
    (x)_n = \frac{x!}{(x - n)!}
}\), undefined for $n < 0$.


\subsection{Primorial}
\label{sec:Primorial}

\( \displaystyle{
    n\# = \prod_{\lbrace i \in \textbf{P} ~:~ i \le n \rbrace} i
}\)
\vspace{1em}

\noindent
\( \displaystyle{
    p_n\# = \prod_{i \in \textbf{P}_{\pi(n)}} i
}\)


\subsection{Gamma function}
\label{sec:Gamma function}

$\Gamma(n) = (n - 1)!$, undefined for $n \le 0$.


\subsection{K-function}
\label{sec:K-function}

\( \displaystyle{
    K(n) = \left \lbrace \begin{array}{ll}
      \displaystyle{\prod_{i = 1}^{n - 1} i^i}  & \textrm{if}~ n \ge 0 \\
      1 & \textrm{if}~ n = -1 \\
      0 & \textrm{otherwise (result is truncated)}
    \end{array} \right .
}\)


\subsection{Binomial coefficient}
\label{sec:Binomial coefficient}

\( \displaystyle{
    \binom{n}{k} = \frac{n!}{k!(n - k)!}
    = \frac{1}{(n - k)!} \prod_{i = k + 1}^n i
    = \frac{1}{k!} \prod_{i = n - k + 1}^n i
}\)


\subsection{Catalan number}
\label{sec:Catalan number}

\( \displaystyle{
    C_n = \left . \binom{2n}{n} \middle / (n + 1) \right .
}\)


\subsection{Fuss–Catalan number}
\label{sec:Fuss-Catalan number} % not en dash

\( \displaystyle{
    A_m(p, r) = \frac{r}{mp + r} \binom{mp + r}{m}
}\)


\newpage
\section{Fibonacci numbers}
\label{sec:Fibonacci numbers}

Fibonacci numbers can be computed efficiently
using the following algorithm:

\begin{alltt}
   static void
   fib_ll(z_t f, z_t g, z_t n)
   \{
       z_t a, k;
       int odd;
       if (zcmpi(n, 1) <= 0) \{
           zseti(f, !zzero(n));
           zseti(f, zzero(n));
           return;
       \}
       zinit(a), zinit(k);
       zrsh(k, n, 1);
       if (zodd(n)) \{
           odd = zodd(k);
           fib_ll(a, g, k);
           zadd(f, a, a);
           zadd(k, f, g);
           zsub(f, f, g);
           zmul(f, f, k);
           zseti(k, odd ? -2 : +2);
           zadd(f, f, k);
           zadd(g, g, g);
           zadd(g, g, a);
           zmul(g, g, a);
       \} else \{
           fib_ll(g, a, k);
           zadd(f, a, a);
           zadd(f, f, g);
           zmul(f, f, g);
           zsqr(a, a);
           zsqr(g, g);
           zadd(g, a);
       \}
       zfree(k), zfree(a);
   \}

   void
   fib(z_t f, z_t n)
   \{
       z_t tmp, k;
       zinit(tmp), zinit(k);
       zset(k, n);
       fib_ll(f, tmp, k);
       zfree(k), zfree(tmp);
   \}
\end{alltt}

\noindent
This algorithm is based on the rules

\vspace{1em}
\( \displaystyle{
    F_{2k + 1} = 4F_k^2 - F_{k - 1}^2 + 2(-1)^k = (2F_k + F_{k-1})(2F_k - F_{k-1}) + 2(-1)^k
}\)
\vspace{1em}

\( \displaystyle{
    F_{2k} = F_k \cdot (F_k + 2F_{k - 1})
}\)
\vspace{1em}

\( \displaystyle{
    F_{2k - 1} = F_k^2 + F_{k - 1}^2
}\)
\vspace{1em}

\noindent
Each call to {\tt fib\_ll} returns $F_n$ and $F_{n - 1}$
for any input $n$. $F_{k}$ is only correctly returned
for $k \ge 0$. $F_n$ and $F_{n - 1}$ is used for
calculating $F_{2n}$ or $F_{2n + 1}$. The algorithm
can be speed up with a larger lookup table than one
covering just the base cases. Alternatively, a naïve
calculation could be used for sufficiently small input.


\newpage
\section{Lucas numbers}
\label{sec:Lucas numbers}

Lucas numbers can be calculated by utilising
{\tt fib\_ll} from \secref{sec:Fibonacci numbers}:

\begin{alltt}
   void
   lucas(z_t l, z_t n)
   \{
       z_t k;
       int odd;
       if (zcmp(n, 1) <= 0) \{
           zset(l, 1 + zzero(n));
           return;
       \}
       zinit(k);
       zrsh(k, n, 1);
       if (zeven(n)) \{
           lucas(l, k);
           zsqr(l, l);
           zseti(k, zodd(k) ? +2 : -2);
           zadd(l, k);
       \} else \{
           odd = zodd(k);
           fib_ll(l, k, k);
           zadd(l, l, l);
           zadd(l, l, k);
           zmul(l, l, k);
           zseti(k, 5);
           zmul(l, l, k);
           zseti(k, odd ? +4 : -4);
           zadd(l, l, k);
       \}
       zfree(k);
   \}
\end{alltt}

\noindent
This algorithm is based on the rules

\vspace{1em}
\( \displaystyle{
    L_{2k} = L_k^2 - 2(-1)^k
}\)
\vspace{1ex}

\( \displaystyle{
    L_{2k + 1} = 5F_{k - 1} \cdot (2F_k + F_{k - 1}) - 4(-1)^k
}\)
\vspace{1em}

\noindent
Alternatively, the function can be implemented
trivially using the rule

\vspace{1em}
\( \displaystyle{
    L_k = F_k + 2F_{k - 1}
}\)


\newpage
\section{Bit operation}
\label{sec:Bit operation unimplemented}

\subsection{Bit scanning}
\label{sec:Bit scanning}

Scanning for the next set or unset bit can be
trivially implemented using {\tt zbtest}. A
more efficient, although not optimally efficient,
implementation would be

\begin{alltt}
   size_t
   bscan(z_t a, size_t whence, int direction, int value)
   \{
       size_t ret;
       z_t t;
       zinit(t);
       value ? zset(t, a) : znot(t, a);
       ret = direction < 0
           ? (ztrunc(t, t, whence + 1), zbits(t) - 1)
           : (zrsh(t, t, whence), zlsb(t) + whence);
       zfree(t);
       return ret;
   \}
\end{alltt}


\subsection{Population count}
\label{sec:Population count}

The following function can be used to compute
the population count, the number of set bits,
in an integer, counting the sign bit:

\begin{alltt}
   size_t
   popcount(z_t a)
   \{
       size_t i, ret = zsignum(a) < 0;
       for (i = 0; i < a->used; i++) \{
           ret += __builtin_popcountll(a->chars[i]);
       \}
       return ret;
   \}
\end{alltt}

\noindent
It requires a compiler extension, if missing,
there are other ways to computer the population
count for a word: manually bit-by-bit, or with
a fully unrolled

\begin{alltt}
   int s;
   for (s = 1; s < 64; s <<= 1)
       w = (w >> s) + w;
\end{alltt}


\subsection{Hamming distance}
\label{sec:Hamming distance}

A simple way to compute the Hamming distance,
the number of differing bits, between two number
is with the function

\begin{alltt}
   size_t
   hammdist(z_t a, z_t b)
   \{
       size_t ret;
       z_t t;
       zinit(t);
       zxor(t, a, b);
       ret = popcount(t);
       zfree(t);
       return ret;
   \}
\end{alltt}

\noindent
The performance of this function could
be improve by comparing character by
character manually with using {\tt zxor}.


\subsection{Character retrieval}
\label{sec:Character retrieval}

\begin{alltt}
uint64_t
getu(z_t a)
\{
    return zzero(a) ? 0 : a->chars[0];
\}
\end{alltt}