# How to compute the values of this function ? ( Fabius function )

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/

-
Related, containing multiple answers: math.stackexchange.com/q/240687/49989 - in fact, my question is your answer and sort of vice-versa. – kram1032 Nov 25 '12 at 1:17

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:

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

-