I think the answers to your first two questions are yes, there is always a codifferential; and no, there is not always a Hodge decomposition. I don't know anything about the Einstein manifold case.
Let's say $E$ is a vector bundle over $M$ with metric and compatible connection $\nabla$. ($E$ could be some tensor bundle with $\nabla$ induced from the Levi-Civita connection, for example.) As you say, $\nabla$ gives an exterior derivative $d^\nabla$ on $E$-valued differential forms that obeys the Leibniz rule
$$ d^\nabla (\omega \otimes e) = d\omega \otimes e + (-1)^{\text{deg} \omega} \omega \wedge \nabla e ,$$
where $\omega$ is a homogeneous form and $e$ is a section of $E$, and $\omega \wedge \nabla e$ means $\sum_i (\omega \wedge dx^i) \otimes \nabla_{\partial_i} e$.
For each degree of forms $k$, you can view $d^\nabla$ as a map $d^\nabla: C^\infty( \Lambda^k \otimes E) \to C^\infty( \Lambda^{k+1} \otimes E)$. Using the metrics, each of those spaces of sections are endowed with $L^2$ inner products. Then $d^\nabla$ always has a formal $L^2$-adjoint $\delta^\nabla$ going the other way: $\delta^\nabla: C^\infty( \Lambda^{k+1} \otimes E) \to C^\infty( \Lambda^k \otimes E)$.
(Note the space of smooth sections can be completed to the space of $L^2$ sections, but $d^\nabla$ and $\delta^\nabla$ are unbounded operators and can generally only be defined on some dense subspace of that $L^2$ space.)
So you can always define $d^\nabla$ and its formal adjoint $\delta^\nabla$. The issue is that $d^\nabla$ squares to zero if and only if the connection $\nabla$ on $E$ is flat, meaning its curvature is zero. (Flat vector bundles $E$, i.e., bundles on which there exists a flat connection, are relatively scarce.) We need $d^\nabla$ to square to zero to even attempt to define the de Rham cohomology of $E$-valued forms, and to have nice properties like the Hodge decomposition. (I want to say that the "Dirac operator" $D = d^\nabla + \delta^\nabla$ and the "Laplacian" $\Delta = D^2$ are not elliptic unless $d^\nabla$ squares to zero, but I need to think about that.)
I don't have a good reference for this at hand, although the Wikipedia article on bundle-valued forms might be useful.
I just started learning this stuff involving vector bundles, but it makes the most sense to me to simply avoid using the $\otimes$ symbol, and rather define everything in terms of $\wedge$ appropriately so that both the rules you state become equivalent.
To do so, we fix a vector bundle $(E,\pi,X)$ and for each integer $q\geq 0$, we shall view an $E$-valued $q$-form on $X$ to be a smooth map over the identity, $s: (TX)^{\oplus q}\to E$ which is fiberwise alternating and multilinear. This is of course nothing but a slight shift in perspective to considering $s$ to be a smooth section of the vector bundle $(\bigwedge^qT^*X)\otimes E$.
Now, given a $p$-form $\omega$ on $X$ (in the usual sense) we can define a wedge product $s\wedge \omega$, which will be an $E$-valued $(p+q)$-form on $X$. The definition is just as we expect: given tangent vectors $h_{1,x},\dots h_{p+q,x}\in T_xX$, we define
\begin{align}
(s\wedge \omega)(h_{1,x},\dots h_{p+q,x}):=\frac{1}{p!q!}\sum_{\sigma\in S_{p+q}}s(h_{\sigma(1),x},\cdots, h_{\sigma(p),x})\cdot \omega_x(h_{\sigma(p+1),x},\cdots, h_{\sigma(p+q),x})
\end{align}
So, this is really just the same definition of the wedge product of usual forms, except now the addition and scalar multiplication on the RHS are taking place in the vector space $E_x$. Similarly, we can define $\omega\wedge s$, and it will turn out as usual that $s\wedge \omega=(-1)^{pq}\omega\wedge s$ (just follow the usual proof, and note that if $c$ is a scalar and $v$ is a vector in any vector space, it doesn't matter whether we write $cv$ or $vc$).
So, now with this understanding of $\wedge$ (both on the left and the right), the two Leibniz rules can be written as (for a usual $p$-form $\omega$ and $E$-valued $q$-form $s$)
- $d_{D}(\omega\wedge s)=(d\omega) \wedge s+(-1)^p\omega \wedge (d_{D}s)$
- $d_{D}(s\wedge \omega)=(d_{D}s) \wedge \omega + (-1)^{q} s \wedge d\omega$
To get from one to the other simply multiply throughout by $(-1)^{pq}$ and use the "anti-commutativity" of the wedge product.
Best Answer
On 1-dimensional manifolds, this is just the holonomy (for circles) or parallel transport (for intervals). One integrates the connection with the path-ordered exponential. It seems that there is a Stokes' theorem for higher gauge theory: There is a notion of 2-holonomy of surfaces for 2-connections, see for example An Invitation to Higher Gauge Theory or Nonabelian Multiplicative Integration on Surfaces.