Fabius function

The Fabius function on the unit interval.

In mathematics, the Fabius function is a smooth cumulative distribution function $F\colon \mathbb {R} \to [0,1]$ that is strictly increasing on $[0,1]$ but nowhere real-analytic. It is the distribution function of the random series $$X=\sum _{n=1}^{\infty }2^{-n}U_{n},$$ where the $U_{n}$ are independent random variables uniformly distributed on $[0,1]$.

Equivalently, it is characterized by $$F(0)=0,\;F(1-x)=1-F(x),\;F'(x)=2F(2x)\,\left(0\leq x\leq {\tfrac {1}{2}}\right).$$ Although every positive-order derivative of $F$ vanishes at both endpoints, $F$ is nonconstant on every subinterval of $[0,1]$.

The same object occurs, after a change of variables, as Rvachev's up function $u$, an even $C^{\infty }$ bump function supported on $[-1,1]$. Its Fourier transform is an infinite product of sinc functions, or equivalently a product of powers of cosines. The function also has an unusual arithmetic side: every value $F(a/2^{n})$ at a dyadic rational is rational, and the numerators, denominators and 2-adic valuations of these values obey strong divisibility laws.

The underlying distribution already appeared in work of Børge Jessen and Aurel Wintner in 1935. Jaap Fabius independently introduced it probabilistically in 1966 as an example of a nowhere-analytic smooth function, and V. A. Rvachev independently introduced the associated compactly supported solution of a functional–differential equation. Near zero, $F$ is smaller than every power of its argument; its exact rate of decay is governed by the lower branch of the Lambert W function together with a small correction that repeats each time the argument is halved. That correction alters $F$ by only about two parts in a million, but it does not fade away as the argument tends to zero, so no asymptotic formula that omits it can be exact.

Notation

The literature uses several incompatible conventions. In this article, $F$ denotes the bounded cumulative distribution function, extended by constants outside $[0,1]$; $u$ denotes Rvachev's compactly supported up function; and ${\mathcal {F}}$ denotes the signed global extension satisfying ${\mathcal {F}}'(x)=2{\mathcal {F}}(2x)$ on the whole real line. Logarithms are natural unless a base is written explicitly, $s_{2}(n)$ is the sum of the binary digits of $n$, and $\nu _{2}(n)$ is the exponent of 2 in $n$.

Definition and elementary properties

There is a unique $C^{\infty }$ function $F\colon \mathbb {R} \to [0,1]$ such that $$F(x)=0\quad (x\leq 0),\qquad F(x)=1\quad (x\geq 1),$$ and, for $0\leq x\leq 1$, $$F(1-x)=1-F(x),$$ while $$F'(x)=2F(2x)\qquad \left(0\leq x\leq {\tfrac {1}{2}}\right).$$ The constant continuation outside the unit interval is smooth because all positive-order derivatives vanish at the endpoints.

Using the symmetry, the differential equation can be written on the whole unit interval as $$F'(x)={\begin{cases}2F(2x),&0\leq x\leq {\tfrac {1}{2}},\\[3pt]2F(2-2x),&{\tfrac {1}{2}}\leq x\leq 1.\end{cases}}$$ The derivative is positive on $(0,1)$, so $F$ is strictly increasing there. In particular, $$F\!\left({\tfrac {1}{2}}\right)={\tfrac {1}{2}},\qquad F'\!\left({\tfrac {1}{2}}\right)=2,$$ and $$F'(1-x)=F'(x),\qquad F'(x)+F'\!\left({\tfrac {1}{2}}-x\right)=2\quad \left(0\leq x\leq {\tfrac {1}{2}}\right).$$ The endpoint flatness is $$F^{(m)}(0)=F^{(m)}(1)=0\qquad (m\geq 1).$$

Integrating the differential equation gives a self-integral equation, which together with the symmetry and the normalization characterizes $F$: $$F(x)=\int _{0}^{2x}F(t)\,dt\qquad \left(0\leq x\leq {\tfrac {1}{2}}\right).$$

