How to compute the values of this function ? ( Fabius function ) It is said not to be analytic but $C^\infty$ everywhere. But I do not even know how to compute its values. Im confused. Here is the link : http://www.math.osu.edu/~edgar.2/selfdiff/


The Fabius function assumes rational values at dyadic rational arguments. We will give an explicit formula for those values. Because dyadic rationals form a dense subset of the reals, and the Fabius function is continuous, its value at any real point can be found as a limit of its values on a sequence of dyadic rationals converging to that point.

The following is based on the book Rvachev V. L., Rvachev V. A., "Non-classical methods of the approximation theory in boundary value problems", Naukova Dumka, Kiev (1979), pp. 116-126$^{[1]}$. The book is in Russian, and as far as I know, no English translation exists. But, fortunately, it mostly consists of math formulae.

Let's define $t_n$ by the recurrence $$t_0 = 1, \quad t_n = (-1)^n \, t_{\lfloor n/2\rfloor}\tag1$$ It is easy to see that $|t_n|=1$, and the signs follow the same pattern as the Thue–Morse sequence. The same sequence is given by a non-recurrent, but rather clumsy formula $$t_n = 1-\frac83 \, \sin^2\left(\frac\pi3 \cdot \sum_{m=1}^n \left[2 + (-1)^{\binom n m}\right]\right)\tag2$$

Let's define $c_n$ by the recurrence $$c_0 = 1, \quad c_n = \frac1{(4^n - 1)(2n+1)} \, \sum_{m\ge1} \binom{2n+1}{2m+1} \, c_{n-m},\tag3$$ so we have $$c_1 = \frac19, \quad c_2 = \frac{19}{675}, \quad c_3 = \frac{583}{59535}, \quad c_4 = \frac{132809}{32531625}, \quad \text{etc.}\tag4$$ Note that only finite number of terms in the sum are non-zero. I am not aware of any non-recurrent formula for this sequence, but it would be very nice to find one.

The values of the Fabius function at dyadic rationals are given by $$F\!\left(\tfrac{s}{2^n}\right) = \frac1{n! \, 2^{\binom{n+1}2}} \, \sum_{m\ge0}\binom n {2m} \, c_m \sum_{1 \le r \le s}(2r-1)^{n-2m} \, t_{s-r}.\tag5$$ Again, note that only finite number of terms in each sum are non-zero.

  • 1
    $\begingroup$ An efficient algorithm to compute values of the Fabius function (for all arguments $0\le x\le1$) implemented in the Wolfram Language can be found here. It gives exact rational values for dyadic rational arguments and arbitrary precision values for other arguments. $\endgroup$ – Vladimir Reshetnikov Feb 14 '17 at 18:18
  • $\begingroup$ A paper explaining some background for this algorithm: J. Arias de Reyna, Arithmetic of the Fabius Function. $\endgroup$ – Vladimir Reshetnikov Jun 14 '18 at 5:20

For the Fabius function we have $$\tag{1}F(0)=0\quad F(1)=1$$ $$F'(x)=\begin{cases} 2F(2x) & 0\leq x\leq 1/2 \\ 2F(2(1-x)) & 1/2\leq x\leq1 \end{cases}\tag{2} $$

Now for $h\in[0, 2^{-i}]$ consider the functions

$$P_i(h)=F(2^{-i}+h) - (-1)^i F(2^{-i}-h)$$

We have $$P_1'(h)=F'(1/2 + h) - F'(1/2 - h) \overset{(2)}{=} 0$$

hence $P_1(h)$ is a constant, which equals $1$ because of $(1)$.

For $i\geq 2$ using the first case of $(2)$ yields

$$P_i'(h)=2F(2^{1-i}+2h) - 2(-1)^{i-1} F(2^{1-i}-2h) = 2P_{i-1}(2h)$$

That is, $P_{i}(h)$ is just the integral of $2P_{i-1}(2h)$, i.e. a degree $i-1$ polynomial with coefficients that are simple to find, except for the additive constant from the integration. For the even $i$s, the constant can be found by solving $P_i(0) = 0$, and for the uneven $i$s, one needs to solve $P_i(0) = 2F(2^{-i})$.

This means that for $i < n$, $P_{i}(h)$ is easily found when knowing $F(2^{-i})$ for $i \in \{1, 3,...,N\}$

Knowing $P_i(h)$ allows us to solve $F(2^{1-i})$ solely in terms of $F(1)$ and $F(2^{-j})$ for $j\in\{1,3,...,i-2\}$. Doing so for the even $i\leq N$ turns out to make up a lower triangular linear system of equations with $U =2\lfloor N/2 \rfloor$ unknowns and rank $U-1$, because the first equation is $F(1/2) = F(1/2)$. However, the value of $F(1/2)$ is obviously $1/2$, by means of which forward substitution solves the rest.

Now that you know $P_i(h)$ for however many $i$s you'd like, you can approximate the value of $F(z)$ for any $z\in[0,1]$ by reflecting the $z$-value about the relevant symmetry-axes among $x=1/2$, $x=1/4$, $x=1/8$, etc. while solving the relevant $P_i(h)$-equation in terms of another Fabius function value which converge super fast towards 0, if you reflect your $z$ towards $x=0$.

On runtime the needed operations increases quadratically w.r.t. amount of reflections. Using the method for $z=2/3$ causes the need of a reflection about all of the axes because it's binary representation is $0.101010101$... That makes it a good choice to see the significant figures relative to the amount of reflections, which looks very much like a quadratic function as well:

enter image description here

In case anyone is interested, I created the method in Mathematica:

enter image description here

  • 4
    $\begingroup$ I re-posted the Mathematica code from the screenshot in plain text format here. $\endgroup$ – Vladimir Reshetnikov Feb 12 '17 at 1:31
  • $\begingroup$ I upvoted this answer. Hesitated to accept. $\endgroup$ – mick Nov 12 '18 at 10:05

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.