I have found a recursive solution to a problem which I am seeking a closed formula. The problem is equivalent to obtaining a function of two variables $P(a,b)$ with the following conditions:
$$0 \le a \le b \le 1 \\ P(a,a) = 0 \\ P(0,b) = b \\ P(a,b)=b-a+\frac{1}{a} \int_0^a P(a-x,b-x)dx$$
Is there a way to derive the close formula for $P(a,b)$?
Alternatively I would be satisfied with obtaining a closed form for $P(a,1)$
Let $b=c+a$ and $P(a,b)=P(a,a+c)=f(a)$. Then (in the last equality I use a change of variables), $$ f(a)=c+1/a\int_0^af(a-s)ds= c+1/a\int_0^af(s)ds. $$ Derivating, we get $$ f’(a)=\frac{1}{a}f(a)-\frac{1}{a^2}\int_0^af(s)ds. $$ Combining the two equations, we get $$ f’(a)= \frac{1}{a}f(a)-\frac{1}{a}f(a)+c/a=\frac{c}{a}. $$ Integrating, we get $$ f(a)=c\log(a)+g $$ for some constant $g$ which could depend on $c$. So, a general formula for $P(a,b)$ is $$ P(a,b)=(b-a)\log(a)+g(a-b), $$ where $g$ is some function (we said it would depend on $c$, that is, $b-a$). The problem now is, for any function $g$ you can never get the second boundary condition because the logarithm diverges there. So either you interpret your equation in a weaker sense, or the problem has no solution. The one I have found, though, is an explicit formula for the integral equation as it is (without boundary conditions). Note that the first condition $P(a,a)=0$ is satisfied by imposing $g(0)=0$.
Edit: I From your comments, if you consider $$ f(a)=c+a\int_0^af(a-s)ds= c+a\int_0^af(s)ds \qquad (1)$$ then yes, you get $$ f’(a)=\frac{(1+a^2)}{a}f(a)-c/a. $$ Formally, it is possible to solve this equation using the Duhamel formula. Calling $g(a)= \frac{(1+a^2)}{a}$, $h(a)=-c/a$: $$ f(a)=e^{\int_{a_0}^a g(r)dr}\left[ f(a_0)+\int_{a_0}^a e^{-\int_{a_0}^sg(r)dr}h(s)ds \right] = $$ $$ =e^{\int_{a_0}^a g(r)dr} f(a_0)+\int_{a_0}^a e^{\int_{s}^ag(r)dr}h(s)ds. $$ It seems there is still a problem if you choose $a_0=0$ since $g$ is not integrable near $0$. But it is likely that this time there are some cancellations when you try to compute the integrals, assuming $f(0)=c$. You could try to do so and check whether “some infinities cancel out” in those integrals.
Edit 2: I have found the following. Given $z$ a primitive of the function $l(a)=-\frac{e^{-a^2/2}}{a^2}$, the solution is given by $$ f(a)=cae^{a^2/2}z(a). $$ This is not very nice, since I think $z$ cannot be espressed in terms of ordinary functions. Such a function $z$ goes to infinity at $a=0$, but it gets balanced by the factor $a$ in the product. So, this seems to be as far as one can get.
Please note that there is also another problem: $z$ is defined only up to a constant, and to determine this constant you cannot rely on the values of $f$ at $a=0$: in fact, you can clearly see from the above formula that adding a constant to $z$ does not change the limit of the function as $a$ goes to zero. My guess at this point is that $f$ is not uniquely determined by equation $(1)$ and by the boundary conditions ($P(0,b)=b$ simply becomes $f(0)=c$, which is probably always satisfies by this family of solutions), and you still miss some hypotheses to determine $P$. I hope this is helpful for you.