Time evolution of harmonic oscillator ground state with a linear interaction Consider the Hamiltonian:
$$\hat{H} = \frac{\hat{p}^{2}}{2m} + \frac{\hat{k} x^{2}}{2} + a \hat{x}$$
My initial state is the ground state of the free harmonic oscillator ($H_{0}$): $\psi_{0}(x,0)$
I want to find the time-evolved state, given as:
$$\psi(x,t) = e^{-i \hat{H}t}\psi_{0}(x,0) = e^{-i(\hat{H_{0}} + a\hat{x}) \ t}\psi_{0}(x,0)$$
Using the fact that:
$$(\hat{H_{0}} + a \hat{x}) \ \psi_{0}(x) = \left(\frac{\omega}{2} + ax\right) \psi_{0}(x)$$
I can now Taylor expand my exponential operator and obtain, by the above action, the following:
$$\psi(x,t) = \left[\exp\left(-i\left(ax+\frac{\omega}{2}\right)t\right)\right]\psi_0(x,0)$$
Is this correct? This doesn't seem true since, normally, this would only be the case if:
$$[\hat{H_{0}}, a \hat{x}]=0$$
 A: First of all, remember that the action of position operator $\hat{x}$ on a wavefunction is multiplication by $x$ because a wavefunction is the position representation of your state ket. The $x$ comes from the action of position operator on its eigenbra
$$\hat{x}\psi(x):=\langle x|\hat{x}|\psi\rangle=x\langle x|\psi\rangle=x\psi(x)$$
In other words, your mistake was assuming that $\psi(x)$ is an eigenfunction of position operator$^1$: multiplication by $x$ makes it a completely different function. After you acknowledge this, the next steps are clearly incorrect for they are based on wrong premises.
Coming to your problem, you can complete the square in your Hamiltonian as follows
$$\hat{H}=\frac{\hat{p}^2}{2m}+\frac{k}{2}\left(\hat{x}^2+\frac{2a}{k}\hat{x}\right)=\frac{\hat{p}^2}{2m}+\frac{k}{2}\left(\hat{x}^2+\frac{2a}{k}\hat{x}+\frac{a^2}{k^2}-\frac{a^2}{k^2}\right)=\\
=\frac{\hat{p}^2}{2m}+\frac{k}{2}\left(\hat{x}+\frac{a}{k}\right)^2-\frac{a^2}{2k}$$
that is a shifted simple harmonic harmonic oscillator.
The time evolution operator is
$$\hat{U}(t)=\exp\left[\frac{\hat{p}^2}{2m}+\frac{k}{2}\left(\hat{x}+\frac{a}{k}\right)^2\right]\exp\left(-\frac{a^2}{2k}\right)$$
The second part is irrelevant to the evolution of the system and it depends on the zero point of potential energy.
Now, your initial state is $|\phi_0\rangle$, the ground state of the harmonic oscillator centered in zero. if we label by $|\phi'_0\rangle$ the eigenstates of the above displaced harmonic oscillator, we can write the time evolution as
$$\hat{U}(t)=\sum_{n=0}^{+\infty}|\phi'_n\rangle\langle\phi'_n|\exp i\left[\omega\left(n+\frac{1}{2}\right)-\frac{a^2}{2k}\right]$$
So, you have for your state at time t
$$|\psi(t)\rangle=\hat{U}(t)|\psi(0)\rangle=\sum_{n=0}^{+\infty}|\phi'_n\rangle\exp i\left[\omega\left(n+\frac{1}{2}\right)-\frac{a^2}{2k}\right]\langle\phi'_n|\phi_0\rangle.$$
To find an explicit answer one should calculate all the $\langle\phi'_n|\phi_0\rangle$ and sum the series.

$^1$ Eigenfunctions of position operator are not even ordinary function, instead they are distributions $\langle x|x'\rangle=\delta(x-x').$
