Let us slightly generalize your functional equation as $u(x,y)=v(x-y)$, $$u(x,y)=\sum_{j=1}^3f_j(x)g_j(y).$$ In your case, $f_3=1$. According to a theorem of Rubel and Gauchman, for a sufficiently smooth function $u$, the necessary and sufficient condition for such representation of sums of products of smooth functions is that the determinant $$\det\left(\frac{\partial^{i+j}u}{\partial x^i\partial y^j}\right)_{i,j=0,\ldots,3}\equiv0.$$ Using $u(x,y)=v(x-y)$ we obtain $u_x=v'$ and $u_y=-v'$, and in general, $$\frac{\partial^{i+j}}{\partial x^i\partial y^j}=(-1)^jv^{(i+j)}.$$ Substituting this to our determinant we obtain the Wrosnkian $$W(v,-v',v'',-v''')\equiv 0.$$ This means that $v,v',v'',v'''$ are linearly dependent, therefore $v$ satisfies a linear differential equation of at most $3$-d order with constant coefficients. So is a generalized exponential sum: $$v(x)=c_1e^{\lambda_1 x}+c_2e^{\lambda_2 x}+c_3e^{\lambda_3 x},$$ if the roots of the characteristic equation are distinct, or $$v(x)=c_1e^{\lambda_1 x}+(c_2+c_3)e^{\lambda_2 x},$$ if there are two roots, or $$v(x)=(c_1+c_2x+c_3x^2)e^{\lambda x},$$ if there is only one root. These are all possible forms of $v$. I leave it to you to determine exactly which of these functions satisfy your specific equation (with $f_3=1$). Solutions that you found correspond to case 3, with $\lambda=0$. Solution that I gave in the comment belongs to case 1 with $c_3=0$ or to case 2 with $c_3=0$.

Ref. MR1024482 Gauchman, Hillel; Rubel, Lee A. Sums of products of functions of x times functions of y. Linear Algebra Appl. 125 (1989), 19–63.

Let us slightly generalize your functional equation as $u(x,y)=v(x-y)$, $$u(x,y)=\sum_{j=1}^3f_j(x)g_j(y).$$ In your case, $f_3=1$. According to a theorem of Rubel and Gauchman, for a sufficiently smooth function $u$, the necessary and sufficient condition for such representation of sums of products of smooth functions is that the determinant $$\det\left(\frac{\partial^{i+j}u}{\partial x^i\partial y^j}\right)_{i,j=0,\ldots,3}\equiv0.$$ Using $u(x,y)=v(x-y)$ we obtain $u_x=v'$ and $u_y=-v'$, and in general, $$\frac{\partial^{i+j}}{\partial x^i\partial y^j}=(-1)^jv^{(i+j)}.$$ Substituting this to our determinant we obtain the Wrosnkian $$W(v,-v',v'',-v''')\equiv 0.$$ This means that $v,v',v'',v'''$ are linearly dependent, therefore $v$ satisfies a linear differential equation of at most $3$-d order with constant coefficients. So is a generalized exponential sum: $$v(x)=c_1e^{\lambda_1 x}+c_2e^{\lambda_2 x}+c_3e^{\lambda_3 x},$$ if the roots of the characteristic equation are distinct, or $$v(x)=c_1e^{\lambda_1 x}+(c_2+c_3)e^{\lambda_2 x},$$ if there are two roots, or $$v(x)=(c_1+c_2x+c_3x^2)e^{\lambda x},$$ if there is only one root. These are all possible forms of $v$. I leave it to you to determine exactly which of these functions satisfy your specific equation (with $f_3=1$) and which of them are concave or whatever other conditions you want to impose. Solutions that you found correspond to case 3, with $\lambda=0$. Solution that I gave in the comment belongs to case 1 with $c_3=0$ or to case 2 with $c_3=0$.

Ref. MR1024482 Gauchman, Hillel; Rubel, Lee A. Sums of products of functions of x times functions of y. Linear Algebra Appl. 125 (1989), 19–63.