Is it possible to find the eigenfunction of the following differential operator? Consider a differential operator of the form
$\sum_{j=-\infty}^\infty(z_{j+1}\partial_{z_j}+z_j\partial_{z_{j+1}}-g z_j-g\partial_{z_j})$
Here g is a constant.
The operator consisting up the first two terms $\sum_j(z_{j+1}\partial_{z_j}+z_j\partial_{z_{j+1}})$ has an eigenfunction of the form $\sum_m z_m e^{ikm}$ wth eigenvalue $2\cos k$ and the operator of the form $g(z_j+\partial_{z_j})$ has an eigenfunction of the form $e^{-z_j^2/2+Ez_j}$ with Eigenvalue $Eg$.
Is there a way to find an eigenfunction of the entire thing? Is there any ansatz that can help me proceed?
 A: I appreciate the answer above because it finds eigenfunctions that are valid in some sense even for an infinite number of variables, however the problem of finding a classification of the eigenfunctions is far from solved.
A characterization can be achieved in many ways, but in the following I will solve the original operator eigenfunction problem as posed with boundary condition for an arbitrary number of variables, first for $g=0$ and then for $g\neq 0$ by employing a transformation that reduces it to the $g=0$ problem.
The operator above for $N$ variables is denoted by $O_N(g)$ henceforth and the problem we are trying to solve is
$$O_{N}(0)f(z_1,..,z_N)=\lambda f(z_1,..., z_N)$$
Note that this is a first order differential equation, that can be solved using the method of characteristics. Trajectories are defined as usual by the characteristic equations
$$\frac{dz_{n}}{dt}=z_{n-1}+z_{n+1}~~~,~~~ z_{N+k}=z_k~,~ k\in \mathbb{Z} $$
$$\frac{df}{dt}=\lambda f$$
Because of the linearity of the equations, $z\propto e^{rt}$ are solutions, and due to the PBC's the candidate $r$'s are eigenvalues of a $N\times N$ circulant matrix. The eigenvalues and eigenvectors of any circulant matrix can be found easily by employing a discrete Fourier transform and it is immediately seen that
$$G_q(t)=\sum_{\ell=0}^{N-1}z_{r+1}(t)e^{\frac{2\pi i qr}{N}}\propto e^{E_q t}~~,~~E_q=2\cos\frac{2\pi q }{N}~~,~~ q=\{0,1,...N-1\}$$
From these expressions it is very easy to construct the complete set of integrals of motion $\{G_0^{E_q}G_q^{-E_0}\}$  and along with the knowledge of the particular solution $f_p=(z_1+...+z_N)^{\lambda/2}$ the most general eigenfunction of $O_N(0)$ is given by
$$f_{\lambda}(z_1,...,z_N)=(G_0)^{\lambda/2}F(G_0^{E_1}G_1^{-E_0},..., G_0^{E_{N-1}}G_{N-1}^{-E_0})=\left(\sum {z_i}\right)^{\lambda/2}F\left(\left(\sum {z_i}\right)^{E_1}\left(\sum {z_i}\omega^i\right)^{-E_0},..., \left(\sum {z_i}\right)^{E_{N-1}}\left(\sum {z_i}\omega^{(N-1)i}\right)^{-E_0}\right)$$
Now to solve the eigenvalue problem for $O_N(g)$, make the change of variables $w_j=z_j-\frac{g}{2}$ and set $f(z_1,...,z_N)=e^{g G_0/2}h(z_1,..., z_N)$, which leads to the equivalent eigenproblem
$$O_N(0)h(w_1,..., w_N)=\left(\lambda+\frac{Ng^2}{2}\right) h(w_1,..., w_N)$$
with the solution
$$f_{\lambda}(z_1,...,z_N)=e^{\frac{g}{2}\sum z_j}\left(\sum {z_i}-\frac{Ng}{2}\right)^{\lambda/2+Ng^2/4}F\left(\left(\sum {z_i}-\frac{Ng}{2}\right)^{E_1}\left(\sum {z_i}\omega^i\right)^{-E_0},..., \left(\sum {z_i}-\frac{Ng}{2}\right)^{E_{N-1}}\left(\sum {z_i}\omega^{(N-1)i}\right)^{-E_0}\right)$$
The reason that the solution to this problem, even though complete, doesn't generalize very well to the infinite variables case is twofold: firstly, the monomials of $z$ get augmented by spurious divergent constants (which wouldn't deter a physicist, but it should deter a mathematician) and secondly the spectrum $E_q$ associated with the characteristic equations becomes dense in the infinite variables limit, which indicates that for a solution to be possible in that regime, the eigenfunction has to be promoted to a functional in a certain way. I will relegate this analysis to a later time for now.
A: Start with replacing the infinite sum $\sum_{j=-\infty}^{\infty}$ by periodic boundaries $\sum_{j=0}^{N-1}$ where $j$ is replaced with $j \mod N$. Next make the ansatz that the eigenfunctions have the form $U(\sum z_i)$. Lets call $Z=\sum z_i$ and note that $\partial_i (\sum z_i) = U'(Z)$ independently of $i$. Thus we need to solve:
$$\sum_{j=0}^{N-1}(z_{j+1}\partial_{z_j}+z_j\partial_{z_{j+1}}-g z_j-g\partial_{z_j})U(Z) = \lambda U(Z)$$
Replacing all $\partial_i U $ with $U'$:
$$\sum_{j=0}^{N-1}(z_{j+1}U'(Z)+z_jU'(Z)-g z_jU(Z)-g U'(Z)) = \lambda U(Z)$$
The summations over the $z_i$ produces $Z$, and the sum over $g$ produces $Ng$, leaving us with
$$(2Z-Ng)U'  = (\lambda + g Z) U$$
hence
$$ (\log U)' = {\lambda + g Z \over (2Z-Ng)} $$
$$\log U ={1 \over 4} (g^2 N + 2 \lambda) \log(2 Z - g N) + (g Z)/2 + constant$$
$$U=(2 Z - g N)^{(g^2 N + 2 \lambda) \over 4} e^{g Z \over 2}$$
Note that this is only a formal solution, as it diverges in $z$ and the solution family is probably non unique. In addition, only solutions with ${(g^2 N + 2 \lambda) \over 4}$ being integer are analytic.
