Consider the second order linear elliptic differential equation

$$Lu=(\sum_{i=1}^d{\partial^2\over\partial\theta_i^2}+{\partial b\over\partial \theta_i}{\partial\over\partial\theta_i})u=exp(i\theta_1)$$

This differential equation is defined on the domain $[0,2\pi]^d$ with periodic boundary condition. The function $b$ are some polynomial in terms of sines and cosines in $\theta_k$, $k=1,..,d$. Is there some algorithm to compute the coefficients of the power series or better yet the Fourier series of the solution $u$?