Solution (source code)

= Solution

Use the usual setting of a <continuously differentiable> <solenoidal vector field> defined on all radial segments from the origin, including the origin. More generally an open <star-shaped set> containing zero suffices. <Differentiation under the integral sign> gives
$$
\nabla\cdot D(x)=\int_0^1t^2(\nabla\cdot B)(tx)\,dt=0.
$$
Apply <divergence and curl of a cross product> with the second field $x$, using $\nabla\cdot x=3$ and $(D\cdot\nabla)x=D$:
$$
\nabla\times(D\times x)=3D-x(\nabla\cdot D)+(x\cdot\nabla)D-D
=2D+(x\cdot\nabla)D.
$$
The scaling identity in part (a) now converts the last expression to a one-variable derivative:
$$
2D+(x\cdot\nabla)D
=\int_0^1\left(2tB(tx)+t^2\frac d{dt}B(tx)\right)dt
=\left[t^2B(tx)\right]_0^1=B(x).
$$
Continuity at the origin ensures the lower endpoint is zero. Hence the <radial vector potential of a solenoidal vector field> is
$$
\boxed{A(x)=D(x)\times x,\qquad \nabla\times A=B.}
$$
The order of the <cross product> matters: $x\times D$ would produce $-B$. This construction is not valid on every punctured or multiply connected domain. For example $B=x/|x|^3$ is solenoidal away from zero, but its radial integral diverges at zero and it has nonzero flux through a surrounding <sphere>. The stated formula therefore implicitly needs the regular radial-domain hypothesis, rather than just local zero <divergence>.