The Lienard Wiechert potential (leaving out the vector potential for simplicity),

$$\phi(\vec{r},t)=\frac{e}{4\pi\epsilon_oR(1-\hat{n}\cdot\vec{\beta})} \Bigg|_{t'=t_{ret}},$$

where the quantities should be evaluated in the retarded time given by $t_{ret}=t-R/c$, should satisfy the wave equation. The wave equation for the potential is derived from Maxwell's equations in the Lorenz gauge and given by

$$\left(\nabla^2-\frac{1}{c^2}\partial_t^2\right)\phi(\vec{r},t)=-\frac{\rho}{\epsilon_o}, $$

where $\rho$ is the charge density and should be $\rho=e\delta(\vec{r}-\vec{r}')$ for the point charge. Plugging $\phi(\vec{r},t)$ into the wave equation we need to compute the second-order spatial and time derivatives with respect to the observer. To assist with the derivation, I have computed the following spatial derivatives in cartesian coordinates ($i=x,y,$ or $z$),

$$\partial_i R=n_i, $$ $$ \partial_i \frac{1}{R} = -\frac{n_i}{R^2},$$ $$\partial_i n_i = \frac{1-n_i^2}{R},$$ $$\partial_i (1-\hat{n}\cdot\vec{\beta}) = \frac{n_i(\hat{n}\cdot\vec{\beta})-\beta_i}{R},$$

and the following observer time derivatives (taking mind to switch from $\partial_t$ to $\partial_{t'}$ using $dt=dt'(1-\hat{n}\cdot\vec{\beta})$ when necessary),

$$\partial_t R=-c\frac{\hat{n}\cdot\vec{\beta}}{(1-\hat{n}\cdot\vec{\beta})}, $$ $$ \partial_t \frac{1}{R} = c\frac{\hat{n}\cdot\vec{\beta}}{R^2(1-\hat{n}\cdot\vec{\beta})},$$ $$\partial_t n_i = \frac{c}{R(1-\hat{n}\cdot\vec{\beta})}\left({n_i(\hat{n}\cdot\vec{\beta})-\beta_i }\right),$$ $$\partial_t (1-\hat{n}\cdot\vec{\beta}) = \frac{c\left({ \vec{\beta}\cdot\vec{\beta}-(\hat{n}\cdot\vec{\beta})^2 }\right)}{R\left({1-\hat{n}\cdot\vec{\beta}}\right)}.$$

Now, when I compute the second-order derivatives I get:

$$\nabla^2\phi(\vec{r},t) = -\frac{e}{2\pi\epsilon_oR^3(1-\hat{n}\cdot\vec{\beta})^3}\left({ 1-(\hat{n}-\vec{\beta})\cdot(\hat{n}-\vec{\beta}) }\right)$$

$$\frac{1}{c^2}\partial_t^2 = \frac{e}{4\pi\epsilon_oR^3(1-\hat{n}\cdot\vec{\beta})^7}\left({ (1-\hat{n}\cdot\vec{\beta})^3 \left({(\hat{n}\cdot\vec{\beta})^2-\vec{\beta}\cdot\vec{\beta} +2(\hat{n}\cdot \vec{\beta})(\hat{n}\cdot\vec{\beta}-\vec{\beta}\cdot\vec{\beta})}\right) + (\hat{n}\cdot\vec{\beta}-\vec{\beta}\cdot\vec{\beta})((\hat{n}\cdot\vec{\beta})^2-\vec{\beta}\cdot\vec{\beta}) }\right) $$

which I dont see any way that their difference will be equal to $-\rho/\epsilon_o$ for a point charge as dictated by the wave equation. Am I doing something wrong here?