After this question popped up again, it seemed to me to scream out for a use of the Fourier Transform (FT). I have decided to post this as an answer, since this approach is transparent and provides a systematic method to obtain all (distributional) solutions.
A simple formal manipulation (to be justified below) shows that the FT $g$ (with independent variable $z$ since, for reasons which will become clear below, we are regarding it as an entire function defined on the complex plane) of a solution satisfies the equation $$(e^{-iy}+1+iy)g(y)=0.$$$$(e^{-iz}+1+iz)g(z)=0.$$
In classical terms, this only provides the zero solution but in the sense of distributions it has many non-trivial ones, in fact, suitable combinations of $\delta$-functions with singularities at the zeros of the function in brackets. This The resulting solution for a simple pole, i.e. the inverse FT of the corresponding delta function is easily translates intoseen to be exactly one of the solutions given in the previous answers. This can be used to give a precise precise description of all possible solutions of the original equation.
[HistoricalIn order to motivate this approach, we begin with a historical aside on the FT for distributions: The problem of extending the FT to distributions, which wasthe latter setting was motivated by applications in physics, was and was studied in detail by many authors at the middle of the last century. Laurent Schwartz famously considered the case of tempered distributions in his seminal monograph. But this but similar constructions can be done forcarried out in many other situations. The basic starting point for such an extension is an initial setting where it isthe FT establishes an isomorphism between two specific l.c. spaces of test functions (that is, functions with good smoothness and/or growth properties). One then uses transposestransposition (in the sense of duality theory for l.c.´s) to translate this into an isomorphism between corresponding dual spaces. The latter consist of distributions(generalised) distributions. Schwartz used the (symmetric) case of the smooth functions of rapid decrease butto extend the FT to tempered distributions but there are many non symmetric cases which are of interest, as is the case here.]
We require the FT’s for arbitrary distributionsconcept and simple properties of the FT applied to ARBITRARY distributions on the line and so used the fact that the FTit is an isomorphism between the smooth functions on the line with compact support and a suitable space of entire functions of exponential growth. (This is the Paley-Wiener-Schwartz theorem, with details details of which can be found in StrichartzStrichartz´ monograph "A Guide to Distribution Theory and Fourier Transforms", where it is Th. 7.2.1 on p. 113). Dualising, we get a version of the FT which works for ALL distributions--it takes its values in a space of analytic functionals (a superspace of the dual space of the Fréchet space of entire functions). This can be described explicitly, but for our purposes it suffices to know that it contains the delta functions on the complex plane, together with suitable combinations. With this apparatus, the above formal computations can be tidied up to give a rigorous treatment which identifies all possible solutions of the above equation.