Solve functional equation $ h(y)+h^{-1}(y)=2y+y^2 $

I was trying to solve a certain physics problem, and encountered the functional equation that contains a function $h$ and its inverse $h^{-1}$: \begin{equation} h(y)+h^{-1}(y)=2y+y^2.\tag{1} \end{equation}

Q: Does $(1)$ have a unique solution and is it possible to find it in closed form?

Equation $(1)$ looks quite simple and probably it is simple to analyze too, but I couldn't figured out how.

Background information: The physics problem I was trying to solve is finding the dependence of the current $J_0$ on the voltage $U_0$ in this infinite chain that contains nonlinear elements with quadratic volt-ampere characteristics $I(V)=\alpha V^2$ and ohmic resistors $R$:

According to dimensional analysis one can write $$ J_0(U_0)=\frac{1}{\alpha R^2}f(\alpha RU_0). $$ Solving simple system of equations I obtained the following functional equation for the unknown function $f$ $$ f(x)=(x-f(x))^2+f(x-f(x)).\tag{2} $$ Now introduce another function $h$ according to $$ x-f(x)=h(x). $$ Then equation $(2)$ becomes $$ x-h(x)=h^2(x)+h(x)-h(h(x)). $$ Let $h(x)=y$, then $x=h^{-1}(y)$ and we get $(1)$.

This is not a textbook problem and I don't even know whether it has a solution. I was studying it out of curiosity.

Solution 1:

Unfortunately, your proposed functional equation: $$ h(y)+h^{-1}(y)=2y+y^2 \tag{1} $$ has no differentiable solutions near $\,y=0.\,$ Suppose that $$ h(y) = 0 + a_1 y + O(y^2). \tag{2} $$ Then its inverse function is $$ h^{-1}(y) = 0 + \frac1{a_1} y + O(y^2). \tag{3} $$ Substituting equations $(2)$ and $(3)$ into $(1)$ gives $$ h(y)\!+\!h^{-1}(y) = \frac{a_1^2+1}{a_1}y \!+\! O(y^2) = 2y \!+\! O(y^2). \tag{4}$$ The only solutions for $\,a_1\,$ are $\,a_1=1\,$ and $\,a_1=-1.\,$ However, we are given that $$ 0<h(x)<x \tag{5}$$ and thus the only possible value is $\,a_1=1.\,$ Suppose we want more terms in the power series $$ h(y) = 0 + y + a_2y^2 + O(y^3). \tag{6}$$ The inverse function is now $$ h^{-1}(y) = 0 + y - a_2y^2 + O(y^3). \tag{7}$$ Adding these two equations together is the equation $$ h(y) + h^{-1}(y) = 0 + 2y + 0y^2 + O(y^3). \tag{8}$$ This implies that equation $(1)$ can not be satisfied.

There is a possibility that there exists an exponent $\,e\,$ not an integer such that $\, h(y) = 0 + y + c y^e + \cdots \,$ which perhaps should be studied according to comments.

In fact, define $$ g(x) : = \sqrt{h(x^2)}. \tag{9}$$ Then solving for a power series expansion gives $$ g(x) = x - \frac{x^2}{\sqrt{6}} + \frac{x^3}6 + O(x^4) \tag{10}$$ which implies that $$ h(x^2) = (x^2+x^4/2)+\frac1{\sqrt{6}}f(x) \tag{11} $$ where $$ f(x) \!=\! -2 x^3 \!-\!\frac{11x^5}{24} \!+\! \frac{117x^7}{1280} \!-\!\frac{5491x^9}{110592} \!+\! \frac{156538363x^{11}}{3715891200} \!+\! O(x)^{13}. \tag{12}$$

A Wolfram Language code to calculate $\,g(x)\,$ is

ClearAll[x, g, gx];
gx[3] = x - 1/Sqrt[6]*x^2 + O[x]^3;
Do[g = Normal[gx[n]] + O[x]^(2+n); gx[n+1] = Simplify[
    g + (x^2 + (Normal[g]/.x -> g)^2 - g^4 - 2*g^2) * 3 /
    ((4+n)*x^2*Sqrt[6])], 
  {n, 3, 6}] 

