Codifferential / divergence of differential form under conformal metric change I have a question related to this and a second post. I want to calculate the codifferential under a conformal metric change, $g_\psi = e^{2\psi} g$. By Besse's book on Einstein manifolds, or an obvious calculation, we have that (acting on $p$-forms)
$$*_{g_\psi} = e^{(n-2p)\psi}  *_g.$$
The codifferential is then given by
$$\delta^{g_\psi} = e^{-2\psi} \left( \delta - (n-2p) \iota_{\nabla \psi}\right).$$
But I don't end up with the correct factor in front of the brackets.
By definition, $\delta = *^{-1} \mathbf d \ *$. So, for some $\alpha \in \Omega^p(M)$ locally (consider a fixed multi-index for the moment),
\begin{align*}
\delta^{g_\psi} \alpha &= *_{g_\psi}^{-1} \mathbf d *_{g\psi} \left(\alpha_i \mathrm dx^i\right)\\
&= *_{g_\psi}^{-1} \mathbf d \left( e^{(n-2p)\psi} \alpha_{i} \mathrm dx^{1} \wedge ... \wedge \widehat{\mathrm dx^{i}} \wedge ... \wedge \mathrm dx^{p} \right)\\
&= *_{g_\psi}^{-1} \Big( (n-2p)e^{(n-2p)\psi} \alpha_{i} \mathrm d\psi  \wedge \mathrm dx^{1} \wedge ... \wedge \widehat{\mathrm dx^{i}} \wedge ... \wedge \mathrm dx^{p}\\
&\qquad+ e^{(n-2p)\psi} (\partial_j\alpha_i) \mathrm dx^j \wedge \mathrm dx^{1} \wedge ... \wedge \widehat{\mathrm dx^{i}} \wedge ... \wedge \mathrm dx^{p} \Big) \\
&\overset?= \underbrace{*_{g_\psi}^{-1}}_{= e^{-(n-2p)\psi} *_g} \Big[ e^{(n-2p)\psi} (n-2p) \mathrm d\psi \wedge \left(\alpha_{i} \mathrm dx^{1} \wedge ... \wedge \widehat{\mathrm dx^{i}} \wedge ... \wedge \mathrm dx^{p}\right)\\
&\qquad + (\partial_i\alpha_{i}) \mathrm dx^{1} \wedge ... \wedge \mathrm dx^{p} \Big] \\
&= -(n-2p) \ \iota_{\nabla\psi}\alpha + \delta, \\
\end{align*}
where $\widehat{\mathrm dx^{i}}$ denotes the element that was left out and using that $*(X^\flat \wedge \alpha) = (-1)^p \iota_X(*\alpha)$. 
(1) Where is the error?
(2) Is there an easy argument to deduce the formula from (cf. here)
$$\text{div}_{g_\psi} X = \text{div}_{g} X + n e^{-2\psi} X(f)?$$
 A: (1) The error indeed appears to occur at the equality marked with a $?$: If we denote $$\hat\ast := \ast_{g_\psi},$$
then in the calculation $\hat\ast^{-1}$ is applied to $d[e^{(n - 2 p) \psi} \ast \alpha]$, which is an $(n - p + 1)$-form. Thus, the inverse Hodge star operators $\ast^{-1}, \hat\ast^{-1}$ are related by
$$\hat\ast^{-1} = e^{[n - (2 (n - p + 1))] \psi} \ast^{-1} = e^{(-n + 2p - 2) \psi} \ast^{-1}.$$
Thus, the total conformal factor that appears is $$e^{(-n + 2p - 2) \psi} \cdot e^{(n - 2p) \psi} = e^{-2 \psi}$$
as claimed.
There's no need here, by the way, to resort to a coordinate computation here. Using the transformation rule for $\ast$, we can write the first few steps of the computation, again using hats $\hat\cdot$ to denote objects constructed from $\hat g := e^{2 \psi} g$, as
$$\hat\delta \alpha
= \hat\ast^{-1} d \hat\ast \alpha
= \hat\ast^{-1} d [ e^{(n - 2 p) \psi} \ast \alpha] .$$
(2) The formula for the change of the divergence of a vector field under a conformal transformation is essentially equivalent to the specialization to $p = 1$ of the formula in this question: The vector field formula gives
$$\widehat{\operatorname{div}} \,{\alpha^{\hat\sharp}}
= \operatorname{div} \alpha^{\hat\sharp} + n \alpha^{\hat\sharp} \cdot f .$$
Unwinding definitions gives that the index-raising operators $\sharp, \hat\sharp$ are related by $$\alpha^{\hat\sharp} = e^{-2 \psi} \alpha^\sharp.$$ Then, substituting, applying the Leibniz rule for $\operatorname{div}$, rearranging, and replacing $\operatorname{div} \alpha^{\sharp}$ with $\delta \alpha$ and $\widehat{\operatorname{div}} \,{\alpha^{\hat\sharp}}$ with $\hat\delta \alpha$ recovers the $p = 1$ formula.