Probabilistic and convolution constructions

Let $U_{1},U_{2},\ldots$ be independent random variables, each uniformly distributed on $[0,1]$. The series $$X=\sum _{n=1}^{\infty }2^{-n}U_{n}$$ converges almost surely and takes values in $[0,1]$. The Fabius function is $$F(x)=\Pr(X\leq x).$$ The decomposition $$X\ {\stackrel {d}{=}}\ {\frac {U+X'}{2}},$$ where $U$ is uniform on $[0,1]$ and $X'$ is an independent copy of $X$, gives $$F(x)=\int _{0}^{1}F(2x-t)\,dt.$$ For $0\leq x\leq {\tfrac {1}{2}}$, the constant continuation $F(t)=0$ for $t\leq 0$ reduces this to the self-integral equation above.

The elementary moments of $X$ are $$\operatorname {E} X={\frac {1}{2}},\qquad \operatorname {Var} (X)={\frac {1}{36}}.$$ Set $$Y=2X-1=\sum _{n=1}^{\infty }2^{-n}V_{n},\qquad V_{n}=2U_{n}-1.$$ The variables $V_{n}$ are independent and uniform on $[-1,1]$. The random variable $Y$ has density $u$, and $$F'(x)=2u(2x-1)\qquad (0\leq x\leq 1).$$ Thus $u$ is the infinite convolution $$u=h_{1}*h_{2}*h_{3}*\cdots ,\qquad h_{n}(t)=2^{n-1}\mathbf {1} _{[-2^{-n},\,2^{-n}]}(t),$$ where $\mathbf {1} _{A}$ denotes the indicator of $A$. This construction accounts for both the compact support and the product formula for the Fourier transform.

Rvachev's up function

The up function associated with $F$ is $$u(t)={\begin{cases}F(t+1),&-1\leq t\leq 0,\\F(1-t),&0\leq t\leq 1,\\0,&|t|>1.\end{cases}}$$ Equivalently, $$u(t)=F(1-|t|)\qquad (|t|\leq 1).$$ It is even, nonnegative, positive on $(-1,1)$, satisfies $u(0)=1$, has integral one, and has topological support exactly $[-1,1]$. Up to a constant factor it is the unique smooth compactly supported solution of the delay differential equation $$u'(t)=2{\bigl (}u(2t+1)-u(2t-1){\bigr )},$$ and the normalization $\int _{\mathbb {R} }u(t)\,dt=1$ fixes that factor. Every derivative of $u$ vanishes at $t=\pm 1$.

Fourier transform

Use the Fourier-transform convention $${\widehat {u}}(z)=\int _{\mathbb {R} }u(t)e^{-2\pi izt}\,dt.$$ Then ${\widehat {u}}$ is an entire even function with ${\widehat {u}}(0)=1$, and has the equivalent product representations $${\begin{aligned}{\widehat {u}}(z)&=\prod _{h=0}^{\infty }{\frac {\sin(\pi z/2^{h})}{\pi z/2^{h}}}\\&=\prod _{m=1}^{\infty }\left(\cos {\frac {\pi z}{2^{m}}}\right)^{m}\\&=\prod _{m=1}^{\infty }\left(1-{\frac {z^{2}}{m^{2}}}\right)^{1+\nu _{2}(m)}.\end{aligned}}$$ Consequently, the zeros of ${\widehat {u}}$ are precisely the nonzero integers, and the zero at $m\in \mathbb {Z} \setminus \{0\}$ has multiplicity $1+\nu _{2}(|m|)$.

Fourier inversion gives $$u(x)=\int _{\mathbb {R} }{\widehat {u}}(t)e^{2\pi itx}\,dt.$$ Because $u$ is smooth and compactly supported, its Fourier transform decreases faster than every inverse power on the real axis.

Partitions of unity and Fourier series

The integer translates of $u$ form a partition of unity: $$\sum _{k\in \mathbb {Z} }u(x+k)=1.$$ More generally, for every positive integer $n$, $$\sum _{k\in \mathbb {Z} }u\!\left(x+{\frac {k}{n}}\right)=n.$$ Each sum is locally finite because $u$ is compactly supported. These identities also follow from the zeros of ${\widehat {u}}$ and the Poisson summation formula.

The 2-periodic extension of $u|_{[-1,1]}$ has the absolutely and uniformly convergent cosine series $$u(x)={\frac {1}{2}}+\sum _{j=0}^{\infty }{\widehat {u}}\!\left({\frac {2j+1}{2}}\right)\cos {\bigl (}(2j+1)\pi x{\bigr )},\qquad -1\leq x\leq 1.$$ Only odd harmonics occur because ${\widehat {u}}(k)=0$ for every nonzero integer $k$.

Fourier–Legendre expansion

Because $u$ is even, only even Legendre polynomials occur in its Fourier–Legendre series: $$u(x)=\sum _{n=0}^{\infty }a_{n}P_{2n}(x),\qquad -1\leq x\leq 1,$$ where $$a_{n}={\frac {4n+1}{2}}\int _{-1}^{1}u(x)P_{2n}(x)\,dx.$$

A formula for the coefficients in terms of reciprocal-power values of $F$ is $${\begin{aligned}a_{n}={}&4^{-n}(4n+1)\sum _{k=0}^{n}(-1)^{n+k}{\binom {2n}{n+k}}{\binom {2n+2k}{2n}}(2k)!\,2^{\binom {2k+1}{2}}F\!\left(2^{-2k-1}\right).\end{aligned}}$$ The series converges absolutely and uniformly on $[-1,1]$, including at the endpoints.

Signed global extension and non-analyticity

The signed global extension on the nonnegative real axis. The signs of successive bumps follow the Thue–Morse sequence.

There is a unique extension ${\mathcal {F}}$ of $F$ to the whole real line that satisfies the self-differential equation everywhere: $${\mathcal {F}}'(x)=2{\mathcal {F}}(2x)\qquad (x\in \mathbb {R} ).$$ It is given by ${\mathcal {F}}(x)=0$ for $x\leq 0$, by ${\mathcal {F}}(x+1)=1-{\mathcal {F}}(x)$ for $0\leq x\leq 1$, and by $${\mathcal {F}}(x+2^{r})=-{\mathcal {F}}(x)\qquad (0\leq x\leq 2^{r})$$ for every positive integer $r$. Equivalently, it is the locally finite sum $${\mathcal {F}}(x)=\sum _{n=0}^{\infty }(-1)^{s_{2}(n)}u(x-2n-1).$$ It agrees with the bounded Fabius function on $[0,1]$, and on the interval $[2n,2n+2]$ is a translated copy of $u$ with sign $(-1)^{s_{2}(n)}$. Its sequence of signs is therefore the Thue–Morse sequence.

Repeated differentiation gives $${\mathcal {F}}^{(k)}(x)=2^{\binom {k+1}{2}}{\mathcal {F}}(2^{k}x)\qquad (k\geq 0).$$ In particular, every positive-order derivative vanishes at every integer.

If $q=a/2^{n}>0$ is a positive dyadic number, with $a$ a positive integer, then $${\mathcal {F}}^{(k)}(q)=2^{\binom {k+1}{2}}{\mathcal {F}}(2^{k}q).$$ For $k>n$, the argument $2^{k}q$ is an even integer, so these derivatives vanish. Hence the Taylor series of ${\mathcal {F}}$ at a dyadic point is a polynomial of degree at most $n$; when $a$ is odd, the degree is exactly $n$. The function does not agree locally with this polynomial, and the dyadic rationals are dense, so $F$ is nowhere real-analytic on $[0,1]$.

Moments and generating functions

Let $$c_{n}=\int _{-1}^{1}t^{2n}u(t)\,dt=\operatorname {E} (Y^{2n})\qquad (n\geq 0).$$ All odd moments of $u$ vanish. The even moments are rational and satisfy $$c_{0}=1,\qquad (2n+1)(2^{2n}-1)c_{n}=\sum _{k=0}^{n-1}{\binom {2n+1}{2k}}c_{k}\quad (n\geq 1).$$ The moment expansion of the Fourier transform is $${\widehat {u}}(z)=\sum _{n=0}^{\infty }(-1)^{n}{\frac {c_{n}}{(2n)!}}(2\pi z)^{2n}.$$ The moments are also related to Fabius values by $$c_{n}=2(2n)!\,2^{\binom {2n+1}{2}}F\!\left(2^{-2n-1}\right).$$

The moment-generating function of $X$ is $${\begin{aligned}D(z)&=\operatorname {E} (e^{zX})=\sum _{n=0}^{\infty }d_{n}{\frac {z^{n}}{n!}}\\&=\prod _{j=1}^{\infty }{\frac {e^{z/2^{j}}-1}{z/2^{j}}}\\&=e^{z/2}{\widehat {u}}\!\left({\frac {iz}{4\pi }}\right)\\&=1+z\int _{0}^{1}e^{zt}u(t)\,dt.\end{aligned}}$$ It obeys the dilation equation $$D(2z)={\frac {e^{z}-1}{z}}D(z),$$ and consequently $$d_{0}=1,\qquad (n+1)(2^{n}-1)d_{n}=\sum _{k=0}^{n-1}{\binom {n+1}{k}}d_{k}\quad (n\geq 1).$$ The coefficients are the raw moments $d_{n}=\operatorname {E} (X^{n})$ and satisfy the identity $$d_{n}=n!\,2^{\binom {n}{2}}F(2^{-n}),$$ which yields the values $F(2^{-n})$ in exact rational arithmetic.

Dyadic values and exact computation

Every value of $F$ at a dyadic rational is rational. A closed formula can be stated using the even moments $c_{k}$. For integers $n\geq 0$ and $0\leq a\leq 2^{n}$, $${\begin{aligned}F\!\left({\frac {a}{2^{n}}}\right)={}&{\frac {2^{-{\binom {n+1}{2}}}}{n!}}\sum _{h=0}^{a-1}(-1)^{s_{2}(h)}\sum _{k=0}^{\lfloor n/2\rfloor }{\binom {n}{2k}}(2a-2h-1)^{n-2k}c_{k}.\end{aligned}}$$ The appearance of the binary digit sum $s_{2}(h)$ reflects the Thue–Morse signs of the global extension.

Values at reciprocal powers of two

$n$	$F(2^{-n})$	$d_{n}=n!\,2^{\binom {n}{2}}F(2^{-n})$
0	$1$	$1$
1	${\tfrac {1}{2}}$	${\tfrac {1}{2}}$
2	${\tfrac {5}{72}}$	${\tfrac {5}{18}}$
3	${\tfrac {1}{288}}$	${\tfrac {1}{6}}$
4	${\tfrac {143}{2073600}}$	${\tfrac {143}{1350}}$
5	${\tfrac {19}{33177600}}$	${\tfrac {19}{270}}$
6	${\tfrac {1153}{561842749440}}$	${\tfrac {1153}{23814}}$
7	${\tfrac {583}{179789679820800}}$	${\tfrac {583}{17010}}$

There is also an efficient recursive evaluator. Arias de Reyna attributes the method to V. A. Rvachev and reports having reconstructed it from a Mathematica program posted on Mathematica Stack Exchange. Given $0<x\leq 1$, choose the unique integer $n\geq 0$ such that $$2^{-n}\leq x<2^{-n+1},$$ and put $y=x-2^{-n}$. Then $$F(x)=-F(y)+\sum _{k=0}^{n}{\frac {2^{\binom {k+1}{2}}F(2^{k-n})}{k!}}\,y^{k}.$$ For a dyadic input, replacing $x$ by $y$ removes its leading nonzero binary digit, so the iteration terminates after as many steps as there are 1s in the binary numerator, and every step is exact rational arithmetic.

Asymptotic behaviour near zero

The flatness at zero implies $$F(x)=o(x^{N})\qquad (x\to 0^{+})$$ for every fixed $N>0$, and more precisely $$\log F(x)\sim -{\frac {\log ^{2}x}{2\log 2}}.$$ A sharp description requires both the Lambert W function and a logarithmically periodic correction.

Sharp Lambert-W form

Put $$L=\log 2,\qquad C_{0}={\frac {6\gamma ^{2}+12\gamma _{1}-\pi ^{2}}{12}},$$ where $\gamma$ is the Euler–Mascheroni constant and $\gamma _{1}$ is the first Stieltjes constant, normalized by $\zeta (s)={\tfrac {1}{s-1}}+\gamma -\gamma _{1}(s-1)+\cdots$. For sufficiently small positive $x$, put $$W=W_{-1}(-Lx),\qquad \lambda (x)=-{\frac {W}{L}},$$ where $W_{-1}$ is the lower real branch of the Lambert $W$ function; equivalently, $\lambda (x)$ is the large positive solution of $\lambda \,2^{-\lambda }=x$. The sharp main term is $$M(x)=-{\frac {7L}{12}}-{\frac {1}{2}}\log(\pi x)+{\frac {C_{0}-W-{\tfrac {1}{2}}W^{2}}{L}}+\Psi \!{\bigl (}\lambda (x){\bigr )},$$ where $\Psi$ is the one-periodic function described below, and, as $x\to 0^{+}$, $$\log F(x)=M(x)+O\!\left({\frac {1}{-\log x}}\right).$$

The periodic correction

The function $\Psi$ is real-valued, smooth and one-periodic, and is normalized to have mean zero over a period. With $\chi _{k}=2\pi ik/L$, its absolutely convergent Fourier series is $$\Psi (t)=\sum _{k\in \mathbb {Z} \setminus \{0\}}\left(-{\frac {\Gamma (-\chi _{k})\zeta (1-\chi _{k})}{L}}\right)e^{2\pi ikt},$$ where $\Gamma$ is the gamma function and $\zeta$ is the Riemann zeta function. The coefficient of index $k$ is nonzero for every $k\neq 0$, because the gamma function has no zeros and the zeta function has no zeros on the line $\operatorname {Re} s=1$, on which $1-\chi _{k}$ lies. Hence $\Psi$ is nonconstant, and deleting $\Psi (\lambda (x))$ from $M(x)$ leaves a bounded oscillation in $\log F(x)$ that the $O(1/(-\log x))$ remainder cannot absorb. The oscillation is small: $\Psi$ has amplitude about $2.1\times 10^{-6}$.

Expanded elementary form

The Lambert phase can itself be expanded in powers of $-\log x$ and $\log(-\log x)$. Writing $$\ell =\log x,\qquad a=\log(-\ell ),\qquad b=\log L,$$ the sharp form becomes $${\begin{aligned}\log F(x)={}&-{\frac {\ell ^{2}}{2L}}+{\frac {\ell a}{L}}-\left({\frac {1}{2}}+{\frac {1+b}{L}}\right)\ell -{\frac {a^{2}}{2L}}+{\frac {ba}{L}}\\&+\left({\frac {6\gamma ^{2}+12\gamma _{1}-\pi ^{2}-6b^{2}}{12L}}-{\frac {7L}{12}}-{\frac {\log \pi }{2}}\right)\\&+{\frac {a^{2}}{2L\ell }}-{\frac {ba}{L\ell }}+\Psi \!{\bigl (}\lambda (x){\bigr )}+O\!\left({\frac {1}{-\ell }}\right).\end{aligned}}$$ At this precision the periodic term is evaluated at the exact phase $\lambda (x)$.

For reciprocal powers of two, let $\lambda _{n}=\lambda (2^{-n})$, so that $\lambda _{n}-\log _{2}\lambda _{n}=n$, and put $$C_{\mathrm {sharp} }={\frac {6\gamma ^{2}+12\gamma _{1}-\pi ^{2}}{12L}}-{\frac {7L}{12}}-{\frac {\log \pi }{2}}.$$ Then $${\begin{aligned}\log F(2^{-n})={}&-{\frac {L}{2}}n^{2}-n\log n+\left(1+{\frac {L}{2}}\right)n-{\frac {\log ^{2}n}{2L}}+C_{\mathrm {sharp} }\\&-{\frac {\log ^{2}n}{2L^{2}n}}+\Psi (\lambda _{n})+O\!\left({\frac {1}{n}}\right).\end{aligned}}$$ The oscillation is periodic in the exact Lambert phase, not in the integer $n$ itself.

Origin of the logarithmic oscillation

For $s>0$ the moment-generating function has the product $$D(-s)=\prod _{j=1}^{\infty }{\frac {1-e^{-s/2^{j}}}{s/2^{j}}},$$ whose logarithm has an exact dyadic-scale decomposition $$\log D(-s)=-{\frac {\log ^{2}s}{2L}}+{\frac {1}{2}}\log s+R(\log _{2}s)+E(s),$$ with $R$ one-periodic and $|E(s)|\leq e^{-s}/(1-e^{-s})^{2}$. Rescaling by a factor of 2 therefore contributes an additive period-one term in the logarithmic variable, and saddle-point inversion carries the centred form of that term to $\Psi (\lambda (x))$ in the asymptotic of $F(x)$.

Full saddle-point expansion

There are continuous one-periodic functions $P_{1},P_{2},\ldots$ such that, for every integer $N\geq 1$, $$\log F(x)=M(x)+\sum _{j=1}^{N-1}{\frac {P_{j}(\lambda (x))}{\lambda (x)^{j}}}+O\!\left(\lambda (x)^{-N}\right),$$ where $M$ is the sharp main term. The first correction is $$P_{1}(t)={\frac {1}{24}}+{\frac {1}{2L}}-{\frac {\Psi ''(t)+(\Psi '(t))^{2}}{2L^{2}}}.$$

History

Arias de Reyna's historical survey records at least six independent appearances of the function or an equivalent one.

• In 1935, Jessen and Wintner considered the associated distribution and obtained its smoothness through the infinite Fourier product.
• In 1966, Fabius gave the probabilistic construction as an example of a smooth nowhere-analytic function.
• In 1971, V. A. Rvachev introduced the compactly supported solution now called the up function through its functional–differential equation; the work is recorded in a monograph by V. L. Rvachev and V. A. Rvachev, and a later survey by V. A. Rvachev developed a broader theory of compactly supported solutions of such equations.
• G. Kh. Kirov and G. A. Totkov gave another independent construction in 1981.
• Arias de Reyna developed the up function systematically in 1982, proving uniqueness, the probability interpretation, partition-of-unity identities, derivative formulas, rationality at dyadic points and exact computation formulas. An English translation appeared in 2017.
• R. Schnabl gave a further independent construction in 1985.

Arias de Reyna's 2018 paper developed the arithmetic of the dyadic values, including the exact 2-adic valuation, denominator bounds and a terminating exact evaluator. Haugland gave another explicit treatment of exact evaluation at arbitrary dyadic arguments.

See also

• Bump function
• Thue–Morse sequence
• Functional differential equation
• Minkowski's question-mark function
• Lambert W function

Further reading

• Rvachev, V. L.; Rvachev, V. A. (1979). Non-classical methods of approximation theory in boundary-value problems (in Russian). Kiev: Naukova Dumka.
• Dimitrov, Youri (2006). Polynomially-divided solutions of bipartite self-differential functional equations (PhD thesis). Ohio University.

External links

• Lean 4 formalization and exact evaluator
• FabiusF at the Wolfram Function Repository
• OEIS: A272755 (numerators of $F(2^{-n})$) and OEIS: A272757 (denominators)
• MathOverflow question on the reciprocal-power values