  In one dimension (1D), a general decaying long-range interaction can be fit
to a sum of exponential interactions $e^{-\lambda r_{ij}}$ with varying
exponents $\lambda$, each of which can be represented by a simple matrix
product operator (MPO) with bond dimension $D=3$. Using this technique,
efficient and accurate simulations of 1D quantum systems with long-range
interactions can be performed using matrix product states (MPS). However, the
extension of this construction to higher dimensions is not obvious. We report
how to generalize the exponential basis to 2D and 3D by defining the basis
functions as the Green's functions of the discretized Helmholtz equation for
different Helmholtz parameters $\lambda$, a construction which is valid for
lattices of any spatial dimension. Compact tensor network representations can
then be found for the discretized Green's functions, by expressing them as
correlation functions of auxiliary fermionic fields with nearest neighbor
interactions via Grassmann Gaussian integration. Interestingly, this analytic
construction in 3D yields a $D=4$ tensor network representation of correlation
functions which (asymptotically) decay as the inverse distance ($r^{-1}_{ij}$),
thus generating the (screened) Coulomb potential on a cubic lattice. These
techniques will be useful in tensor network simulations of realistic materials.
