Hodge decomposition in Minkowski space
There are no "rapid decaying harmonic 2-forms" in Minkowski space.
Consider the expression $$ 0 = \mathrm{d} \star \mathrm{d} F = \partial^i \partial_{[i}F_{jk]} = \frac13 (\Box F_{jk} + \partial_j \partial^i F_{ki} + \partial_k \partial^i F_{ij})$$ The second and third terms in the parentheses vanish by $\mathrm{d}\star F = 0$, so we see that over Minkowski space, Maxwell's equation implies that the components of the Faraday tensor must solve a linear wave equation.
(In the case you have source, you also see that the equations are equivalent to $\Box F_{jk} = J_{jk}$ where $J_{jk}$ is built out of the magnetic and electric currents in the appropriate way.)
Now, the linear decay estimates for the wave equation (that $|F| \sim \frac{1}{1+|t|}$) is sharp. In fact, by the finite speed of propagation + strong Huygens principle, if $F$ decays too rapidly in time, conservation of energy immediately implies that $F$ must vanish identically!
(If you don't want to go through the wave equation, you can also argue through the Maxwell equation by considering the energy integrals.)
Conservation of energy also implies that over Minkowski space $\mathbb{R}^{1,3}$ there cannot exist a non-vanishing "harmonic" form with finite space-time $L^2$ norm. So harmonic forms cannot exist, even in a class that is not necessary "rapidly" decaying.
Now, under the influence of external sources $j_e$ and $j_m$, you can pose the ansatz (by linearity of the equations) $F = G + H$, where $$ \begin{align} \mathrm{d} G & = j_e & \mathrm{d} H & = 0 \newline \mathrm{d} \star G & = 0 & \mathrm{d} \star H &= j_m \end{align} $$ Then the cohomology result you mentioned implies that $H = \mathrm{d}A$ and $\star G = \mathrm{d}C$ for $A,C$. The difference of the true solution with your ansatz is $\tilde{F} - F$ a harmonic two form, which can be completely absorbed into the electric potential $A$ if you want...
So the claimed decomposition is available whenever it is possible to solve the equations (possibly nonuniquely) for $G$ and $H$. For compatible source terms (you need $\mathrm{d} j_e =0 = \mathrm{d} j_m$ as a necessary condition) which are $L^1$ in space-time, the existence of a solution follows from the well-posedness of the corresponding Cauchy problem, for example.
Much of "Hodge theory" is false for "compact semi-Riemannian manifolds" in general. For example, there is no good connection between the cohomology and the number of harmonic forms. (Simplest example: on any compact Riemannian manifold the only harmonic functions are constants. But take, for example, the torus $\mathbb{T}^2 = \mathbb{S}^1\times\mathbb{S}^1$ with the Lorentzian metric $\mathrm{d}t^2 - \mathrm{d}x^2$. Then for any $k$ the function $\sin (kx) \sin(kt)$ is "harmonic", and they are linearly independent over $L^2$. )
Number 4 above also implies that the Hodge decomposition is not true for general compact semi-Riemannian manifolds. Consider the exact same $\mathbb{T}^2$ as above. Let $\alpha$ be the one form $\sin(t) \cos(x) \mathrm{d}t$. We have that $$ \mathrm{d}\alpha = \sin(t) \sin(x) \mathrm{d}t\wedge \mathrm{d}x \qquad \mathrm{d}\star \alpha = \cos(t) \cos(x) \mathrm{d}t \wedge \mathrm{d}x$$ To allow a Hodge decomposition requires that there exists $A$ and $C$, scalars, such that $$\mathrm{d}\star\mathrm{d}C = \mathrm{d}\alpha \implies \Box C = \pm \sin(t)\sin(x) $$ and $$ \mathrm{d}\star \mathrm{d}A = \mathrm{d}\star\alpha \implies \Box A = \pm \cos(t) \cos(x) $$ where $\Box$ is the wave operator $\partial_t^2 - \partial_x^2$.
Now taking Fourier transform (since we would desire $C$ and $A$ be in $L^2$ at least), we see that we have a problem. If $f\in L^2(\mathbb{T}^2)$ we have $$ \widehat{\Box f} = -(\tau^2 - \xi^2)\widehat{f} = 0 \qquad\text{ whenever } |\tau| = |\xi| $$ But the Fourier support of $\sin(t)\sin(x)$ and $\cos(t)\cos(x)$ are precisely when $\tau = \xi = 1$. So there in fact does not exist a Hodge-like decomposition for the $C^\infty$ form $\alpha$.
The above just goes to say that
- You are not going to find anything in the literature concerning Hodge decomposition for general compact semi-Riemannian manifolds, as such as theory does not exist
- You are also unlikely to find literature that studies the Hodge decomposition in the non-compact case in general: on non-compact manifolds causality violations (closed time-like curves, for example) can, in principle, also cause problems with solving the wave equation (both homogeneous and inhomogeneous), which will prevent an analogue for Hodge decomposition to hold.
- Your best bet, in so far as two-forms are concerned, are precisely the mathematical-physics literature concerning Maxwell equations over globally hyperbolic Lorentzian manifolds. For these manifolds, at least in principle, some version of the Hodge decomposition can be had by solving the appropriate initial value problems (something like $\Box C = j_e, \Box A = j_m, \mathrm{d}F_0 = 0 = \mathrm{d}\star F_0$ where $F_0$ is the harmonic part, with initial data prescribed so that the data for $C$ and $A$ vanishes on the initial surface and $F_0$ is equal to $F$ on the initial surface).
Willie's answer is of course correct. But the part that I think is most interesting in the physical context is only briefly mentioned in point 3. of the last paragraph. Let me expand on that.
The relevant condition that make everything work nicely is that the background Lorentzian spacetime be globally hyperbolic. There are two important spaces of forms, $\Omega^p_0$ and $\Omega^p_{SC}$, $p$-forms with compact and spacelike compact supports, respectively. A set $X$ is spacelike compact if it is contained in the causal influence set of a compact set $Y$, $X\subseteq J(Y)$.
On Lorentzian manifolds, the D'Alembertian operator $\square = d\delta + \delta d$ is the analog of the Laplacian operator $\Delta$ in Riemannian signature. Note that I'm using the notation $\delta = (-1)^p{\star d \star}$. It is not elliptic, but it is hyperbolic with very nice properties. Instead of the nice analytical properties of the elliptic Laplacian, we have the following exact sequence (the image of each map is equal to the kernel of the next one): $$ 0 \to \Omega^p_0 \stackrel{\square}{\to} \Omega^p_0 \stackrel{G}{\to} \Omega^p_{SC} \stackrel{\square}{\to} \Omega^p_{SC} \to 0 , $$ where $G=G^+-G^-$ is the causal Green function, defined as the difference of the retarded ($G^+$) and advanced ($G^-$) Green functions. The retarded and advanced Green functions of the D'Alembertian always exist on Globally hyperbolic manifolds (see, for instance, here or here).
A version of the result that you want is this: if $F\in \Omega^2_{SC}$, $d F = 0$ and $\delta F = 0$, then there exist two forms $A\in \Omega^1_{SC}$ and $B\in \Omega^3_{SC}$, with $\delta A = 0$ and $d B = 0$, such that $F = d A + \delta B$. Moreover, there exist $\alpha \in \Omega^1_0$ and $\beta \in \Omega^3_0$ such that $A=G\alpha$ and $B=G\beta$. This result can be found, for instance, as Prop.2.2 here. I think you can straight forwardly adapt the proof to the case with non-vanishing currents $j_e$ or $j_m$.