$f(x)$, $g(x)$ and $h(x)$ are known continuous functions on a bounded domain (specifically, probability density functions with known parameters defined on a finite interval). It is required to solve for the smooth positive monotonic function $p(x)$ such that $$p(2x)f(x)+p(x)g(x)=h(x)+g(x)$$
In principle this could be solved by a direct iterative numerical procedure -- any other ideas/suggestions?
If you know what $p(0)$ is, you might be able to express it as a series. Write the functional equation as $$p(x)=\frac{h(x/2)+g(x/2)}{f(x/2)}-\frac{g(x/2)}{f(x/2)}p(x/2)$$ and then plug this into itself over and over again, giving us $$p(x)=-\frac{g(x/2^n)}{f(x/2^n)}p(x/2^n)+\sum_{k=1}^{n}(-1)^{k+1}\frac{h(x/2^k)+g(x/2^k)}{f(x/2^k)}\prod_{i=1}^{k-1}\frac{g(x/2^i)}{f(x/2^i)}$$ where the empty product $\prod_{i=1}^0=1$ according to convention. Letting $n\to\infty$ and taking advantage of continuity, we have $$p(x)=-\frac{g(0)p(0)}{f(0)}+\sum_{k=1}^{\infty}(-1)^{k+1}\frac{h(x/2^k)+g(x/2^k)}{f(x/2^k)}\prod_{i=1}^{k-1}\frac{g(x/2^i)}{f(x/2^i)}$$ This is an infinite series representation for $p(x)$ (which may or may not converge depending on the functions $f,g,h$), and it looks pretty ugly but could simplify nicely depending on what functions you choose for $f,g,h$.
You could do a similar process (but backwards) which deals with $2^n$ rather than $2^{-n}$ and depends on $p(\infty)$.