Non-local Bi-linear term

I have been trying to understand the effects from non-locality with the following "generalization" of the harmonic potential: \[S\left[q(t), \, J(t) \right] = \int dt \left[-\frac{1}{2}m\dot{q}^{2} - q \cdot J(t)\right] + \frac{1}{2}m\omega^{2} \int dt \int ds \left[ q(t) \cdot K(t, \, s) \cdot q(s) \right] \] where the kernel is \[K_{jk}(t, \, s) = \delta_{jk}\sqrt{\frac{1}{2 \pi \epsilon}} \exp{\left(-\frac{(t-s)^{2}}{2\epsilon}\right)}\] I am interested in this kernel since it becomes the Dirac delta kernel in the limit $ \epsilon \rightarrow 0 $. The kernel $K(t, \, s)$ also appears in the Weierstrass transform of $q(t)$: \[w(s) = \int dt \left[ K(s, \, t) \cdot q(t) \right]\] Since taking $\epsilon \rightarrow 0$ leads to the Dirac delta kernel, this limit identifies the function $w$ with $q$.

The equations of motion for $ q(t) $ read \[m\ddot{q} + m\omega^{2} \int ds \left[ K(t, \, s) \cdot q(s) \right] = J(t)\] This is now an integro-differential equation. It is still a linear equation, though. So then given two linearly-independent solutions $q_{1}(t)$ and $q_{2}(t)$ we can have the linear combination of them again be a solution: \[q_{cl}(t) = A q_{1}(t) + B q_{2}(t)\] We should explicitly include the dependence on the deformation parameter $\epsilon$. This parameter has dimension of time-squared. It deforms the solution for the harmonic potential. We expect \[\lim_{\epsilon \rightarrow 0}\; q_{1}(t, \, \epsilon) = \cos{(\omega t)} \qquad \lim_{\epsilon \rightarrow 0}\; q_{2}(t, \, \epsilon) = \sin{(\omega t)}\] The kernel can be rewritten in terms of the generating function for Hermite polynomials: \[K(t, \, s) = G(t, \, s) \rho(s) \qquad G(t, \, s) = \exp{\left(\frac{ts}{\epsilon} - \frac{t^{2}}{2\epsilon}\right)} \qquad \rho(s) = \sqrt{\frac{1}{2\pi \epsilon}}\exp{\left(- \frac{s^{2}}{2\epsilon}\right)} \] The generating function $G(t, \, s)$ satisfies \[G(t, \, s) = \sum_{n = 0}^{\infty} \frac{1}{n!}\left(\frac{t}{\sqrt{2\epsilon}}\right)^{n} H_{n}\left(\frac{s}{\sqrt{2 \epsilon}}\right)\] We can introduce the inner product: \[\left\langle f | g \right\rangle = \int dt \left[\rho(t)f(t)g(t)\right]\] The non-local term can be rewritten as $\left\langle G(t, \, s) | q(s) \right\rangle$.

The equation of motion for $q(t)$ can be rewritten as \[q(t) + \omega^{2}\int dv \, \mathcal{K}(t, \, v) \cdot q(v) = \mathcal{J}(t)\] where \[\mathcal{K}(t, \, v) = \int^{t}ds \int^{s}du \, K(u, \, v) \qquad \mathcal{J}(t) = \frac{1}{m} \int^{t}ds \int^{s}du \, J(u) + q(t_{A}) + \dot{q}(t_{A}) (t - t_{A})\] In this form it can be recognized as an inhomogeneous Fredholm equation of the second type.

We can work on a finite time interval: \[t_{A} < t < t_{Z}\] Then, the integrals of the kernel can be found: \[\mathcal{K}(t, \, v) = \epsilon \left( K(t, \, v) - K(t_{A}, \, v) \right) + (t - v) \int^{t} du \, K(u, \, v)\] Note that $\mathcal{K}(t_{A}, \, v) = 0$. The solution to the Fredholm equation can be written in terms of a Liouville-Neumann series \[q(t) = \sum_{n = 0}^{\infty} (-1)^{n}\omega^{2n} q_{n}(t)\] where \[q_{0}(t) = \mathcal{J}(t)\] \[q_{1}(t) = \int dv_{1} \, \mathcal{K}(t, \, v_{1}) \cdot \mathcal{J}(v_{1})\] \[q_{2}(t) = \int dv_{2} dv_{1} \, \mathcal{K}(t, \, v_{2}) \cdot \mathcal{K}(v_{2}, \, v_{1}) \cdot \mathcal{J}(v_{1})\] \[\vdots\] \[q_{n}(t) = \int dv_{n} \ldots dv_{1} \, \mathcal{K}(t, \, v_{n}) \cdot \mathcal{K}(v_{n}, \, v_{n-1}) \cdots \mathcal{K}(v_{3}, \, v_{2}) \cdot \mathcal{K}(v_{2}, \, v_{1}) \cdot \mathcal{J}(v_{1})\] When $J = 0$ we obtain the solution to the homogeneous case.