Find the approximation of the inverse of binary entropy function

I need to find a good approximation for the binary entropy function:

$H_2(x) = -x \log_2(x) -(1-x) \log_2(1-x)$ where x $\in [0,1].$

I tried to find this using matematica and matlab but these tools gives to me an error. Since it seems that this function is not invertible.

As you said, I suppose that the function$$h(x) = -x \log_{2}(x) -(1-x) \log_{2}(1-x)\tag1$$ is not invertible.

However, if you are not too concerned by the values close to the boundaries, $h(x)$ can be approximated by a $[4,4]$ Padé approximant built at $x=\frac 12$. This would give $$h(x) \sim \frac{1- \left(\frac{124}{49}+\frac{2}{\log (2)}\right)\left(x-\frac{1}{2}\right)^2+ \left(\frac{152}{245}+\frac{548}{147 \log (2)}\right)\left(x-\frac{1}{2}\right)^4 } {1-\frac{124}{49} \left(x-\frac{1}{2}\right)^2+\frac{152}{245} \left(x-\frac{1}{2}\right)^4 }\tag2$$ which reduces to a quadratic equation in $\left(x-\frac{1}{2}\right)^2$.

For sure, this could be improved using the $[6,6]$ Padé approximant which would reduces to a cubic equation in $\left(x-\frac{1}{2}\right)^2$ which would be less pleasant.

In order to check the quality of the approximation, give $x$ a value in $(1)$ to get $h(x)$ and for this value, solve $(2)$ to get the solution. $$\left( \begin{array}{cc} x_{\text{given}} & x_{\text{recomputed}} \\ 0.05 & 0.048894 \\ 0.10 & 0.099763 \\ 0.15 & 0.149949 \\ 0.20 & 0.199990 \\ 0.25 & 0.249998 \\ 0.30 & 0.300000 \\ 0.35 & 0.350000 \\ 0.40 & 0.400000 \\ 0.45 & 0.450000 \\ 0.50 & 0.500000 \\ 0.55 & 0.550000 \\ 0.60 & 0.600000 \\ 0.65 & 0.650000 \\ 0.70 & 0.700000 \\ 0.75 & 0.750002 \\ 0.80 & 0.800010 \\ 0.85 & 0.850051 \\ 0.90 & 0.900237 \\ 0.95 & 0.951106 \end{array} \right)$$

Edit

The approximation can be improved using again

$$h(x) \sim f(x)=\frac{1+ a\left(x-\frac{1}{2}\right)^2+ b\left(x-\frac{1}{2}\right)^4 } {1+c \left(x-\frac{1}{2}\right)^2+d \left(x-\frac{1}{2}\right)^4 }$$ In order to respect the boundary conditions, we need to fix $\color{red}{b=-4 (a+4)}$. Now, considering the norm $$\Phi(a,c,d)=\int_0^1 \big( h(x)-f(x) \big) ^2\,dx $$ its numerical minimization with respect to $a$, $c$ and $d$ (this is equivalent to a nonlinear regression for an infinte number of data points) leads to $$\{a=-7.14483221448716,\,\, c=-4.28551864360839,\,\,d=2.63456254971723\}$$ At this point $\Phi(a,c,d)=3.44\times 10^{-7}$ while using the Padé approximant the norm would be $3.94\times 10^{-5}$ (that is to say $114$ times larger).

Update

If we take the problem for any base $a$,we have $$h_a(x) = -x \log_{a}(x) -(1-x) \log_{a}(1-x)$$ that is to say $$h_a(x)\log_e(a)=-x \log_{e}(x) -(1-x) \log_{e}(1-x)$$ and the $[4,4]$ Padé approximant of the rhs is

$$\frac{\log_e(2) + a\left(x-\frac{1}{2}\right)^2+ b\left(x-\frac{1}{2}\right)^4 }{{1-\frac{124}{49} \left(x-\frac{1}{2}\right)^2+\frac{152}{245} \left(x-\frac{1}{2}\right)^4 } } $$ where $$a=-\left(2+\frac{124 }{49}\log_e (2)\right)\qquad \text{and} \qquad b=\frac{548}{147}+\frac{152 }{245}\log_e (2)$$ and the norm is equal to $1.89\times 10^{-5}$.

Update

Answering this old question, I proposed two rather good approximations of the inverse $$x_0=\frac{1}{2} \left(1-\sqrt{1-h^{4/3}(x)}\right)$$ $$x_1=\frac{\log \left(1+\sqrt{1-h^{4/3}(x)}\right)+(y-1) \log (2)}{2 \tanh ^{-1}\left(\sqrt{1-h^{4/3}(x)}\right)}$$