The is an Ornstein-Uhlebeck process with time dependent drift and you can solve this linear ODE with inhomogenous forcing term formally as (i'm going to do it with unit mass but it is straight forward to adapt to the general case)
$$
\begin{align}
v(t) = e^{-\gamma t} v_0 + \int_0^t e^{-\gamma(t-\tau)} \xi(\tau) \operatorname{d}\tau + F_0 \int_0^te^{-\gamma(t- \tau)} \cos (\omega \tau)\operatorname{d}\tau,
\end{align}
$$
this is just a Gaussian random variable - all of the stochastic behaviour coming from $\int_0^t e^{-\gamma(t-\tau)}\xi(\tau)\operatorname{d\tau}$ which being a linear transformation of a Gaussian random variable is itself a Gaussian random variable - there are existing questions justifying this fact.
Now the moment generating function of a Gaussian random variable is known to be
$$
M(\vartheta)=e^{\mu_t \vartheta + \frac{\vartheta^2\sigma^2_t }{2}},
$$
where
$$
\mu_t := \mathbb{E}\left[v(t)\right], \qquad \sigma_t^2 :=\mathbb{E}\left[ (v(t)-\mu_t)^2\right],
$$
where the subscript indicates the time dependence of these values.
Taking the logarithm of this we get the cumulant generating function and we also see that it terminates after the second order term.
It is straight forward to calculate $\mu_t$ and $\sigma_t^2$ from the explicit solution above, and this is what makes the OU process so tractable in comparison with the those cases where the Gaussian force acts multiplicatively with the state or more general nonlinear diffusions. Doing this we get
$$
\begin{align}
\mu_t &= \mathbb{E}\left[e^{-\gamma t}v_0 + \int_0^t e^{-\gamma(t-\tau)} \xi(\tau) \operatorname{d}\tau + F_0 \int_0^te^{-\gamma(t- \tau)} \cos (\omega \tau)\operatorname{d}\tau\right] \\
&= e^{-\gamma t}v_0 + F_0 \int_0^te^{-\gamma(t- \tau)} \cos (\omega \tau)\operatorname{d}\tau + \mathbb{E}\left[ \int_0^t e^{-\gamma(t-\tau)} \xi(\tau) \operatorname{d}\tau \right] \\
&= e^{-\gamma t}v_0 + F_0 \int_0^te^{-\gamma(t- \tau)} \cos (\omega \tau)\operatorname{d}\tau
\end{align}
$$
and also to calculate $\sigma_t^2$ we have the more general result that the covariance is given by
$$
\begin{align}
\operatorname{Cov}\left\{\int_0^s e^{-\gamma(t-s)}\operatorname{d}\sigma, \int_0^te^{-\gamma(t-\tau)}\operatorname{d}\tau\right\} &= \int_0^s \int_0^te^{-\gamma(t+s - (\tau + \sigma)}\langle \xi(\tau)\xi(\sigma)\rangle \operatorname{d}\tau \operatorname{d}\sigma \\
&= \frac{D}{\gamma}\left(e^{-\gamma|t-s|}-e^{-\gamma(t+s)} \right),
\end{align}
$$
and in particular
$$
\sigma_t^2 = \frac{D}{\gamma}\left(1-e^{-2\gamma t} \right).
$$
