In General: Suppose we have $f(\alpha) = \alpha,$ while $f'(\alpha) = s$ with $0 < s < 1,$ also $f'(x) > 0$ and $f(x) \neq x$ for $x \neq \alpha,$ $C < x < D,$ for constants $C < \alpha < D.$ The first thing we do is move the fixpoint to the origin, with
$$ g(x) = f(x + \alpha) - \alpha. $$
Then $g'(0) = s.$ As a result, pages 122-123 in KCG or pages 46-47 in Alexander or pages 77-79 in Milnor, we can solve Schroeder's equation
$$ \psi(g(x)) = s \psi(x) $$
with $$ \psi'(0) = 1 $$
by
$$ \psi(x) = \lim_{n \rightarrow \infty} \frac{g^n(x)}{s^n}, $$
where $g^1(x) = g(x)$ and then
$$ g^n(x) = g \left( g^{n-1}(x) \right) $$
From above,
$$ s \psi(x) = \psi(g(x)) $$ so
$$ \psi^{-1} (s \psi(x)) = g(x ) = f(x+\alpha) - \alpha. $$
Substitute $x \mapsto x - \alpha$ to get
$$ \psi^{-1} (s \psi(x - \alpha)) = f(x) - \alpha. $$
At this point, we can use $\sqrt s$ to define the half-iterate,
$$ h(x) = \alpha + \psi^{-1} \left( \sqrt s \; \psi(x - \alpha) \right) $$
We will need the evident
$$ \psi \left( h(x) - \alpha \right) = \sqrt s \; \psi(x - \alpha) $$
We may now simply calculate
$$ h(h(x)) = \alpha + \psi^{-1} \left( \sqrt s \; \psi(h(x) - \alpha) \right) $$
$$ h(h(x)) = \alpha + \psi^{-1} \left( \sqrt s \; \psi( \alpha + \psi^{-1} \left( \sqrt s \; \psi(x - \alpha) \right) - \alpha) \right) $$
$$ h(h(x)) = \alpha + \psi^{-1} \left( \sqrt s \; \psi( \psi^{-1} \left( \sqrt s \; \psi(x - \alpha) \right) ) \right) $$
$$ h(h(x)) = \alpha + \psi^{-1} \left( \sqrt s \; \left( \sqrt s \; \psi(x - \alpha) \right) \right) $$
$$ h(h(x)) = \alpha + \psi^{-1} \left( s \; \psi(x - \alpha) \right) $$
$$ h(h(x)) = \alpha + f(x) - \alpha = f(x). $$
==========================================================================
This one has the advantage that the fixpoint at $\phi = \frac{1 + \sqrt 5}{2}$ has the function's derivative greater than one. This means that there is a real analytic half iterate for, at least, $x > 0.$ All you need to do is solve the Schroeder equation around $\phi.$ I will need to look up some things, but the references you want are Alexander and KCG.
Oh: it is not difficult to prove that there can be no answer for this that extends to the entire complex plane, or even extends to the entire real line. That is, however, not the end of the story.
Alright, pages 46 and 47 in Alexander, A History of Complex Dynamics. We need to use the inverse function to get derivative below $1,$ so we define
$$ h(x) = \sqrt {1+x}. $$ Or pages 122-123 in Kuczma, Choczewski, and Ger, Iterative Functional Equations.
I need some more time to make an actual computer program. I have done the more difficult case with derivative $1,$ including program, http://mathoverflow.net/questions/45608/formal-power-series-convergence/46765#46765 and later http://math.stackexchange.com/questions/911818/how-to-obtain-fx-if-it-is-known-that-ffx-x2x/912324#912324
Anyway, as I said, your problem is conceptually far simpler.
It is necessary to move the fixpoint to the origin, so
$$ g(x) = \sqrt {\phi^2 + x \;} - \phi. $$
We will then get the Schroeder function
$$ \psi(x) = \lim_{n \rightarrow \infty} (2 \phi)^n g^n(x), $$
where $g^1 (x) = g(x)$ and then $$g^n(x) = g (g^{n-1}(x))$$
Here is an example with $x=1,$ showing $\psi(1) \approx 0.88523$
1 0.2840790438404122 0.9192990968507166
2 0.08552494260243937 0.8956288264762806
3 0.02621627545939931 0.888431399703
4 0.008081094572713665 0.8862183392045823
5 0.002495271498931473 0.8855355193262247
6 0.0007708976551425994 0.8853246168742123
7 0.0002382029425400667 0.8852594540251384
8 7.360708310000241e-05 0.8852393185124893
9 2.274567970927954e-05 0.8852330963800049
10 7.028786312091029e-06 0.8852311736375063
11 2.172012962375902e-06 0.8852305794564203
12 6.711887781118975e-07 0.8852303957800051
13 2.074087255277135e-07 0.8852303388330724
14 6.40928197181978e-08 0.8852303215475368
15 1.98057703570953e-08 0.8852303146984961
16 6.120319584468348e-09 0.8852303085424645
17 1.891282686017348e-09 0.885230272694271
18 5.844384975972616e-10 0.8852302824650963
19 1.806013116834038e-10 0.8852297127013157
20 5.580891304646229e-11 0.8852303238549764
===========================================================
int main()
{
double phi = ( 1 + sqrt(5) ) / 2;
double phi2 = 1 + phi;
cout.precision(16);
cout << endl << phi << " " << phi2 << endl;
double x = 1.0;
double mult = 1;
for(int n = 1; n <= 20; ++n)
{
mult *= 2 * phi;
x = sqrt( phi2 + x) - phi;
cout << n << " " << x << " " << mult * x << endl;
}
return 0;
}
==========================================================
