Helmholtz 3d/4d theorem on electrodynamics is a special case of a general theorem:
1: let $A = \sum A_i dx^i $ a 1-form with $d A = \partial_k A_i dx^k\wedge dx^i =0 $ ($ R^3 : \nabla \times A=0 $) everywhere in R^n, then all path integrals from the origin 0 to any point x is independent of the path, integrals over closed paths vanish.
2: The path integral defines a potential $\Phi$ such that $A = d \Phi =\sum \partial_k \Phi dx^k$.
3: Now let $R^n$ be a space with an euclidean scalar product $< dx^i, dx^k > \text{-->} \delta^{ik}$
and let $ dV = dx^1 \wedge ... \wedge dx^n $ be the antisymmetric volume n-form.
Then one can define the metric dual form of A by taking the scalar product over the wedge product, commuting all factors to the left and contract
$ < dV, \sum A_i dx^i > == \sum < dx^1 .. dx^n , A_i dx^i > = \sum (-1)^k A_k \partial_{dx^k}dV $
where the partial derivative is short for "delete_that" in monomials.
So from the 1-form, conventionally seen in physics as a vector field by its component vector, we get a n-1-form, that forms a flow over surfaces. Its differential n-form is its production density. Again taking its dual, it is a scalar function, the sum over the second derivatives.
4:$ < dV, d < dV, A > > = \nabla . A = \Delta \Phi$
The n-form derivative $<dV, A>$ is called the divergence of A by the genereral formula
$ \nabla . B = \lim_{dV -> 0} 1/dV \int_{\partial_x^k } B . dO^k $,
that is the flow outward of the surface equals the volume integral over the production volume density. Taking the point limit, its the density. By this definition, B is not supposed to be differentiable; only the surface integrals must exist and have a point limit.
With these definitions the theorem reads:
Let A be a closed differential form in R^n everywhere: $partial_i A_k == \partial_k A_i =0$, then it is exact, ie. the derivative of a 2-times continuously differentiable function $\Phi$.
Any 2-times differentiable function $\Phi$ from $ R^n \to C$ has $\partial_{i,k}\Phi==\partial_{k,i}\Phi$. $\nabla \times A=0$ is only a very weak smoothness condition to guarantee the existence of $\Phi$.
This is important, because in the physics of Maxwell, potential functions are distributions, generally, having their sources in its singular points, lines and surfaces.
So the question reduces to: if a n-1-form over surfaces is the dual of the derivative of a 2-differentiable potential, then it is the density of a potential with $\Delta \Phi ==0$.
5: The algebraic real kernel of $\Delta$ as a differential map are the linear functions.
Liouville: A harmonic bounded function is constant.
Hilbert spaces: The dimension of the space of square integrable complex functions with $\Delta \psi=0$ is zero.
But: In distributional spaces of continous functions eg, a surface of a jump of the derivative produces a delta-density of the Laplacian.
Standard example the capacitor with
$\Phi={{ 0,x<0}, {x, x<1}, {1, x>=1} }$
5: On a flat torus, that is a unit cube in $R^n$ with opposite sides identified, all this does not apply because of toplogy: there are n linear independent constant vector fields along the n- axis.
This fact is a common case in electrodynamics: On a the superconducting torus you have two situations: The current as a a coil, field inside, or current around the axis, field outside. Both have rot A=0 div A=0, as well as rot j=0, div j=0.
- The complex case is especially important, because of its intimate relation to complex analysis.
$\partial_{(x+ I y)} f(x + I y) == f'(x+ I y)$
if $f$ is independent of $x-I y$
$ \forall x,y \partial_{(x+ I y)} f(x + I y) = 0$
Liouville: A complex differentiable bounded function is constant.
Proof: $\partial_{x+ I y,{x-Iy}} f(x+I y) =\Delta f(x+I Y) = 0$
Caveat:$\Delta Exp[ k x + I k y] ==0$ and so
$\Delta Exp[ k x + I k y] Exp[-(x^2+y^2)/s^2]/s $
can be made arbitrarily small in a very large volume