As some comments suggest, it seems that the power series for $\,f(x)\,$ has zero radius of convergence. This makes finding properties of it difficult. Perhaps we need to let $\,x\,$ approach $\,\infty.$ In this case we find that $\,h(x) \approx \sqrt{x}\,$ with infinite number of other terms. Using equation $(1)$ we can find the expansion $$ h(x) = x^{1/2} - 1 + \frac12 x^{-1/4} + \frac1x s(x) \tag{13} $$ where $$ s(x) := \sum_{k=2}^\infty 2^{-k}(x^{-\frac32 2^{-k}} - x^{-2^{-k}} ). \tag{14} $$ This gives moderately good approximations as $\,x\,$ gets large and down to $1.$

One method that leads to equations $(13)$ and $(14)$ is as follows. From equation $(1)$ we immediately get $$ h^{-1}(y) = y^2 + 2y - h(y) \tag{15} $$ and if $\,y=h(x),\,$ then $$ x = h(x)^2 + 2h(x) - h(h(x)). \tag{16} $$ We start with an approximation and try to find what additional term will satisfy equation $(16)$. So we guess $$ h(x) = x^{1/2} + cx^e + \cdots \tag{17} $$ where $\,\cdots\,$ denotes terms with smaller exponents. Substitute equation $(17)$ in equation $(16)$ to get $$ x \!=\! (x \!+\! 2cx^{1/2+e} \!+\! \cdots) \!+\! 2x^{1/2} \!+\! \cdots \!=\! x \!+\! 2x^{1/2}(cx^e \!+\! 1) + \cdots. \tag{18} $$ This implies $\,0 = cx^e + 1.\,$ Solving this we get our next guess as $$ \,h(x) = x^{1/2} - 1 + cx^e + \cdots. \tag{19} $$ Repeating this process leads to equations $(13)$ and $(14)$.

The series in $(14)$ appears to converge but I have no proof, only numerical evidence. I also have no proof that $(13)$ satisfies the functional equation $(1)$.

An answer to this question by 'Semiclasical' contains sequence recursions $$ u_n = u_{n-1} - j_{n-1} \quad \text{ and } \quad j_n = j_{n-1} - u_n^2. \tag{20} $$ with the property that $\, u_n = h(u_{n-1})\,$ and that $\,u_n \to 0^+.\,$ I found that for suitable starting values of $\,u_1\,$ and $\,j_1\,$ we have $$ u_n = 6n^{-2} - \frac{15}2n^{-4} + \frac{663}{40} n^{-6} - \frac{43647}{800}n^{-8} + \cdots. \tag{21} $$

An answer to this question by 'Cesareo' uses recursion to construct a sequence of functions $\,h_n(x)\,$ which seems to converge to a global solution. It would be nice to give a proof of this.

Solution 2:

From

$$ x = 2h(x)+h^2(x)-h(h(x))\Rightarrow h(x) = \pm\sqrt{h(h(x))+x+1}-1 $$

now using an iterative approximation over the positive leaf such as

$$ h_{k+1}(x) = \sqrt{h_k(h_k(x))+x+1}-1 $$

or according to the MATHEMATICA script which follows, (inspired after a fruitful discussion with Semiclassical)

Clear[h]
h[x_, 1] := x
h[x_, n_] := h[x, n] = Sqrt[h[h[x, n - 1], n - 1] + x + 1] - 1

we obtain beginning with $h_1(x) = x$ the successive approximations in red, and black for $h_5(x)$ showing a good convergence for $x > 0$.

Also the accomplishment for $h^{-1}(x)+h(x) = 2x+x^2$ can be analyzed in the following plot

err5 = y - h[2 y + y^2 - h[y, 5], 5];
Plot[Abs[err5], {y, 0, 4}, PlotStyle -> {Thick, Black}]