# Calculus/Inverting vector calculus operators

In the chapter on vector calculus, the differential operator of the gradient ($\nabla f$), the divergence ($\nabla \cdot \mathbf {F}$), and the curl ($\nabla \times \mathbf {F}$) were introduced. This chapter will focus on inverting these differential operators.
The gradient, divergence, and curl operators are all "linear", meaning that given arbitrary scalar fields $f_{1},f_{2}:\mathbb {R} ^{3}\to \mathbb {R}$, vector fields $\mathbf {F} _{1},\mathbf {F} _{2}:\mathbb {R} ^{3}\to \mathbb {R} ^{3}$, and scalars $\alpha _{1},\alpha _{2}\in \mathbb {R}$, that:
$\nabla (\alpha _{1}f_{1}+\alpha _{2}f_{2})=\alpha _{1}(\nabla f_{1})+\alpha _{2}(\nabla f_{2})$
$\nabla \cdot (\alpha _{1}\mathbf {F} _{1}+\alpha _{2}\mathbf {F} _{2})=\alpha _{1}(\nabla \cdot \mathbf {F} _{1})+\alpha _{2}(\nabla \cdot \mathbf {F} _{2})$
$\nabla \times (\alpha _{1}\mathbf {F} _{1}+\alpha _{2}\mathbf {F} _{2})=\alpha _{1}(\nabla \times \mathbf {F} _{1})+\alpha _{2}(\nabla \times \mathbf {F} _{2})$
More generally given a family of scalar fields $f[\mathbf {q} ']$, vector fields $\mathbf {F} [\mathbf {q} ']$, and coefficients $\alpha [\mathbf {q} ']$ (each $\mathbf {q} '$ corresponds to a scalar field, a vector field, and a scalar coefficient),
$\nabla \left(\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\alpha [\mathbf {q} ']f[\mathbf {q} ']dV'\right)=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\alpha [\mathbf {q} '](\nabla f[\mathbf {q} '])dV'$
$\nabla \cdot \left(\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\alpha [\mathbf {q} ']\mathbf {F} [\mathbf {q} ']dV'\right)=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\alpha [\mathbf {q} '](\nabla \cdot \mathbf {F} [\mathbf {q} '])dV'$
$\nabla \times \left(\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\alpha [\mathbf {q} ']\mathbf {F} [\mathbf {q} ']dV'\right)=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\alpha [\mathbf {q} '](\nabla \times \mathbf {F} [\mathbf {q} '])dV'$

## Inverting Linear Operators

Recall from linear algebra that when given a bijective linear operator ${\mathcal {L}}:\mathbb {R} ^{n}\to \mathbb {R} ^{n}$, that an inverse can be created by computing a solution to each of ${\mathcal {L}}(\mathbf {x} _{i})=\mathbf {e} _{i}$ where $\mathbf {e} _{i}$ is the $i^{\text{th}}$ elementary basis vector for each $i=1,2,\dots ,n$. When attempting to find an $\mathbf {x} \in \mathbb {R} ^{n}$ that solves ${\mathcal {L}}(\mathbf {x} )=\mathbf {y}$ for an arbitrary $\mathbf {y} =\sum _{i=1}^{n}y_{i}\mathbf {e} _{i}$, the solutions to ${\mathcal {L}}(\mathbf {x} _{i})=\mathbf {e} _{i}$ can be stacked in a linear manner to get $\mathbf {x} =\sum _{i=1}^{n}y_{i}\mathbf {x} _{i}$ as a possible solution. This yields an approach to inverting ${\mathcal {L}}$.
This same approach will be used to compute "Green's functions" for each vector calculus differential operator.

## The Dirac delta function

The Dirac delta function is a hypothetical function $\delta :\mathbb {R} \to \mathbb {R}$ that returns $0$ for all $x\neq 0$ and $+\infty$ for $x=0$. The Dirac delta is not meant to be evaluated at $x=0$, but instead assumed to satisfy the following integral property: $\int _{x=-\infty }^{+\infty }\delta (x)dx=1$. More generally, given any interval $[a,b]\subseteq \mathbf {R}$ that strictly contains 0: $a<0<b$, then $\int _{x=a}^{b}\delta (x)dx=1$. $\delta (x)$ is effectively a density function that describes an infinitely dense total mass of 1 at $x=0$. For an arbitrary $x'\in \mathbb {R}$, $\delta (x-x')$ describes a density function that describes an infinitely dense total mass of 1 at $x=x'$.
Even though $\delta (x)=\left\{{\begin{array}{cc}+\infty &(x=0)\\0&(x\neq 0)\end{array}}\right.$, it is not the case that $\delta (2x)=\delta (x)$, or $2\delta (x)=\delta (x)$. The integrals are: $\int _{x=-\infty }^{+\infty }\delta (2x)dx={\frac {1}{2}}\int _{x=-\infty }^{+\infty }\delta (2x)(2dx)={\frac {1}{2}}\int _{y=-\infty }^{+\infty }\delta (y)dy={\frac {1}{2}}$ and $\int _{x=-\infty }^{+\infty }2\delta (x)dx=2\int _{x=-\infty }^{+\infty }\delta (x)dx=2$. In general, for any $a\neq 0$, $\delta (ax)={\frac {1}{|a|}}\delta (x)$.

### Addressing non-differentiability

Consider the piece-wise function $f:\mathbb {R} \to \mathbb {R}$ defined by: $f(x)=\left\{{\begin{array}{cc}-1-{\frac {1}{2}}x&(x<0)\\3-x&(x\geq 0)\end{array}}\right.$. It is common to accept that $f(x)$ is not differentiable at $x=0$ and that $f'(x)=\left\{{\begin{array}{cc}-{\frac {1}{2}}&(x<0)\\-1&(x>0)\end{array}}\right.$. With the Dirac delta function, the derivative of $f$ can be expressed as $f'(x)=4\delta (x)+\left\{{\begin{array}{cc}-{\frac {1}{2}}&(x<0)\\0&(x=0)\\-1&(x>0)\end{array}}\right.$. With this derivative, part II of the fundamental theorem of calculus holds even for intervals that contain $x=0$.
In this chapter, it will generally be assumed that all scalar fields $f:\mathbb {R} ^{3}\to \mathbb {R}$ and vector fields $\mathbf {F} :\mathbb {R} ^{3}\to \mathbb {R} ^{3}$ are continuous and differentiable everywhere. However if this is not the case, the Dirac delta function will be used to model the derivative operators at points of non-differentiability.

### Inverting linear differential operators

Given an arbitrary function $f:\mathbb {R} \to \mathbb {R}$, $f$ can be expressed as the linear combination $f(x)=\int _{x'=-\infty }^{+\infty }f(x')\delta (x-x')dx'$ of Dirac delta functions. Let ${\mathcal {D}}$ be a linear differential operator that takes a real valued single variable function $f(x)$, and returns another real valued single variable function $g(x)$. If a solution $h(x;x')$ exists to ${\mathcal {D}}(h(x;x'))=\delta (x-x')$ for all $x'\in \mathbb {R}$, then given an arbitrary $g(x)$, a solution $f(x)$ exists to ${\mathcal {D}}(f(x))=g(x)$ which is $f(x)=\int _{x'=-\infty }^{+\infty }g(x')h(x;x')dx'$. The family of functions $h(x;x')$ are referred to as "Green's functions".

### Dirac delta function variants

In this chapter, given an arbitrary point $\mathbf {q} '\in \mathbb {R} ^{3}$, the function $\delta (\mathbf {q} ;\mathbf {q} ')$ will denote a 3-dimensional variant of the Dirac delta function. The key integral property is that for any volume $\Omega \subseteq \mathbb {R} ^{3}$, that $\iiint _{\mathbf {q} \in \Omega }\delta (\mathbf {q} ;\mathbf {q} ')dV=\left\{{\begin{array}{cc}1&(\mathbf {q} '\in \Omega )\\0&(\mathbf {q} '\notin \Omega )\end{array}}\right.$. Note that $\delta (\mathbf {q} ;\mathbf {q} ')=\delta (\mathbf {q} ';\mathbf {q} )$. It should also be noted that $\delta (\mathbf {q} ;\mathbf {q} ')$ carries with it the dimensions ${\text{length}}^{-3}$.
Given an arbitrary oriented curve $C$, ${\vec {\delta }}(\mathbf {q} ;C)$ will denote another variant of the Dirac delta function which returns vectors. ${\vec {\delta }}(\mathbf {q} ;C)=\mathbf {0}$ for all $\mathbf {q} \notin C$, and is infinite in the direction of $C$ for all points from $C$. The key integral property is that for any oriented surface $\sigma$, that $\iint _{\mathbf {q} \in \sigma }{\vec {\delta }}(\mathbf {q} ;C)\cdot \mathbf {dS} =N$ where $N$ is the net number of times $C$ passes through $\sigma$ in the preferred direction ($C$ passing through $\sigma$ in the reverse direction reduces $N$ by 1). It should also be noted that ${\vec {\delta }}(\mathbf {q} ;C)$ carries with it the dimensions ${\text{length}}^{-2}$.
Given an arbitrary oriented surface $\sigma$, ${\vec {\delta }}(\mathbf {q} ;\sigma )$ will denote another variant of the Dirac delta function which returns vectors. ${\vec {\delta }}(\mathbf {q} ;\sigma )=\mathbf {0}$ for all $\mathbf {q} \notin \sigma$, and is infinite in the direction of the oriented normals of $\sigma$ for all points from $\sigma$. The key integral property is that for any oriented curve $C$, that $\int _{\mathbf {q} \in C}{\vec {\delta }}(\mathbf {q} ;\sigma )\cdot d\mathbf {q} =N$ where $N$ is the net number of times $C$ passes through $\sigma$ in the preferred direction ($C$ passing through $\sigma$ in the reverse direction reduces $N$ by 1). It should also be noted that ${\vec {\delta }}(\mathbf {q} ;\sigma )$ carries with it the dimensions ${\text{length}}^{-1}$.

## Multi-paths and multi-surfaces

Given an oriented curve $C$, then $C$ can be denoted by the Dirac delta vector field ${\vec {\delta }}(\mathbf {q} ;C)$. If $C$ is continuous and starts at $\mathbf {q} _{0}$ and ends at $\mathbf {q} _{1}$, then it can be proven using Gauss's Divergence Theorem that $\nabla \cdot {\vec {\delta }}(\mathbf {q} ;C)=\delta (\mathbf {q} ;\mathbf {q} _{0})-\delta (\mathbf {q} ;\mathbf {q} _{1})$.
Given a collection $C_{1},C_{2},\dots ,C_{k}$ of oriented paths, then the vector field ${\vec {\delta }}(\mathbf {q} ;C_{1})+{\vec {\delta }}(\mathbf {q} ;C_{2})+\dots +{\vec {\delta }}(\mathbf {q} ;C_{k})$ effectively denotes the "superposition" of $C_{1},C_{2},\dots ,C_{k}$. This superposition is referred to as a "multi-path". Not all paths have to have a weight of 1. With multi-path ${\frac {1}{2}}{\vec {\delta }}(\mathbf {q} ;C_{1})+{\frac {1}{2}}{\vec {\delta }}(\mathbf {q} ;C_{2})$, the weights on $C_{1}$ and $C_{2}$ are both $0.5$. This multi-path is an even 50%/50% superposition between $C_{1}$ and $C_{2}$.
Any vector field $\mathbf {F}$ can be envisioned as a superposition of a possibly infinite number of paths. Each path may have an infinitesimal weight. When a vector field is envisioned as a multi-path, the decomposition into individual paths is not unique. When vector field $\mathbf {F}$ denotes a multi-path, $\nabla \cdot \mathbf {F}$ is the net density of path origin points minus the density of destination points: the path starting point density.
When $\nabla \cdot \mathbf {F} =0$ everywhere, $\mathbf {F}$ can be envisioned as a superposition of a possibly infinite number of paths that are either closed or extend to infinity. If it is also the case that $\mathbf {F} (\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{F}})$ for some $\alpha _{F}>2$, then all of the paths have to close and $\mathbf {F}$ is effectively a "multi-loop". ($\mathbf {F} (\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{F}})$ if and only if there exists some threshold $c_{1}>0$ and factor $c_{2}>0$ such that $\forall \mathbf {q} \in \mathbb {R} ^{3}:|\mathbf {q} |>c_{1}\implies |\mathbf {F} (\mathbf {q} )|<c_{2}(1/|\mathbf {q} |^{\alpha _{F}})$)
In the image to the right, a divergence free vector field that denotes flow density is decomposed into the superposition of multiple simple loops. 2 dimensional space is depicted as a lattice of infinitely small squares. The vector field is the top-left section. The flow along each horizontal edge, denoted by the direction and number of arrows, is the horizontal component of the vector field at the current edge (or neighboring vertex). The flow along each vertical edge, denoted by the direction and number of arrows, is the vertical component of the vector field at the current edge (or neighboring vertex). The remaining 3 sections show 3 simple loops whose superposition forms the vector field.
Given an oriented surface $\sigma$, then $\sigma$ can be denoted by the Dirac delta vector field ${\vec {\delta }}(\mathbf {q} ;\sigma )$. If $\sigma$ has the counter-clockwise oriented boundary $C$, then it can be proven using Stokes' Theorem that $\nabla \times {\vec {\delta }}(\mathbf {q} ;\sigma )={\vec {\delta }}(\mathbf {q} ;C)$.
Given a collection $\sigma _{1},\sigma _{2},\dots ,\sigma _{k}$ of oriented surfaces, then the vector field ${\vec {\delta }}(\mathbf {q} ;\sigma _{1})+{\vec {\delta }}(\mathbf {q} ;\sigma _{2})+\dots +{\vec {\delta }}(\mathbf {q} ;\sigma _{k})$ effectively denotes the "superposition" of $\sigma _{1},\sigma _{2},\dots ,\sigma _{k}$. This superposition is referred to as a "multi-surface". Not all surfaces have to have a weight of 1. With multi-surface ${\frac {1}{2}}{\vec {\delta }}(\mathbf {q} ;\sigma _{1})+{\frac {1}{2}}{\vec {\delta }}(\mathbf {q} ;\sigma _{2})$, the weights on $\sigma _{1}$ and $\sigma _{2}$ are both $0.5$. This multi-surface is an even 50%/50% superposition between $\sigma _{1}$ and $\sigma _{2}$.
Any vector field $\mathbf {F}$ can be envisioned as a superposition of a possibly infinite number of surfaces. Each surface may have an infinitesimal weight. When a vector field is envisioned as a multi-surface, the decomposition into individual surfaces is not unique. When vector field $\mathbf {F}$ denotes a multi-surface, $\nabla \times \mathbf {F}$ is the multi-loop that is the counter-clockwise oriented boundary of the multi-surface denoted by $\mathbf {F}$.
When $\nabla \times \mathbf {F} =\mathbf {0}$ everywhere, $\mathbf {F}$ can be envisioned as a superposition of a possibly infinite number of surfaces that are either closed with no boundaries or extend to infinity. If it is also the case that $\mathbf {F} (\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{F}})$ for some $\alpha _{F}>1$, then all of the surfaces have to close without extending to infinity and $\mathbf {F}$ is effectively a "multi closed surface".
Given an oriented curve $C$ and a vector field $\mathbf {F}$, then the path integral $\int _{\mathbf {q} \in C}\mathbf {F} (\mathbf {q} )\cdot d\mathbf {q}$ is equivalent to the volume integral $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;C))dV$. This statement, while intuitive to a certain degree given the definition of ${\vec {\delta }}(\mathbf {q} ;C)$, is proven in the box below:
To start, note that for all $\mathbf {q} \notin C$ that $\mathbf {F} (\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;C)=0$. Positions $\mathbf {q} \notin C$ that do not lie in $C$ have no contribution to the integral $\int _{\mathbf {q} \in C}\mathbf {F} (\mathbf {q} )\cdot d\mathbf {q}$ nor have any contribution to the integral $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;C))dV$.
Vector field ${\vec {\delta }}(\mathbf {q} ;C)$ denotes a flow density field where a total flow of 1 is moving along curve $C$. Imagining $C$ as an infinitesimally narrow "pipe", since a flow of 1 is moving through this pipe, the density of the flow is infinite. To properly model ${\vec {\delta }}(\mathbf {q} ;C)$, path $C$ will be assumed to have a non-zero infinitesimally small cross-sectional area of $a$ so the flow density is ${\frac {1}{a}}$. Path $C$ with a non-zero thickness of $a$ will be denoted by $C'(a)$. For any point $\mathbf {q} \in C'(a)$, $\mathbf {u} (\mathbf {q} )$ will denote a unit length vector that is parallel to the direction of $C$. ${\vec {\delta }}(\mathbf {q} ;C)$ will be approximated by ${\vec {\Delta }}(\mathbf {q} ;C'(a))=\left\{{\begin{array}{cc}{\frac {1}{a}}\mathbf {u} (\mathbf {q} )&(\mathbf {q} \in C'(a))\\\mathbf {0} &(\mathbf {q} \notin C'(a))\end{array}}\right.$. The volume integral becomes: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;C))dV=\lim _{a\to 0^{+}}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\Delta }}(\mathbf {q} ;C'(a)))dV$
Now consider an arbitrary infinitesimal displacement $d\mathbf {q}$ along $C$ that is located at $\mathbf {q} \in C$. The corresponding volume segment from $C'(a)$ has a volume of $dV=a|d\mathbf {q} |$. The contribution to the volume integral $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\Delta }}(\mathbf {q} ;C'(a)))dV$ is approximately $(\mathbf {F} (\mathbf {q} )\cdot {\vec {\Delta }}(\mathbf {q} ;C'(a)))dV=(\mathbf {F} (\mathbf {q} )\cdot {\frac {1}{a}}\mathbf {u} (\mathbf {q} ))(a|d\mathbf {q} |)$ $=\mathbf {F} (\mathbf {q} )\cdot (\mathbf {u} (\mathbf {q} )|d\mathbf {q} |)$ $=\mathbf {F} (\mathbf {q} )\cdot d\mathbf {q}$ which is equivalent to the contribution to the path integral $\int _{\mathbf {q} \in C}\mathbf {F} (\mathbf {q} )\cdot d\mathbf {q}$. Integrating over all path segments $d\mathbf {q}$ and letting $a\to 0^{+}$ yields: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;C))dV=\int _{\mathbf {q} \in C}\mathbf {F} (\mathbf {q} )\cdot d\mathbf {q}$
In a general sense, the infinitesimal displacement $d\mathbf {q}$ from a path integral along $C$ can be replaced by ${\vec {\delta }}(\mathbf {q} ;C)dV$ when the integral is converted to a volume integral over $\mathbb {R} ^{3}$. For example, given a continuous curve $C$ which begins at $\mathbf {q} _{0}$ and ends at $\mathbf {q} _{1}$, it can be derived that: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}{\vec {\delta }}(\mathbf {q} ;C)dV=\int _{\mathbf {q} \in C}d\mathbf {q} =\mathbf {q} _{1}-\mathbf {q} _{0}$.
Given an oriented surface $\sigma$ and a vector field $\mathbf {F}$, then the surface integral $\iint _{\mathbf {q} \in \sigma }\mathbf {F} (\mathbf {q} )\cdot \mathbf {dS}$ is equivalent to the volume integral $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;\sigma ))dV$. This statement, while intuitive to a certain degree given the definition of ${\vec {\delta }}(\mathbf {q} ;\sigma )$, is proven in the box below:
To start, note that for all $\mathbf {q} \notin \sigma$ that $\mathbf {F} (\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;\sigma )=0$. Positions $\mathbf {q} \notin \sigma$ that do not lie in $\sigma$ have no contribution to the integral $\iint _{\mathbf {q} \in \sigma }\mathbf {F} (\mathbf {q} )\cdot \mathbf {dS}$ nor have any contribution to the integral $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;\sigma ))dV$.
Vector field ${\vec {\delta }}(\mathbf {q} ;\sigma )$ denotes a vector field where a total gain of 1 is acquired when a path integral $\int _{\mathbf {q} \in C}{\vec {\delta }}(\mathbf {q} ;\sigma )\cdot d\mathbf {q}$ crosses $\sigma$ in the preferred direction. Imagining $\sigma$ as an infinitesimally thin "sheet", since a gain of 1 is acquired moving through this sheet in the preferred direction, the rate of gain is infinite. To properly model ${\vec {\delta }}(\mathbf {q} ;\sigma )$, surface $\sigma$ will be assumed to have a non-zero infinitesimally small thickness of $a$ so the rate of gain is ${\frac {1}{a}}$. Surface $\sigma$ with a non-zero thickness of $a$ will be denoted by $\sigma '(a)$. For any point $\mathbf {q} \in \sigma '(a)$, $\mathbf {n} (\mathbf {q} )$ will denote a unit length vector that is normal to $\sigma$ and points in the preferred direction. ${\vec {\delta }}(\mathbf {q} ;\sigma )$ will be approximated by ${\vec {\Delta }}(\mathbf {q} ;\sigma '(a))=\left\{{\begin{array}{cc}{\frac {1}{a}}\mathbf {n} (\mathbf {q} )&(\mathbf {q} \in \sigma '(a))\\\mathbf {0} &(\mathbf {q} \notin \sigma '(a))\end{array}}\right.$. The volume integral becomes: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;\sigma ))dV=\lim _{a\to 0^{+}}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\Delta }}(\mathbf {q} ;\sigma '(a)))dV$
Now consider an arbitrary surface vector $\mathbf {dS}$ that denotes an infinitesimal portion of $\sigma$ that is located at $\mathbf {q} \in \sigma$. The corresponding volume segment from $\sigma '(a)$ has a volume of $dV=a|\mathbf {dS} |$. The contribution to the volume integral $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\Delta }}(\mathbf {q} ;\sigma '(a)))dV$ is approximately $(\mathbf {F} (\mathbf {q} )\cdot {\vec {\Delta }}(\mathbf {q} ;\sigma '(a)))dV=(\mathbf {F} (\mathbf {q} )\cdot {\frac {1}{a}}\mathbf {n} (\mathbf {q} ))(a|\mathbf {dS} |)$ $=\mathbf {F} (\mathbf {q} )\cdot (\mathbf {n} (\mathbf {q} )|\mathbf {dS} |)$ $=\mathbf {F} (\mathbf {q} )\cdot \mathbf {dS}$ which is equivalent to the contribution to the surface integral $\iint _{\mathbf {q} \in \sigma }\mathbf {F} (\mathbf {q} )\cdot \mathbf {dS}$. Integrating over all surface segments $\mathbf {dS}$ and letting $a\to 0^{+}$ yields: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;\sigma ))dV=\iint _{\mathbf {q} \in \sigma }\mathbf {F} (\mathbf {q} )\cdot \mathbf {dS}$
In a general sense, the infinitesimal surface vector $\mathbf {dS}$ from a surface integral over $\sigma$ can be replaced by ${\vec {\delta }}(\mathbf {q} ;\sigma )dV$ when the integral is converted to a volume integral over $\mathbb {R} ^{3}$. For example, given an oriented surface $\sigma$ which has a total surface vector of $\mathbf {S}$, it can be derived that: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}{\vec {\delta }}(\mathbf {q} ;\sigma )dV=\iint _{\mathbf {q} \in \sigma }\mathbf {dS} =\mathbf {S}$.
Given a curve $C$ and a surface $\sigma$, the net number of times $C$ passes through $\sigma$ in the preferred direction, denoted by $N$, is given by:
$N=\int _{\mathbf {q} \in C}{\vec {\delta }}(\mathbf {q} ;\sigma )\cdot d\mathbf {q} =\iint _{\mathbf {q} \in \sigma }{\vec {\delta }}(\mathbf {q} ;C)\cdot \mathbf {dS} =\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}{\vec {\delta }}(\mathbf {q} ;C)\cdot {\vec {\delta }}(\mathbf {q} ;\sigma )dV$
Given vector fields $\mathbf {F}$ and $\mathbf {G}$ where $\mathbf {F}$ denotes a multi-path and $\mathbf {G}$ denotes a multi-surface, then $N=\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot \mathbf {G} (\mathbf {q} ))dV$ denotes the total "flux" of multi-path $\mathbf {F}$ through multi-surface $\mathbf {G}$.
If $C$ is closed loop and $\sigma$ is a closed surface then $N=0$ since every time $C$ passes through $\sigma$, $C$ must pass through $\sigma$ in the opposite direction in order to close itself. More generally, if vector field $\mathbf {F}$ denotes a multi-loop, which means that $\nabla \cdot \mathbf {F} =0$ and $\mathbf {F} (\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{F}})$ for some $\alpha _{F}>2$, and vector field $\mathbf {G}$ denotes a multi closed surface, which means that $\nabla \times \mathbf {G} =\mathbf {0}$ and $\mathbf {G} (\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{G}})$ for some $\alpha _{G}>1$, then $N=\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot \mathbf {G} (\mathbf {q} ))dV=0$. An algebraic proof (which is not truly necessary) is given in the box below:
To start, let $\mathbf {F} (\mathbf {q} )=F_{x}(\mathbf {q} )\mathbf {i} +F_{y}(\mathbf {q} )\mathbf {j} +F_{z}(\mathbf {q} )\mathbf {k}$ and $\mathbf {G} (\mathbf {q} )=G_{x}(\mathbf {q} )\mathbf {i} +G_{y}(\mathbf {q} )\mathbf {j} +G_{z}(\mathbf {q} )\mathbf {k}$.
$\nabla \cdot \mathbf {F} =0\implies {\frac {\partial F_{z}}{\partial z}}=-\left({\frac {\partial F_{x}}{\partial x}}+{\frac {\partial F_{y}}{\partial y}}\right)$
$\implies F_{z}(x,y,z)=-\int _{z'=-\infty }^{z}\left({\frac {\partial F_{x}}{\partial x}}{\bigg |}_{(x,y,z')}+{\frac {\partial F_{y}}{\partial y}}{\bigg |}_{(x,y,z')}\right)dz'$
$\implies F_{z}=-({\frac {\partial I_{F,x}}{\partial x}}+{\frac {\partial I_{F,y}}{\partial y}})$ where $I_{F,x}(x,y,z)=\int _{z'=-\infty }^{z}F_{x}(x,y,z')dz'$ and $I_{F,y}(x,y,z)=\int _{z'=-\infty }^{z}F_{y}(x,y,z')dz'$
$\nabla \times \mathbf {G} =\mathbf {0} \implies {\frac {\partial G_{x}}{\partial z}}={\frac {\partial G_{z}}{\partial x}}\quad {\text{and}}\quad {\frac {\partial G_{y}}{\partial z}}={\frac {\partial G_{z}}{\partial y}}$
$\implies G_{x}(x,y,z)=\int _{z'=-\infty }^{z}{\frac {\partial G_{z}}{\partial x}}{\bigg |}_{(x,y,z')}dz'\quad {\text{and}}\quad G_{y}(x,y,z)=\int _{z'=-\infty }^{z}{\frac {\partial G_{z}}{\partial y}}{\bigg |}_{(x,y,z')}dz'$
$\implies G_{x}={\frac {\partial I_{G,z}}{\partial x}}\quad {\text{and}}\quad G_{y}={\frac {\partial I_{G,z}}{\partial y}}$ where $I_{G,z}(x,y,z)=\int _{z'=-\infty }^{z}G_{z}(x,y,z')dz'$
It is important to note that:
    As $x\to \pm \infty$, $y\to \pm \infty$, or $z\to \pm \infty$, that $\mathbf {F} ,\mathbf {G} \to \mathbf {0}$.
    As $x\to \pm \infty$, $y\to \pm \infty$, or $z\to -\infty$, that $I_{F,x},I_{F,y},I_{G,z}\to 0$.
    Since $F_{z},G_{x},G_{y}\to 0$ as $z\to \pm \infty$, then ${\frac {\partial I_{F,x}}{\partial x}}+{\frac {\partial I_{F,y}}{\partial y}},{\frac {\partial I_{G,z}}{\partial x}},{\frac {\partial I_{G,z}}{\partial y}}\to 0$ as $z\to \pm \infty$.
Evaluating the integral gives:
$\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {F} (\mathbf {q} )\cdot \mathbf {G} (\mathbf {q} ))dV=\iiint _{(x,y,z)\in \mathbb {R} ^{3}}(F_{x}G_{x}+F_{y}G_{y}+F_{z}G_{z})dxdydz$
$=\iiint _{(x,y,z)\in \mathbb {R} ^{3}}{\biggF_{x}{\frac {\partial I_{G,z}}{\partial x}}+F_{y}{\frac {\partial I_{G,z}}{\partial y}}-G_{z}({\frac {\partial I_{F,x}}{\partial x}}+{\frac {\partial I_{F,y}}{\partial y}}){\bigg )}dxdydz$
$=\iiint _{(x,y,z)\in \mathbb {R} ^{3}}{\bigg{\frac {\partial I_{F,x}}{\partial z}}{\frac {\partial I_{G,z}}{\partial x}}+{\frac {\partial I_{F,y}}{\partial z}}{\frac {\partial I_{G,z}}{\partial y}}-{\frac {\partial I_{G,z}}{\partial z}}({\frac {\partial I_{F,x}}{\partial x}}+{\frac {\partial I_{F,y}}{\partial y}}){\bigg )}dxdydz$
$=\iiint _{(x,y,z)\in \mathbb {R} ^{3}}{\bigg({\frac {\partial I_{F,x}}{\partial z}}{\frac {\partial I_{G,z}}{\partial x}}-{\frac {\partial I_{G,z}}{\partial z}}{\frac {\partial I_{F,x}}{\partial x}})+({\frac {\partial I_{F,y}}{\partial z}}{\frac {\partial I_{G,z}}{\partial y}}-{\frac {\partial I_{G,z}}{\partial z}}{\frac {\partial I_{F,y}}{\partial y}}){\bigg )}dxdydz$
$=\iiint _{(x,y,z)\in \mathbb {R} ^{3}}{\bigg({\frac {\partial }{\partial z}}(I_{F,x}{\frac {\partial I_{G,z}}{\partial x}})-{\frac {\partial }{\partial x}}(I_{F,x}{\frac {\partial I_{G,z}}{\partial z}}))+({\frac {\partial }{\partial z}}(I_{F,y}{\frac {\partial I_{G,z}}{\partial y}})-{\frac {\partial }{\partial y}}(I_{F,y}{\frac {\partial I_{G,z}}{\partial z}})){\bigg )}dxdydz$
$=\iiint _{(x,y,z)\in \mathbb {R} ^{3}}{\bigg{\frac {\partial }{\partial z}}(I_{F,x}{\frac {\partial I_{G,z}}{\partial x}}+I_{F,y}{\frac {\partial I_{G,z}}{\partial y}})-{\frac {\partial }{\partial x}}(I_{F,x}G_{z})-{\frac {\partial }{\partial y}}(I_{F,y}G_{z}){\bigg )}dxdydz$
$=\iint _{(x,y)\in \mathbb {R} ^{2}}(I_{F,x}{\frac {\partial I_{G,z}}{\partial x}}+I_{F,y}{\frac {\partial I_{G,z}}{\partial y}}){\bigg |}_{z=-\infty }^{+\infty }dxdy-\iint _{(y,z)\in \mathbb {R} ^{2}}I_{F,x}G_{z}{\bigg |}_{x=-\infty }^{+\infty }dydz-\iint _{(x,z)\in \mathbb {R} ^{2}}I_{F,y}G_{z}{\bigg |}_{y=-\infty }^{+\infty }dxdz$
$=0-0-0=0$
An example application of multi-paths and multi-surfaces is given in the box below:
Given a closed loop $C$ that is carrying an electric current of $I$, the magnetic field $\mathbf {B}$ generated by this closed loop obeys the following equations:
$\nabla \cdot \mathbf {B} =0$
$\nabla \times \mathbf {B} =\mu _{0}I{\vec {\delta }}(\mathbf {q} ;C)$ (Ampere's Law in magnetostatics)
where $\mu _{0}=4\pi \times 10^{-7}{\text{N}}/{\text{A}}^{2}$ is magnetic permeability of free space. It is also required that $\mathbf {B} (\mathbf {q} )\to \mathbf {0}$ as $|\mathbf {q} |\to +\infty$.
Let $C_{1}$ and $C_{2}$ be two closed loops which are the counterclockwise boundaries of surfaces $\sigma _{1}$ and $\sigma _{2}$. Let $\mathbf {B} _{1}$ denote the magnetic field generated by running an electric current of $I_{1}$ around $C_{1}$. Note that $\mathbf {B} _{1}$ is proportional to $I_{1}$. Now let $\Phi _{B,2}$ denote the total flux of $\mathbf {B} _{1}$ through $\sigma _{2}$ in the preferred direction: $\Phi _{B,2}=\iint _{\mathbf {q} \in \sigma _{2}}\mathbf {B} _{1}(\mathbf {q} )\cdot \mathbf {dS}$. Note that since $\nabla \cdot \mathbf {B} _{1}=0$, that $\Phi _{B,2}$ is a function of $C_{2}$ instead of $\sigma _{2}$. It is the case that $\Phi _{B,2}$ is proportional to $I_{1}$, and the constant of proportionality $M(C_{1},C_{2})={\frac {\Phi _{B,2}}{I_{1}}}$ is the "mutual inductance" from $C_{1}$ to $C_{2}$.
The mutual inductance $M(C_{1},C_{2})$ is purely a function of $C_{1}$ and $C_{2}$. It is also the case that the mutual inductance is symmetric: $M(C_{1},C_{2})=M(C_{2},C_{1})$. The symmetry of the mutual inductance is not obvious, and while the symmetry is apparent from explicit formulas for the mutual inductance such as the "Neumann formula", the symmetry can be made clear by interpreting the magnetic field as both a multi-loop and a multi-surface.
The vector field $\mathbf {b} _{1}={\frac {\mathbf {B} _{1}}{\mu _{0}I_{1}}}$ satisfies both $\nabla \cdot \mathbf {b} _{1}=0$ and $\nabla \times \mathbf {b} _{1}={\vec {\delta }}(\mathbf {q} ;C_{1})$. While not proven here, the magnetic field generated by a closed loop of current that does not extend to infinity is $O(1/|\mathbf {q} |^{3})$ so $\mathbf {b} _{1}$ is both a multi-loop and a multi-surface with counterclockwise boundary $C_{1}$. Doing the same for loop $C_{2}$ gives $\mathbf {b} _{2}={\frac {\mathbf {B} _{2}}{\mu _{0}I_{2}}}$ where $\mathbf {B} _{2}$ is the magnetic field generated solely by running a current of $I_{2}$ around $C_{2}$. $\mathbf {b} _{2}$ is both a multi-loop and a multi-surface with counterclockwise boundary $C_{2}$.
The flux of $\mathbf {B} _{1}$ through $\sigma _{2}$ is:
$\Phi _{B,2}=\iint _{\mathbf {q} \in \sigma _{2}}\mathbf {B} _{1}(\mathbf {q} )\cdot \mathbf {dS}$ $=I_{1}\mu _{0}\iint _{\mathbf {q} \in \sigma _{2}}\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {dS}$ $=I_{1}\mu _{0}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;\sigma _{2}))dV$ $=I_{1}\mu _{0}\left(\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {b} _{2}(\mathbf {q} ))dV+\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot ({\vec {\delta }}(\mathbf {q} ;\sigma _{2})-\mathbf {b} _{2}(\mathbf {q} )))dV\right)$
Since $\nabla \times ({\vec {\delta }}(\mathbf {q} ;\sigma _{2})-\mathbf {b} _{2}(\mathbf {q} ))={\vec {\delta }}(\mathbf {q} ;C_{2})-{\vec {\delta }}(\mathbf {q} ;C_{2})=\mathbf {0}$, the vector field ${\vec {\delta }}(\mathbf {q} ;\sigma _{2})-\mathbf {b} _{2}(\mathbf {q} )$ is a multi closed surface, and the total flux of multi-loop $\mathbf {b} _{1}$ through this multi closed surface is 0: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot ({\vec {\delta }}(\mathbf {q} ;\sigma _{2})-\mathbf {b} _{2}(\mathbf {q} )))dV=0$.
Therefore $\Phi _{B,2}=I_{1}\mu _{0}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {b} _{2}(\mathbf {q} ))dV$ so the mutual inductance is $M(C_{1},C_{2})=\mu _{0}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {b} _{2}(\mathbf {q} ))dV$. The symmetry of the mutual inductance is now apparent from this formula: $M(C_{1},C_{2})=\mu _{0}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {b} _{2}(\mathbf {q} ))dV=\mu _{0}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{2}(\mathbf {q} )\cdot \mathbf {b} _{1}(\mathbf {q} ))dV=M(C_{2},C_{1})$.
Given a closed loop $C$ that is carrying an electric current of $I$, the magnetic field $\mathbf {B}$ and the auxiliary field $\mathbf {H}$ generated by this closed loop obeys the following equations:
$\nabla \cdot \mathbf {B} =0$
$\nabla \times \mathbf {H} =I{\vec {\delta }}(\mathbf {q} ;C)$ (Ampere's Law in a magnetic material)
$\mathbf {H} ={\frac {1}{\mu }}\mathbf {B}$ (Definition of the auxiliary field $\mathbf {H}$)
where $\mu$ is the magnetic permeability of the magnetic material. $\mu$ may vary over space. It is also required that $\mathbf {B} (\mathbf {q} )\to \mathbf {0}$ as $|\mathbf {q} |\to +\infty$.
Let $C_{1}$ and $C_{2}$ be two closed loops which are the counterclockwise boundaries of surfaces $\sigma _{1}$ and $\sigma _{2}$. Let $\mathbf {B} _{1}$ and $\mathbf {H} _{1}$ respectively denote the magnetic field and auxiliary field generated by running an electric current of $I_{1}$ around $C_{1}$. Note that both $\mathbf {B} _{1}$ and $\mathbf {H} _{1}$ are proportional to $I_{1}$. Now let $\Phi _{B,2}$ denote the total flux of $\mathbf {B} _{1}$ through $\sigma _{2}$ in the preferred direction: $\Phi _{B,2}=\iint _{\mathbf {q} \in \sigma _{2}}\mathbf {B} _{1}(\mathbf {q} )\cdot \mathbf {dS}$. Note that since $\nabla \cdot \mathbf {B} _{1}=0$ ($\mathbf {B} _{1}$ is a multi-loop), that $\Phi _{B,2}$ is a function of $C_{2}$ instead of $\sigma _{2}$. It is the case that $\Phi _{B,2}$ is proportional to $I_{1}$, and the constant of proportionality $M(C_{1},C_{2})={\frac {\Phi _{B,2}}{I_{1}}}$ is the "mutual inductance" from $C_{1}$ to $C_{2}$.
The mutual inductance $M(C_{1},C_{2})$ is purely a function of $C_{1}$ and $C_{2}$. It is also the case that the mutual inductance is symmetric: $M(C_{1},C_{2})=M(C_{2},C_{1})$. The symmetry of the mutual inductance is not obvious, and can be made clear by interpreting the magnetic field as a multi-loop and the auxiliary field as a multi-surface.
The vector fields $\mathbf {b} _{1}={\frac {\mathbf {B} _{1}}{I_{1}}}$ and $\mathbf {h} _{1}={\frac {\mathbf {H} _{1}}{I_{1}}}$ satisfies $\nabla \cdot \mathbf {b} _{1}=0$, $\nabla \times \mathbf {h} _{1}={\vec {\delta }}(\mathbf {q} ;C_{1})$, and $\mathbf {h} _{1}={\frac {1}{\mu }}\mathbf {b} _{1}$. While not proven here, the magnetic field generated by a closed loop of current that does not extend to infinity is $O(1/|\mathbf {q} |^{3})$ so $\mathbf {b} _{1}$ is a multi-loop and $\mathbf {h} _{1}$ is a multi-surface with counterclockwise boundary $C_{1}$. Doing the same for loop $C_{2}$ gives $\mathbf {b} _{2}={\frac {\mathbf {B} _{2}}{I_{2}}}$ and $\mathbf {h} _{2}={\frac {\mathbf {H} _{2}}{I_{2}}}$. $\mathbf {B} _{2}$ is the magnetic field and $\mathbf {H} _{2}$ is the auxiliary field generated solely by running a current of $I_{2}$ around $C_{2}$. $\mathbf {b} _{2}$ is a multi-loop and $\mathbf {h} _{2}$ is a multi-surface with counterclockwise boundary $C_{2}$.
The flux of $\mathbf {B} _{1}$ through $\sigma _{2}$ is:
$\Phi _{B,2}=\iint _{\mathbf {q} \in \sigma _{2}}\mathbf {B} _{1}(\mathbf {q} )\cdot \mathbf {dS}$ $=I_{1}\iint _{\mathbf {q} \in \sigma _{2}}\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {dS}$ $=I_{1}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot {\vec {\delta }}(\mathbf {q} ;\sigma _{2}))dV$ $=I_{1}\left(\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {h} _{2}(\mathbf {q} ))dV+\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot ({\vec {\delta }}(\mathbf {q} ;\sigma _{2})-\mathbf {h} _{2}(\mathbf {q} )))dV\right)$
Since $\nabla \times ({\vec {\delta }}(\mathbf {q} ;\sigma _{2})-\mathbf {h} _{2}(\mathbf {q} ))={\vec {\delta }}(\mathbf {q} ;C_{2})-{\vec {\delta }}(\mathbf {q} ;C_{2})=\mathbf {0}$, the vector field ${\vec {\delta }}(\mathbf {q} ;\sigma _{2})-\mathbf {h} _{2}(\mathbf {q} )$ is a multi-closed surface, and the total flux of multi-loop $\mathbf {b} _{1}$ through this multi-closed surface is 0: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot ({\vec {\delta }}(\mathbf {q} ;\sigma _{2})-\mathbf {h} _{2}(\mathbf {q} )))dV=0$.
Therefore $\Phi _{B,2}=I_{1}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {h} _{2}(\mathbf {q} ))dV$ $=I_{1}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}{\frac {1}{\mu (\mathbf {q} )}}(\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {b} _{2}(\mathbf {q} ))dV$ so the mutual inductance is $M(C_{1},C_{2})=\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}{\frac {1}{\mu (\mathbf {q} )}}(\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {b} _{2}(\mathbf {q} ))dV$. The symmetry of the mutual inductance is now apparent from this formula: $M(C_{1},C_{2})=\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}{\frac {1}{\mu (\mathbf {q} )}}(\mathbf {b} _{1}(\mathbf {q} )\cdot \mathbf {b} _{2}(\mathbf {q} ))dV=\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}{\frac {1}{\mu (\mathbf {q} )}}(\mathbf {b} _{2}(\mathbf {q} )\cdot \mathbf {b} _{1}(\mathbf {q} ))dV=M(C_{2},C_{1})$.
Given a continuous oriented curve $C$ that originates from $\mathbf {q} _{0}$ and terminates at $\mathbf {q} _{1}$, it has been noted that $\nabla \cdot {\vec {\delta }}(\mathbf {q} ;C)=\delta (\mathbf {q} ;\mathbf {q} _{0})-\delta (\mathbf {q} ;\mathbf {q} _{1})$. It can be derived that: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}{\vec {\delta }}(\mathbf {q} ;C)dV=\int _{\mathbf {q} \in C}d\mathbf {q}$ $=\mathbf {q} _{1}-\mathbf {q} _{0}$ $=\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\delta (\mathbf {q} ;\mathbf {q} _{1})-\delta (\mathbf {q} ;\mathbf {q} _{0}))\mathbf {q} dV$ $=-\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\nabla \cdot {\vec {\delta }}(\mathbf {q} ;C))|_{\mathbf {q} }\mathbf {q} dV$. Generalizing to a multi-path vector field $\mathbf {F}$, it is the case that: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}\mathbf {F} (\mathbf {q} )dV=-\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}(\nabla \cdot {\vec {\mathbf {F} }})|_{\mathbf {q} }\mathbf {q} dV$. Recall that $\nabla \cdot \mathbf {F}$ is a measure of the "starting point density" of $\mathbf {F}$.
Given an oriented surface $\sigma$ with counterclockwise boundary $C$, it has been noted that $\nabla \times {\vec {\delta }}(\mathbf {q} ;\sigma )=\delta (\mathbf {q} ;C)$. It can be derived that: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}{\vec {\delta }}(\mathbf {q} ;\sigma )dV=\iint _{\mathbf {q} \in \sigma }\mathbf {dS}$ $={\frac {1}{2}}\int _{\mathbf {q} \in C}\mathbf {q} \times d\mathbf {q}$ $={\frac {1}{2}}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}\mathbf {q} \times ({\vec {\delta }}(\mathbf {q} ;C)dV)$ $={\frac {1}{2}}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}\mathbf {q} \times (\nabla \times {\vec {\delta }}(\mathbf {q} ;\sigma ))|_{\mathbf {q} }dV$. Generalizing to a multi-surface vector field $\mathbf {F}$, it is the case that: $\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}\mathbf {F} (\mathbf {q} )dV={\frac {1}{2}}\iiint _{\mathbf {q} \in \mathbb {R} ^{3}}\mathbf {q} \times (\nabla \times {\vec {\mathbf {F} }})|_{\mathbf {q} }dV$. Recall that $\nabla \times \mathbf {F}$ is the counterclockwise boundary density of $\mathbf {F}$.

## Inverting the divergence operator (the inverse square law)

Given an arbitrary scalar field $f:\mathbb {R} ^{3}\to \mathbb {R}$, the problem of interest is that of finding a vector field $\mathbf {F} :\mathbb {R} ^{3}\to \mathbb {R} ^{3}$ that satisfies $\nabla \cdot \mathbf {F} =f$ and $\nabla \times \mathbf {F} =\mathbf {0}$. In other words, we want to find a vector field whose divergence is given by $f$ and is irrotational. For reasons that will soon become apparent, it will be assumed that $f(\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{f}})$ as $|\mathbf {q} |\rightarrow +\infty$ for some $\alpha _{f}>1$. The "big O" means that $\exists c_{1},c_{2}>0:\forall \mathbf {q} \in \mathbb {R} ^{3}:|\mathbf {q} |>c_{1}\implies |f(\mathbf {q} )|<c_{2}(1/|\mathbf {q} |^{\alpha _{f}})$.
$f$ can be expressed as a linear combination of Dirac delta functions: $f(\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}f(\mathbf {q} ')\delta (\mathbf {q} ;\mathbf {q} ')dV'$.
If an irrotational vector field $\mathbf {H}$ can be determined such that $(\nabla \cdot \mathbf {H} )|_{\mathbf {q} }=\delta (\mathbf {q} ;\mathbf {0} )$ (that is divergence-free everywhere except the origin where the divergence is infinite), then for all $\mathbf {q} '\in \mathbb {R} ^{3}$: $(\nabla \cdot \mathbf {H} (\mathbf {q} -\mathbf {q} '))|_{\mathbf {q} }=\delta (\mathbf {q} ;\mathbf {q} ')$. Since the divergence operator distributes over linear combinations, $\mathbf {F} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}f(\mathbf {q} ')\mathbf {H} (\mathbf {q} -\mathbf {q} ')dV'$ is an irrotational vector field that satisfies $\nabla \cdot \mathbf {F} =f$.
An intuitive candidate for an irrotational vector field that is divergence free everywhere except the origin is a radially symmetric vector field $\mathbf {H} (\mathbf {q} )={\frac {k}{r^{2}}}{\hat {\mathbf {r} }}$ where $r=|\mathbf {q} |$ is the distance from the origin, and ${\hat {\mathbf {r} }}={\frac {\mathbf {q} }{|\mathbf {q} |}}$ is the unit vector that points away from the origin (see spherical coordinates). $k$ is unknown at this point. The inverse square indicates that the flow diffuses out over a larger area as the distance from the origin increases. It can easily be checked that $\mathbf {H}$ is irrotational everywhere (including the origin), and is divergence free everywhere except the origin. All that remains is to determine $k$ such that $(\nabla \cdot \mathbf {H} )|_{\mathbf {q} }=\delta (\mathbf {q} ;\mathbf {0} )$. Consider a sphere of radius $R$ centered on the origin. The total outwards flow/flux through the surface of this sphere is ${\frac {k}{R^{2}}}\cdot (4\pi R^{2})=4\pi k$, so the total flow generated inside the sphere is $4\pi k$. Since $(\nabla \cdot \mathbf {H} )|_{\mathbf {q} }=\delta (\mathbf {q} ;\mathbf {0} )$, the total flow generated inside the sphere is $1$. This gives $k={\frac {1}{4\pi }}$. Therefore $\mathbf {H} (\mathbf {q} )={\frac {1}{4\pi r^{2}}}{\hat {\mathbf {r} }}={\frac {\mathbf {q} }{4\pi |\mathbf {q} |^{3}}}$.
In total, $\mathbf {F} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}f(\mathbf {q} '){\frac {\mathbf {q} -\mathbf {q} '}{4\pi |\mathbf {q} -\mathbf {q} '|^{3}}}dV'$ is an irrotational vector field that satisfies $\nabla \cdot \mathbf {F} =f$. The assumption that $f(\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{f}})$ as $|\mathbf {q} |\rightarrow +\infty$ for some $\alpha _{f}>1$ implies that the volume integral does not diverge to infinity. Also note that the apparent singularity at $\mathbf {q} '=\mathbf {q}$ does not impact the integral (to be discussed below).

### Uniqueness

An important question is if the vector field $\mathbf {F} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}f(\mathbf {q} '){\frac {\mathbf {q} -\mathbf {q} '}{4\pi |\mathbf {q} -\mathbf {q} '|^{3}}}dV'$ is the only irrotational vector field that satisfies $\nabla \cdot \mathbf {F} =f$. If $\mathbf {F} '$ is another possible solution, then $\mathbf {F} ''=\mathbf {F} '-\mathbf {F}$ is a vector field that is both irrotational and divergence free at all points: $\nabla \times \mathbf {F} ''=\mathbf {0}$ and $\nabla \cdot \mathbf {F} ''=0$. There is then the following theorem:

#### Divergence free and irrotational vector fields

Given a vector field $\mathbf {F}$ that is divergence free (meaning: For all closed oriented surfaces $\sigma$: $\iint _{\mathbf {q} \in \sigma }\mathbf {F} (\mathbf {q} )\cdot \mathbf {dS} =0$)
and irrotational, (meaning: For all closed oriented curves $C$: $\int _{\mathbf {q} \in C}\mathbf {F} (\mathbf {q} )\cdot d\mathbf {q} =0$)
it must then be the case that $\mathbf {F}$ is everywhere constant: $\mathbf {F} (\mathbf {q} )=\mathbf {C}$
$\mathbf {F} ''=\mathbf {F} '-\mathbf {F}$ is a constant vector $\mathbf {C}$, and therefore an irrotational vector field $\mathbf {F}$ that satisfies $\nabla \cdot \mathbf {F} =f$ is unique up to the addition of a constant vector field.

### About Improper Integrals

The volume integral $\mathbf {F} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}f(\mathbf {q} '){\frac {\mathbf {q} -\mathbf {q} '}{4\pi |\mathbf {q} -\mathbf {q} '|^{3}}}dV'$ has a pole/singularity when $\mathbf {q} '=\mathbf {q}$. Its range also extends to infinity. Both of these irregularities have the potential to result in a divergent (infinite) integral.
To analyse whether or not the integral diverges due to the pole/singularity or infinite range, the volume integral will be expressed as the integral of concentric spherical shells centered on $\mathbf {q}$:
$\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}f(\mathbf {q} '){\frac {\mathbf {q} -\mathbf {q} '}{4\pi |\mathbf {q} -\mathbf {q} '|^{3}}}dV'=\int _{r=0}^{+\infty }\iint _{\mathbf {q} '\in b(\mathbf {q} ;r)}f(\mathbf {q} '){\frac {(\mathbf {q} -\mathbf {q} ')/r}{4\pi r^{2}}}|\mathbf {dS} '|dr$ where $b(\mathbf {q} ;r)=\{\mathbf {q} ':|\mathbf {q} '-\mathbf {q} |=r\}$ denotes the outwards oriented surface of a sphere centered on $\mathbf {q}$ with radius $r$. $B(\mathbf {q} ;r)=\{\mathbf {q} ':|\mathbf {q} '-\mathbf {q} |<r\}$ denotes a solid sphere centered on $\mathbf {q}$ with radius $r$.
The inner surface integral does not present any irregularities. The lower bound of $r=0$ denotes the pole, and the upper bound of $r=+\infty$ denotes the infinite range.
The integral $\mathbf {F} (\mathbf {q} ;\epsilon ;R)=\int _{r=\epsilon }^{+R}\iint _{\mathbf {q} '\in b(\mathbf {q} ;r)}f(\mathbf {q} '){\frac {(\mathbf {q} -\mathbf {q} ')/r}{4\pi r^{2}}}|\mathbf {dS} '|dr$ where the inner radius $\epsilon$ is small and the outer radius $R$ is large is regular. The goal is to analyse the integral's behavior as $\epsilon \rightarrow 0^{+}$ and $R\rightarrow +\infty$.
Assume that $|f(\mathbf {q} ')|<K(\mathbf {q} ;r)$ for all $\mathbf {q} '\in b(\mathbf {q} ;r)$ for each radius $r$. The magnitude of the inner surface integral is now bounded from above by:
$\left|\iint _{\mathbf {q} '\in b(\mathbf {q} ;r)}f(\mathbf {q} '){\frac {(\mathbf {q} -\mathbf {q} ')/r}{4\pi r^{2}}}|\mathbf {dS} '|\right|\leq \iint _{\mathbf {q} '\in b(\mathbf {q} ;r)}|f(\mathbf {q} ')|{\frac {|(\mathbf {q} -\mathbf {q} ')/r|}{4\pi r^{2}}}|\mathbf {dS} '|\leq \iint _{\mathbf {q} '\in b(\mathbf {q} ;r)}K(\mathbf {q} ;r){\frac {1}{4\pi r^{2}}}|\mathbf {dS} '|$
$\quad \quad =K(\mathbf {q} ;r){\frac {1}{4\pi r^{2}}}\iint _{\mathbf {q} '\in b(\mathbf {q} ;r)}|\mathbf {dS} '|=K(\mathbf {q} ;r)$
As can be seen, the surface area of $b(\mathbf {q} ;r)$ cancels out ${\frac {1}{4\pi r^{2}}}$.
The magnitude of $\mathbf {F} (\mathbf {q} ;\epsilon ;R)$ is bounded from above by:
$|\mathbf {F} (\mathbf {q} ;\epsilon ;R)|\leq \int _{r=\epsilon }^{+R}K(\mathbf {q} ;r)dr$
Assuming that $f(\mathbf {q} ')$ is bounded everywhere, $K(\mathbf {q} ;r)$ is bounded everywhere, so it is clear that $|\mathbf {F} (\mathbf {q} ;\epsilon ;R)|$ does not approach infinity as $\epsilon \rightarrow 0^{+}$. This settles the pole/singularity at $\mathbf {q} '=\mathbf {q}$. For the integral to not diverge as $R\rightarrow +\infty$, the condition that $f(\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{f}})$ as $|\mathbf {q} |\rightarrow +\infty$ for some $\alpha _{f}>1$ is sufficient. This condition implies that $K(\mathbf {q} ;r)$ is less than a multiple of $1/r^{\alpha _{f}}$ when $r$ becomes sufficiently large. The integral of $1/r^{\alpha _{f}}$ converges at infinity provided that $\alpha _{f}>1$.

## Inverting the gradient operator

Given a vector field $\mathbf {F} :\mathbb {R} ^{3}\to \mathbb {R} ^{3}$, the problem of interest is finding a scalar field $f:\mathbb {R} ^{3}\to \mathbb {R}$ such that $\nabla f=\mathbf {F}$. It should be noted that $f$ does not exist for most vector fields. The gradient theorem implies that for all closed continuous curves $C$, that $\int _{\mathbf {q} \in C}(\nabla f)|_{\mathbf {q} }\cdot d\mathbf {q} =0$. It must be the case that $\int _{\mathbf {q} \in C}\mathbf {F} (\mathbf {q} )\cdot d\mathbf {q} =0$ for all closed continuous curves $C$: $\mathbf {F}$ must be conservative. Equivalently, $\nabla \times \mathbf {F} =\mathbf {0}$ everywhere: $\mathbf {F}$ should be irrotational.
Given a conservative vector field $\mathbf {F}$, $f$ can be determined by choosing an origin point $\mathbf {q} _{0}$, and a constant $f_{0}$. For each $\mathbf {q} \in \mathbb {R} ^{3}$, a continuous oriented curve $C[\mathbf {q} ]$ that starts at $\mathbf {q} _{0}$ and ends at $\mathbf {q}$ should be generated. $f(\mathbf {q} )$ is assigned: $f(\mathbf {q} )=f_{0}+\int _{\mathbf {q} '\in C[\mathbf {q} ]}\mathbf {F} (\mathbf {q} ')\cdot d\mathbf {q} '$. The choice of curve $C[\mathbf {q} ]$ is irrelevant since $\mathbf {F}$ is conservative/irrotational.
It can be confirmed that $f(\mathbf {q} )=f_{0}+\int _{\mathbf {q} '\in C[\mathbf {q} ]}\mathbf {F} (\mathbf {q} ')\cdot d\mathbf {q} '$ satisfies $\nabla f=\mathbf {F}$ for each point $\mathbf {q} \in \mathbb {R} ^{3}$ by evaluating the directional derivative along arbitrary curves that pass through $\mathbf {q}$. For an arbitrary curve $C'$ parameterized by $t$ that passes through $\mathbf {q}$ at $t=t_{0}$, the directional derivative of $f$ at $t=t_{0}$ is $\mathbf {F} (\mathbf {q} )\cdot {\frac {d\mathbf {q} _{C'}}{dt}}{\bigg |}_{t=t_{0}}$. This confirms that $\nabla f=\mathbf {F}$.

### Spherical volume integral solution

Assume that $\mathbf {F} (\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{F}})$ as $|\mathbf {q} |\to +\infty$ for some $\alpha _{F}>1$. This means that for some $\alpha _{F}>1$, that there exists $c_{1},c_{2}>0$ such that $\forall \mathbf {q} \in \mathbb {R} ^{3}:|\mathbf {q} |>c_{1}\implies |\mathbf {F} (\mathbf {q} )|<c_{2}(1/|\mathbf {q} |^{\alpha _{F}})$.
This property implies that given a sufficiently large $R\gg 0$, that any path integral between any two points outside of the sphere $B(\mathbf {0} ;R)$ that remains outside of the sphere is arbitrarily small. Therefore the origin point $\mathbf {q} _{0}$ can freely shift between points at infinity.
Choose infinity as the origin point and let $f_{0}=0$. Given an arbitrary point $\mathbf {q}$, choose a direction quantified by unit length vector $\mathbf {u}$. The path $C[\mathbf {q} ]$ travels backwards along a ray that starts at $\mathbf {q}$ and points in the direction given by $\mathbf {u}$. $f(\mathbf {q} )=\int _{\mathbf {q} '\in C[\mathbf {q} ]}\mathbf {F} (\mathbf {q} ')\cdot d\mathbf {q} '=\int _{r'=0}^{+\infty }\mathbf {F} (\mathbf {q} +r'\mathbf {u} )\cdot (dr'(-\mathbf {u} ))$. Since this integral does not depend on the direction $\mathbf {u}$, the average over all directions $\mathbf {u}$ is:
$f(\mathbf {q} )={\frac {1}{4\pi }}\iint _{|\mathbf {u} |=1}\left(\int _{r'=0}^{+\infty }\mathbf {F} (\mathbf {q} +r'\mathbf {u} )\cdot (dr'(-\mathbf {u} ))\right)|\mathbf {dS} _{\mathbf {u} }|$ where $\mathbf {dS} _{\mathbf {u} }$ is an infinitesimal surface portion of the unit sphere that constrains $\mathbf {u}$.
$f(\mathbf {q} )={\frac {1}{4\pi }}\iint _{|\mathbf {u} |=1}\int _{r'=0}^{+\infty }{\frac {\mathbf {F} (\mathbf {q} +r'\mathbf {u} )\cdot (-\mathbf {u} )}{r'^{2}}}(dr'\cdot r'^{2}|\mathbf {dS} _{\mathbf {u} }|)$
Letting the position vector variable $\mathbf {q} ''$ be $\mathbf {q} ''=r'\mathbf {u}$, the new volume differential is $dV''=dr'\cdot r'^{2}|\mathbf {dS} _{\mathbf {u} }|$.
$f(\mathbf {q} )={\frac {1}{4\pi }}\iiint _{\mathbf {q} ''\in \mathbb {R} ^{3}}{\frac {\mathbf {F} (\mathbf {q} +\mathbf {q} '')\cdot (-\mathbf {q} ''/|\mathbf {q} ''|)}{|\mathbf {q} ''|^{2}}}dV''$
Lastly, letting $\mathbf {q} '=\mathbf {q} +\mathbf {q} ''$ gives:
$f(\mathbf {q} )={\frac {1}{4\pi }}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {\mathbf {F} (\mathbf {q} ')\cdot (\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}dV'$
This new formula $f(\mathbf {q} )={\frac {1}{4\pi }}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {\mathbf {F} (\mathbf {q} ')\cdot (\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}dV'$ is a volume integral that expresses the potential $f$ as a linear combination of functions that exhibit a degree of spherical symmetry. This formula is similar to the inverse square law for the inverse of the divergence.
Next, the above formula will be derived using a Green's function approach in a manner similar to the inverse square law for the divergence operator.

### Green's function solution

This section will derive a formula identical to the formula above using a Green's function approach. While the derivation will be complicated and the result will not be new, the derivation itself will yield many interesting intermediate results. Again, it will be assumed that $\mathbf {F} (\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{F}})$ as $|\mathbf {q} |\to +\infty$ for some $\alpha _{F}>1$. It will also be assumed that $\mathbf {F}$ is continuous.
Vector field $\mathbf {F} (\mathbf {q} )=F_{x}(\mathbf {q} )\mathbf {i} +F_{y}(\mathbf {q} )\mathbf {j} +F_{z}(\mathbf {q} )\mathbf {k}$ can be expressed as the following linear combination of Dirac delta functions: $\mathbf {F} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}(F_{x}(\mathbf {q} ')(\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {i} )+F_{y}(\mathbf {q} ')(\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {j} )+F_{z}(\mathbf {q} ')(\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {k} ))dV'$. This linear combination however does not facilitate the inversion of the gradient since for each $\mathbf {q} '\in \mathbb {R} ^{3}$, the vector fields $\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {i}$, $\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {j}$, and $\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {k}$ are all not conservative, so there do not exist any scalar fields $h_{x}(\mathbf {q} ;\mathbf {q} ')$, $h_{y}(\mathbf {q} ;\mathbf {q} ')$, and $h_{z}(\mathbf {q} ;\mathbf {q} ')$ such that $\nabla _{\mathbf {q} }h_{x}(\mathbf {q} ;\mathbf {q} ')=\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {i}$, $\nabla _{\mathbf {q} }h_{y}(\mathbf {q} ;\mathbf {q} ')=\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {j}$, and $\nabla _{\mathbf {q} }h_{z}(\mathbf {q} ;\mathbf {q} ')=\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {k}$. This prevents a simple solution of $f(\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}(F_{x}(\mathbf {q} ')h_{x}(\mathbf {q} ;\mathbf {q} ')+F_{y}(\mathbf {q} ')h_{y}(\mathbf {q} ;\mathbf {q} ')+F_{z}(\mathbf {q} ')h_{z}(\mathbf {q} ;\mathbf {q} '))dV'$.
It is first necessary to express $\mathbf {F}$ as a linear combination of vector fields that are conservative, so that for each component/basis vector field, a scalar field exists where the gradient is the component vector field. In mathematical terms $\mathbf {F}$ should be decomposed into the following linear combination:
$\mathbf {F} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}(F_{x}(\mathbf {q} ')\mathbf {G} _{x}(\mathbf {q} ;\mathbf {q} ')+F_{y}(\mathbf {q} ')\mathbf {G} _{y}(\mathbf {q} ;\mathbf {q} ')+F_{z}(\mathbf {q} ')\mathbf {G} _{z}(\mathbf {q} ;\mathbf {q} '))dV'$ where for each $\mathbf {q} '\in \mathbb {R} ^{3}$, vector fields $\mathbf {G} _{x}(\mathbf {q} ;\mathbf {q} ')$, $\mathbf {G} _{y}(\mathbf {q} ;\mathbf {q} ')$, and $\mathbf {G} _{z}(\mathbf {q} ;\mathbf {q} ')$ are all conservative. Specifically a vector valued function $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ')$ is required that is linear with respect to $\mathbf {F} '$: $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',F_{x}'\mathbf {i} +F_{y}'\mathbf {j} +F_{z}'\mathbf {k} )=F_{x}'\mathbf {G} _{x}(\mathbf {q} ;\mathbf {q} ')+F_{y}'\mathbf {G} _{y}(\mathbf {q} ;\mathbf {q} ')+F_{z}'\mathbf {G} _{z}(\mathbf {q} ;\mathbf {q} ')$. In essence, $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ')$ is effectively a linear combination of $\mathbf {G} _{x}(\mathbf {q} ;\mathbf {q} ')$, $\mathbf {G} _{y}(\mathbf {q} ;\mathbf {q} ')$, and $\mathbf {G} _{z}(\mathbf {q} ;\mathbf {q} ')$ where the coefficients are the components of $\mathbf {F} '$. The linear combination that $\mathbf {F}$ must be decomposed into is: $\mathbf {F} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} '))dV'$.
A candidate function that exhibits a degree of spherical symmetry has the form:
$\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ')=a(|\mathbf {q} -\mathbf {q} '|){\textbf {proj}}(\mathbf {F} ';\mathbf {q} -\mathbf {q} ')+b(|\mathbf {q} -\mathbf {q} '|){\textbf {perp}}(\mathbf {F} ';\mathbf {q} -\mathbf {q} ')$
where ${\textbf {proj}}(\mathbf {v} ;\mathbf {u} )={\frac {(\mathbf {v} \cdot \mathbf {u} )\mathbf {u} }{|\mathbf {u} |^{2}}}$ is the projection of $\mathbf {v}$ onto $\mathbf {u}$, and ${\textbf {perp}}(\mathbf {v} ;\mathbf {u} )=\mathbf {v} -{\textbf {proj}}(\mathbf {v} ;\mathbf {u} )$ is the perpendicular component of $\mathbf {v}$ relative to $\mathbf {u}$. $a(r)$ and $b(r)$ are real valued functions that scale the components of $\mathbf {F} '$ relative to the displacement $\mathbf {q} -\mathbf {q} '$. $a(r)$ and $b(r)$ have yet to be determined.
$a(r)$ and $b(r)$ have to be chosen such that $\mathbf {F} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} '))dV'$ provided that $\mathbf {F}$ satisfies all of the required conditions (most importantly $\mathbf {F}$ is conservative) and that $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ')$ is always a conservative vector field with $\mathbf {q}$ as the position parameter.
Using the assumption that $\mathbf {F}$ is continuous, for an arbitrarily small $\epsilon$, it can be assumed that $\mathbf {F} (\mathbf {q} )\approx \mathbf {F} (\mathbf {q} ')$ for $|\mathbf {q} -\mathbf {q} '|<\epsilon$. $\epsilon >0$ can be added as an additional parameter to $\mathbf {G}$, $a$, and $b$ to get:
$\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )=a(|\mathbf {q} -\mathbf {q} '|;\epsilon ){\textbf {proj}}(\mathbf {F} ';\mathbf {q} -\mathbf {q} ')+b(|\mathbf {q} -\mathbf {q} '|;\epsilon ){\textbf {perp}}(\mathbf {F} ';\mathbf {q} -\mathbf {q} ')$
The condition that $\mathbf {F} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} '))dV'$ becomes: $\mathbf {F} (\mathbf {q} )=\lim _{\epsilon \to 0^{+}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'$.
It can be proven with difficulty that if $k(r)$ where $r>0$ is an arbitrary monotone decreasing function and $\lim _{r\to 0^{+}}k(r)=+\infty$ and $\lim _{r\to +\infty }k(r)=0$, then choosing $a(r;\epsilon )={\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}c&(r\leq \epsilon )\\{\frac {(c-1)\epsilon ^{2}}{k(\epsilon )^{2}}}{\frac {k(r)^{2}}{r^{2}}}&(r>\epsilon )\end{array}}\right.$ and $b(r;\epsilon )={\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}c&(r\leq \epsilon )\\{\frac {(c-1)\epsilon ^{2}}{k(\epsilon )^{2}}}{\frac {k(r)}{r}}{\frac {dk}{dr}}{\bigg |}_{r}&(r>\epsilon )\end{array}}\right.$ will satisfy the condition $\mathbf {F} (\mathbf {q} )=\lim _{\epsilon \to 0^{+}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'$.
$\mathbf {F}$ being conservative, continuous, and vanishing at infinity are all important.
Theorem Statement
Let $\mathbf {F}$ be an arbitrary vector field that is:
    Conservative
    Continuous
    $O(1/|\mathbf {q} |^{\alpha _{F}})$ as $|\mathbf {q} |\to +\infty$ for some $\alpha _{F}>1$
Let arbitrary function $k:(0,+\infty )\to (0,+\infty )$ be such that:
    $k$ is monotone decreasing
    $\lim _{r\to 0^{+}}k(r)=+\infty$ and $\lim _{r\to +\infty }k(r)=0$
Let $c$ be an arbitrary real number.
If $a(r;\epsilon )={\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}c&(r\leq \epsilon )\\{\frac {(c-1)\epsilon ^{2}}{k(\epsilon )^{2}}}{\frac {k(r)^{2}}{r^{2}}}&(r>\epsilon )\end{array}}\right.$, $b(r;\epsilon )={\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}c&(r\leq \epsilon )\\{\frac {(c-1)\epsilon ^{2}}{k(\epsilon )^{2}}}{\frac {k(r)}{r}}{\frac {dk}{dr}}{\bigg |}_{r}&(r>\epsilon )\end{array}}\right.$, and $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )=a(|\mathbf {q} -\mathbf {q} '|;\epsilon ){\textbf {proj}}(\mathbf {F} ';\mathbf {q} -\mathbf {q} ')+b(|\mathbf {q} -\mathbf {q} '|;\epsilon ){\textbf {perp}}(\mathbf {F} ';\mathbf {q} -\mathbf {q} ')$ then
$\mathbf {F} (\mathbf {q} )=\lim _{\epsilon \to 0^{+}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'$
Proof:
To begin, fix $\mathbf {q}$. To simplify notation let:
    $r=|\mathbf {q} '-\mathbf {q} |$
    ${\hat {\mathbf {r} }}={\frac {\mathbf {q} '-\mathbf {q} }{r}}$ denotes the unit length vector at position $\mathbf {q} '$ that points away from $\mathbf {q}$
    $\theta \in [0,\pi ]$ denote the angle between $\mathbf {q} '-\mathbf {q}$ and $\mathbf {F} (\mathbf {q} )$
When $r\leq \epsilon$, $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )={\frac {c}{{\frac {4}{3}}\pi \epsilon ^{3}}}\mathbf {F} (\mathbf {q} ')$. $\mathbf {F}$ is continuous, so approximating $\mathbf {F} (\mathbf {q} ')$ with $\mathbf {F} (\mathbf {q} )$ gives: $\iiint _{r\leq \epsilon }\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'$ $=\iiint _{r\leq \epsilon }{\frac {c}{{\frac {4}{3}}\pi \epsilon ^{3}}}\mathbf {F} (\mathbf {q} ')dV'$ $\approx c\mathbf {F} (\mathbf {q} )$. The error relative to $\epsilon$ in the approximation vanishes as $\epsilon \to 0^{+}$.
The challenge is to evaluate the integral when $r>\epsilon$: $\iiint _{r>\epsilon }\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'$ where $\mathbf {F} (\mathbf {q} ')$ can no longer be approximated by $\mathbf {F} (\mathbf {q} )$. To start, hollow out the vector field around $\mathbf {q}$ to get: $\mathbf {F} _{\text{hollow}}(\mathbf {q} ')=I(\mathbf {q} ';\epsilon ){\vec {\delta }}(\mathbf {q} ';s(\mathbf {q} ;\epsilon ))+\left\{{\begin{array}{cc}\mathbf {0} &(r\leq \epsilon )\\\mathbf {F} (\mathbf {q} ')&(r>\epsilon )\end{array}}\right.$ where $s(\mathbf {q} ;\epsilon )=\{\mathbf {q} ':|\mathbf {q} '-\mathbf {q} |=\epsilon \}$ denotes the outwards oriented surface of a sphere centered on $\mathbf {q}$ with radius $\epsilon$. $I(\mathbf {q} ';\epsilon )$ is chosen so that the term $I(\mathbf {q} ';\epsilon ){\vec {\delta }}(\mathbf {q} ';s(\mathbf {q} ;\epsilon ))$ makes $\mathbf {F} _{\text{hollow}}$ conservative. For $\epsilon$ arbitrarily small, $I(\mathbf {q} ';\epsilon )\approx |\mathbf {F} (\mathbf {q} )|\epsilon \cos \theta$. Since $\mathbf {F} (\mathbf {q} ')$ is $O(1/|\mathbf {q} '|^{\alpha _{F}})$ as $|\mathbf {q} '|\to +\infty$ for some $\alpha _{F}>1$, $\mathbf {F} _{\text{hollow}}$ can be viewed as a multi closed surface.
Given a closed surface $\sigma$, consider moving the points according to the following transformation: $\mathbf {q} '\mapsto \mathbf {q} +k(r){\hat {\mathbf {r} }}$. Given infinitesimal surface vector $\mathbf {dS}$ that is parallel to ${\hat {\mathbf {r} }}$, then $\mathbf {dS}$ is scaled by a factor of ${\frac {k(r)^{2}}{r^{2}}}$. If $\mathbf {dS}$ is perpendicular to ${\hat {\mathbf {r} }}$, then $\mathbf {dS}$ is scaled by a factor of ${\frac {k(r)}{r}}{\frac {dk}{dr}}{\bigg |}_{r}$. An arbitrary surface vector $\mathbf {dS}$ is stretched and squeezed to become $\mathbf {g} '(\mathbf {q} ',\mathbf {dS} )={\frac {k(r)^{2}}{r^{2}}}{\textbf {proj}}(\mathbf {dS} ;{\hat {\mathbf {r} }})+{\frac {k(r)}{r}}{\frac {dk}{dr}}{\bigg |}_{r}{\textbf {perp}}(\mathbf {dS} ;{\hat {\mathbf {r} }})$. Replacing each $\mathbf {dS}$ with $\mathbf {g} '(\mathbf {q} ',\mathbf {dS} )$ does not rupture the surface, and the total surface vector remains 0. Hence, $\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {g} '(\mathbf {q} ',{\vec {\delta }}(\mathbf {q} ';\sigma ))dV'=\iint _{\mathbf {q} '\in \sigma }\mathbf {g} '(\mathbf {q} ',\mathbf {dS} ')=0$. Since $\mathbf {F} _{\text{hollow}}$ is a multi closed surface, $\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {g} '(\mathbf {q} ',\mathbf {F} _{\text{hollow}}(\mathbf {q} '))dV'=0$
We now evaluate:
$\iiint _{r>\epsilon }\mathbf {g} '(\mathbf {q} ',\mathbf {F} (\mathbf {q} '))dV'=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {g} '(\mathbf {q} ',\mathbf {F} _{\text{hollow}}(\mathbf {q} '))dV'-\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {g} '(\mathbf {q} ',I(\mathbf {q} ';\epsilon ){\vec {\delta }}(\mathbf {q} ';s(\mathbf {q} ;\epsilon )))dV'$
$=0-\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {k(r)^{2}}{r^{2}}}I(\mathbf {q} ';\epsilon ){\vec {\delta }}(\mathbf {q} ';s(\mathbf {q} ;\epsilon ))dV'$
$\approx -\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {k(r)^{2}}{r^{2}}}|\mathbf {F} (\mathbf {q} )|\epsilon \cos \theta {\vec {\delta }}(\mathbf {q} ';s(\mathbf {q} ;\epsilon ))dV'$
$=-|\mathbf {F} (\mathbf {q} )|{\frac {k(\epsilon )^{2}}{\epsilon }}\iint _{\mathbf {q} '\in s(\mathbf {q} ;\epsilon )}\cos \theta \mathbf {dS} '$
$=-|\mathbf {F} (\mathbf {q} )|{\frac {k(\epsilon )^{2}}{\epsilon }}({\frac {4}{3}}\pi \epsilon ^{2}{\frac {\mathbf {F} (\mathbf {q} )}{|\mathbf {F} (\mathbf {q} )|}})$
$=-{\frac {4}{3}}\pi \epsilon k(\epsilon )^{2}\mathbf {F} (\mathbf {q} )$
It is the case that:
$\iiint _{r>\epsilon }\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'={\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}{\frac {(c-1)\epsilon ^{2}}{k(\epsilon )^{2}}}\iiint _{r>\epsilon }\mathbf {g} '(\mathbf {q} ',\mathbf {F} (\mathbf {q} '))dV'$ $\approx {\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}({\frac {4}{3}}\pi \epsilon ^{3}(1-c)\mathbf {F} (\mathbf {q} ))$ $=(1-c)\mathbf {F} (\mathbf {q} )$
Therefore: $\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'=\iiint _{r\leq \epsilon }\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'+\iiint _{r>\epsilon }\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'\approx c\mathbf {F} (\mathbf {q} )+(1-c)\mathbf {F} (\mathbf {q} )=\mathbf {F} (\mathbf {q} )$ so $\lim _{\epsilon \to 0^{+}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'=\mathbf {F} (\mathbf {q} )$
End of proof
For $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$ to always be a conservative vector field with $\mathbf {q}$ as the position parameter, $a(r;\epsilon )=b(r;\epsilon )+r{\frac {\partial b}{\partial r}}(r;\epsilon )$ for all $r\geq 0$. Choosing $k(r)={\frac {1}{\sqrt {r}}}$ and $c={\frac {1}{3}}$ will make $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$ conservative. Hence $a(r;\epsilon )={\frac {1}{4\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}1&(r\leq \epsilon )\\-{\frac {2\epsilon ^{3}}{r^{3}}}&(r>\epsilon )\end{array}}\right.$ and $b(r;\epsilon )={\frac {1}{4\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}1&(r\leq \epsilon )\\{\frac {\epsilon ^{3}}{r^{3}}}&(r>\epsilon )\end{array}}\right.$
Lastly, a scalar valued function $h(\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$ should be chosen that satisfies $\nabla _{\mathbf {q} }h(\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )=\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$. One possible $h$ is:
$h(\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )={\frac {1}{4\pi }}(\mathbf {F} '\cdot (\mathbf {q} -\mathbf {q} '))\left\{{\begin{array}{cc}{\frac {1}{\epsilon ^{3}}}&(|\mathbf {q} -\mathbf {q} '|\leq \epsilon )\\{\frac {1}{|\mathbf {q} -\mathbf {q} '|^{3}}}&(|\mathbf {q} -\mathbf {q} '|>\epsilon )\end{array}}\right.$
Therefore one possible solution is $f(\mathbf {q} )=\lim _{\epsilon \to 0^{+}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}h(\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'={\frac {1}{4\pi }}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {\mathbf {F} (\mathbf {q} ')\cdot (\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}dV'$. This is identical to the solution given in the previous section.

## Inverting the curl operator (the Boit-Savart law)

Given a vector field $\mathbf {F} :\mathbb {R} ^{3}\to \mathbb {R} ^{3}$, the problem of interest is finding a divergence free vector field $\mathbf {B} :\mathbb {R} ^{3}\to \mathbb {R} ^{3}$ such that $\nabla \times \mathbf {B} =\mathbf {F}$ and $\nabla \cdot \mathbf {B} =0$. Since the curl is always divergence free, $\nabla \cdot (\nabla \times \mathbf {B} )=0$ always, it is required that $\mathbf {F}$ be divergence free everywhere: $\nabla \cdot \mathbf {F} =0$.

### Green's function solution

The Green's function solution that will be presented here will be very similar to the Green's function solution given for inverting the gradient. It will be assumed that $\mathbf {F} (\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{F}})$ as $|\mathbf {q} |\to +\infty$ for some $\alpha _{F}>1$. It will also be assumed that $\mathbf {F}$ is continuous.
Vector field $\mathbf {F} (\mathbf {q} )=F_{x}(\mathbf {q} )\mathbf {i} +F_{y}(\mathbf {q} )\mathbf {j} +F_{z}(\mathbf {q} )\mathbf {k}$ can be expressed as the following linear combination of Dirac delta functions: $\mathbf {F} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}(F_{x}(\mathbf {q} ')(\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {i} )+F_{y}(\mathbf {q} ')(\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {j} )+F_{z}(\mathbf {q} ')(\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {k} ))dV'$. This linear combination however does not facilitate the inversion of the curl since for each $\mathbf {q} '\in \mathbb {R} ^{3}$, the vector fields $\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {i}$, $\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {j}$, and $\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {k}$ are all not divergence free, so there do not exist any vector fields $\mathbf {H} _{x}(\mathbf {q} ;\mathbf {q} ')$, $\mathbf {H} _{y}(\mathbf {q} ;\mathbf {q} ')$, and $\mathbf {H} _{z}(\mathbf {q} ;\mathbf {q} ')$ such that $\nabla _{\mathbf {q} }\times \mathbf {H} _{x}(\mathbf {q} ;\mathbf {q} ')=\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {i}$, $\nabla _{\mathbf {q} }\times \mathbf {H} _{y}(\mathbf {q} ;\mathbf {q} ')=\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {j}$, and $\nabla _{\mathbf {q} }\times \mathbf {H} _{z}(\mathbf {q} ;\mathbf {q} ')=\delta (\mathbf {q} ;\mathbf {q} ')\mathbf {k}$. This prevents a simple solution of $\mathbf {B} (\mathbf {q} )=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}(F_{x}(\mathbf {q} ')\mathbf {H} _{x}(\mathbf {q} ;\mathbf {q} ')+F_{y}(\mathbf {q} ')\mathbf {H} _{y}(\mathbf {q} ;\mathbf {q} ')+F_{z}(\mathbf {q} ')\mathbf {H} _{z}(\mathbf {q} ;\mathbf {q} '))dV'$.
Using identical reasoning as was used for inverting the gradient, a vector valued function $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$ is required such that:
* $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$ is linear with respect to $\mathbf {F} '$ so that
$\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )=\mathbf {G} (\mathbf {q} ;\mathbf {q} ',F_{x}(\mathbf {q} ')\mathbf {i} +F_{y}(\mathbf {q} ')\mathbf {j} +F_{z}(\mathbf {q} ')\mathbf {k} ;\epsilon )=F_{x}(\mathbf {q} ')\mathbf {G} _{x}(\mathbf {q} ;\mathbf {q} ';\epsilon )+F_{y}(\mathbf {q} ')\mathbf {G} _{y}(\mathbf {q} ;\mathbf {q} ';\epsilon )+F_{z}(\mathbf {q} ')\mathbf {G} _{z}(\mathbf {q} ;\mathbf {q} ';\epsilon )$
* $\mathbf {F} (\mathbf {q} )=\lim _{\epsilon \to 0^{+}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'$.
* $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$ is divergence free with respect to $\mathbf {q}$ as the input parameter.
It can again be proven with difficulty that if $k(r)$ where $r>0$ is an arbitrary monotone decreasing function and $\lim _{r\to 0^{+}}k(r)=+\infty$, and $\lim _{r\to +\infty }k(r)=0$, then choosing $a(r;\epsilon )={\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}c&(r\leq \epsilon )\\{\frac {(c-1)\epsilon }{k(\epsilon )}}{\frac {dk}{dr}}{\bigg |}_{r}&(r>\epsilon )\end{array}}\right.$ and $b(r;\epsilon )={\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}c&(r\leq \epsilon )\\{\frac {(c-1)\epsilon }{k(\epsilon )}}{\frac {k(r)}{r}}&(r>\epsilon )\end{array}}\right.$ will satisfy the condition $\mathbf {F} (\mathbf {q} )=\lim _{\epsilon \to 0^{+}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'$.
$\mathbf {F}$ being divergence free, continuous, and vanishing at infinity are all important.
Theorem Statement
Let $\mathbf {F}$ be an arbitrary vector field that is:
    Divergence Free
    Continuous
    $O(1/|\mathbf {q} |^{\alpha _{F}})$ as $|\mathbf {q} |\to +\infty$ for some $\alpha _{F}>1$
Let arbitrary function $k:(0,+\infty )\to (0,+\infty )$ be such that:
    $k$ is monotone decreasing
    $\lim _{r\to 0^{+}}k(r)=+\infty$ and $\lim _{r\to +\infty }k(r)=0$
Let $c$ be an arbitrary real number.
If $a(r;\epsilon )={\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}c&(r\leq \epsilon )\\{\frac {(c-1)\epsilon }{k(\epsilon )}}{\frac {dk}{dr}}{\bigg |}_{r}&(r>\epsilon )\end{array}}\right.$, $b(r;\epsilon )={\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}c&(r\leq \epsilon )\\{\frac {(c-1)\epsilon }{k(\epsilon )}}{\frac {k(r)}{r}}&(r>\epsilon )\end{array}}\right.$, and $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )=a(|\mathbf {q} -\mathbf {q} '|;\epsilon ){\textbf {proj}}(\mathbf {F} ';\mathbf {q} -\mathbf {q} ')+b(|\mathbf {q} -\mathbf {q} '|;\epsilon ){\textbf {perp}}(\mathbf {F} ';\mathbf {q} -\mathbf {q} ')$ then
$\mathbf {F} (\mathbf {q} )=\lim _{\epsilon \to 0^{+}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'$
Proof:
To begin, fix $\mathbf {q}$. To simplify notation let:
    $r=|\mathbf {q} '-\mathbf {q} |$
    ${\hat {\mathbf {r} }}={\frac {\mathbf {q} '-\mathbf {q} }{r}}$ denotes the unit length vector at position $\mathbf {q} '$ that points away from $\mathbf {q}$
    $\theta \in [0,\pi ]$ denote the angle between $\mathbf {q} '-\mathbf {q}$ and $\mathbf {F} (\mathbf {q} )$
    ${\hat {\mathbf {\phi } }}={\frac {\mathbf {F} (\mathbf {q} )\times (\mathbf {q} '-\mathbf {q} )}{|\mathbf {F} (\mathbf {q} )\times (\mathbf {q} '-\mathbf {q} )|}}$ denote the unit length vector that points in a counterclockwise direction around $\mathbf {F} (\mathbf {q} )$
When $r\leq \epsilon$, $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )={\frac {c}{{\frac {4}{3}}\pi \epsilon ^{3}}}\mathbf {F} (\mathbf {q} ')$. $\mathbf {F}$ is continuous, so approximating $\mathbf {F} (\mathbf {q} ')$ with $\mathbf {F} (\mathbf {q} )$ gives: $\iiint _{r\leq \epsilon }\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'$ $=\iiint _{r\leq \epsilon }{\frac {c}{{\frac {4}{3}}\pi \epsilon ^{3}}}\mathbf {F} (\mathbf {q} ')dV'$ $\approx c\mathbf {F} (\mathbf {q} )$. The error relative to $\epsilon$ in the approximation vanishes as $\epsilon \to 0^{+}$.
The challenge is to evaluate the integral when $r>\epsilon$: $\iiint _{r>\epsilon }\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'$ where $\mathbf {F} (\mathbf {q} ')$ can no longer be approximated by $\mathbf {F} (\mathbf {q} )$. To start, hollow out the vector field around $\mathbf {q}$ to get: $\mathbf {F} _{\text{hollow}}(\mathbf {q} ')=I(\mathbf {q} ';\epsilon )({\vec {\delta }}(\mathbf {q} ';s(\mathbf {q} ;\epsilon ))\times {\hat {\mathbf {\phi } }})+\left\{{\begin{array}{cc}\mathbf {0} &(r\leq \epsilon )\\\mathbf {F} (\mathbf {q} ')&(r>\epsilon )\end{array}}\right.$ where $s(\mathbf {q} ;\epsilon )=\{\mathbf {q} ':|\mathbf {q} '-\mathbf {q} |=\epsilon \}$ denotes the outwards oriented surface of a sphere centered on $\mathbf {q}$ with radius $\epsilon$. $I(\mathbf {q} ';\epsilon )$ is chosen so that the term $I(\mathbf {q} ';\epsilon )({\vec {\delta }}(\mathbf {q} ';s(\mathbf {q} ;\epsilon ))\times {\hat {\mathbf {\phi } }})$ makes $\mathbf {F} _{\text{hollow}}$ divergence free. For $\epsilon$ arbitrarily small, $I(\mathbf {q} ';\epsilon )\approx {\frac {1}{2}}|\mathbf {F} (\mathbf {q} )|\epsilon \sin \theta$. Vector field $\mathbf {F} _{\text{hollow}}$ can be viewed as a superposition of paths that are either closed or extend to infinity. The fact that some paths may extend to infinity will become moot.
Given a path $C$ with no endpoints (either closed or extends to infinity) and that does not pass through $\mathbf {q}$, consider moving the points according to the following transformation: $\mathbf {q} '\mapsto \mathbf {q} +k(r){\hat {\mathbf {r} }}$. Given infinitesimal displacement vector $d\mathbf {q}$ that is parallel to ${\hat {\mathbf {r} }}$, then $d\mathbf {q}$ is scaled by a factor of ${\frac {dk}{dr}}$. If $d\mathbf {q}$ is perpendicular to ${\hat {\mathbf {r} }}$, then $d\mathbf {q}$ is scaled by a factor of ${\frac {k(r)}{r}}$. An arbitrary displacement vector $d\mathbf {q}$ is stretched and squeezed to become $\mathbf {g} '(\mathbf {q} ',d\mathbf {q} )={\frac {dk}{dr}}{\bigg |}_{r}{\textbf {proj}}(d\mathbf {q} ;{\hat {\mathbf {r} }})+{\frac {k(r)}{r}}{\textbf {perp}}(d\mathbf {q} ;{\hat {\mathbf {r} }})$. Replacing each $d\mathbf {q}$ with $\mathbf {g} '(\mathbf {q} ',d\mathbf {q} )$ does not rupture the curve. Moreover the inverting nature of $k$ closes curves that extend to infinity. Hence, $\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {g} '(\mathbf {q} ',{\vec {\delta }}(\mathbf {q} ';C))dV'=\int _{\mathbf {q} '\in C}\mathbf {g} '(\mathbf {q} ',d\mathbf {q} ')=0$. Since $\mathbf {F} _{\text{hollow}}$ is a superposition of paths with no endpoints and that do not pass through $\mathbf {q}$, $\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {g} '(\mathbf {q} ',\mathbf {F} _{\text{hollow}}(\mathbf {q} '))dV'=0$
We now evaluate:
$\iiint _{r>\epsilon }\mathbf {g} '(\mathbf {q} ',\mathbf {F} (\mathbf {q} '))dV'$ $=\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {g} '(\mathbf {q} ',\mathbf {F} _{\text{hollow}}(\mathbf {q} '))dV'-\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {g} '(\mathbf {q} ',I(\mathbf {q} ';\epsilon )({\vec {\delta }}(\mathbf {q} ';s(\mathbf {q} ;\epsilon ))\times {\hat {\mathbf {\phi } }}))dV'$
$=0-\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {k(r)}{r}}I(\mathbf {q} ';\epsilon )({\vec {\delta }}(\mathbf {q} ';s(\mathbf {q} ;\epsilon ))\times {\hat {\mathbf {\phi } }})dV'$
$\approx -\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {k(r)}{r}}{\frac {1}{2}}|\mathbf {F} (\mathbf {q} )|\epsilon \sin \theta ({\vec {\delta }}(\mathbf {q} ';s(\mathbf {q} ;\epsilon ))\times {\hat {\mathbf {\phi } }})dV'$
$=-{\frac {1}{2}}|\mathbf {F} (\mathbf {q} )|k(\epsilon )\iint _{\mathbf {q} '\in s(\mathbf {q} ;\epsilon )}\sin \theta (\mathbf {dS} '\times {\hat {\mathbf {\phi } }})$
$=-{\frac {1}{2}}|\mathbf {F} (\mathbf {q} )|k(\epsilon )({\frac {8}{3}}\pi \epsilon ^{2}{\frac {\mathbf {F} (\mathbf {q} )}{|\mathbf {F} (\mathbf {q} )|}})$
$=-{\frac {4}{3}}\pi k(\epsilon )\epsilon ^{2}\mathbf {F} (\mathbf {q} )$
It is the case that:
$\iiint _{r>\epsilon }\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'={\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}{\frac {\epsilon (c-1)}{k(\epsilon )}}\iiint _{r>\epsilon }\mathbf {g} '(\mathbf {q} ',\mathbf {F} (\mathbf {q} '))dV'$ $\approx {\frac {1}{{\frac {4}{3}}\pi \epsilon ^{3}}}({\frac {4}{3}}\pi \epsilon ^{3}(1-c)\mathbf {F} (\mathbf {q} ))$ $=(1-c)\mathbf {F} (\mathbf {q} )$
Therefore: $\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'=\iiint _{r\leq \epsilon }\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'+\iiint _{r>\epsilon }\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'\approx c\mathbf {F} (\mathbf {q} )+(1-c)\mathbf {F} (\mathbf {q} )=\mathbf {F} (\mathbf {q} )$ so $\lim _{\epsilon \to 0^{+}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'=\mathbf {F} (\mathbf {q} )$
End of proof
For $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$ to always be a divergence free vector field with $\mathbf {q}$ as the position parameter, $2(a(r;\epsilon )-b(r;\epsilon ))+r{\frac {\partial a}{\partial r}}(r;\epsilon )$ for all $r\geq 0$. Choosing $k(r)={\frac {1}{r^{2}}}$ and $c={\frac {2}{3}}$ will make $\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$ divergence free. Hence $a(r;\epsilon )={\frac {1}{2\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}1&(r\leq \epsilon )\\{\frac {\epsilon ^{3}}{r^{3}}}&(r>\epsilon )\end{array}}\right.$ and $b(r;\epsilon )={\frac {1}{2\pi \epsilon ^{3}}}\left\{{\begin{array}{cc}1&(r\leq \epsilon )\\-{\frac {\epsilon ^{3}}{2r^{3}}}&(r>\epsilon )\end{array}}\right.$
Lastly, a vector valued function $\mathbf {H} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$ should be chosen that satisfies $\nabla _{\mathbf {q} }\times \mathbf {H} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )=\mathbf {G} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )$ and $\nabla _{\mathbf {q} }\cdot \mathbf {H} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )=0$. One possible $\mathbf {H}$ is:
$\mathbf {H} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} ';\epsilon )={\frac {1}{4\pi }}(\mathbf {F} '\times (\mathbf {q} -\mathbf {q} '))\left\{{\begin{array}{cc}{\frac {1}{\epsilon ^{3}}}&(|\mathbf {q} -\mathbf {q} '|\leq \epsilon )\\{\frac {1}{|\mathbf {q} -\mathbf {q} '|^{3}}}&(|\mathbf {q} -\mathbf {q} '|>\epsilon )\end{array}}\right.$
Therefore one possible solution is $\mathbf {B} (\mathbf {q} )=\lim _{\epsilon \to 0^{+}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}\mathbf {H} (\mathbf {q} ;\mathbf {q} ',\mathbf {F} (\mathbf {q} ');\epsilon )dV'={\frac {1}{4\pi }}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {\mathbf {F} (\mathbf {q} ')\times (\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}dV'$.

## Summary

Let $f$ be a continuous scalar field such that $f(\mathbf {q} )$ is $O(1/|\mathbf {q} |^{\alpha _{f}})$ as $|\mathbf {q} |\to +\infty$ for some $\alpha _{f}>1$. An irrotational vector field $\mathbf {F}$ that satisfies $\nabla \cdot \mathbf {F} =f$ and $\nabla \times \mathbf {F} =\mathbf {0}$ is
    $\mathbf {F} (\mathbf {q} )={\frac {1}{4\pi }}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {f(\mathbf {q} ')(\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}dV'$
Let $\mathbf {F}$ be a continuous conservative ($\nabla \times \mathbf {F} =\mathbf {0}$) vector field that is $O(1/|\mathbf {q} |^{\alpha _{F}})$ as $|\mathbf {q} |\to +\infty$ for some $\alpha _{F}>1$. A scalar field $f$ that satisfies $\nabla f=\mathbf {F}$ is
    $f(\mathbf {q} )={\frac {1}{4\pi }}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {\mathbf {F} (\mathbf {q} ')\cdot (\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}dV'$
Let $\mathbf {F}$ be a continuous divergence free ($\nabla \cdot \mathbf {F} =0$) vector field that is $O(1/|\mathbf {q} |^{\alpha _{F}})$ as $|\mathbf {q} |\to +\infty$ for some $\alpha _{F}>1$. A vector field $\mathbf {B}$ that satisfies $\nabla \times \mathbf {B} =\mathbf {F}$ and $\nabla \cdot \mathbf {B} =0$ is
    $\mathbf {B} (\mathbf {q} )={\frac {1}{4\pi }}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {\mathbf {F} (\mathbf {q} ')\times (\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}dV'$

## Example Applications
### Electromagnetism

Maxwell's Equations for electromagnetism are:
* Gauss's Law of Electric Fields: $\nabla \cdot \mathbf {E} ={\frac {\rho }{\epsilon _{0}}}$
* Faraday's Law: $\nabla \times \mathbf {E} =-{\frac {\partial \mathbf {B} }{\partial t}}$
* Gauss's Law of Magnetic Fields: $\nabla \cdot \mathbf {B} =0$
* Ampere's Law with Maxwell's correction: $\nabla \times \mathbf {B} =\mu _{0}(\mathbf {J} +\epsilon _{0}{\frac {\partial \mathbf {E} }{\partial t}})$
$\mathbf {E}$ and $\mathbf {B}$ are vector fields that denote the electric and magnetic field respectively. $\rho$ is a scalar field that denotes the charge density. $\mathbf {J}$ is a vector field that denotes the current density. $\epsilon _{0}=8.85\times 10^{-12}{\text{C}}^{2}/({\text{Nm}}^{2})$ is a constant that denotes the electric permittivity of free space, and $\mu _{0}=4\pi \times 10^{-7}{\text{N}}/{\text{A}}^{2}$ is a constant that denotes the magnetic permeability of free space.
When there is no changes in $\mathbf {E}$, $\mathbf {B}$, $\rho$, and $\mathbf {J}$ with respect to time, Maxwell's equations are reduced to:
* $\nabla \cdot \mathbf {E} ={\frac {\rho }{\epsilon _{0}}}$
* $\nabla \times \mathbf {E} =\mathbf {0}$
* $\nabla \cdot \mathbf {B} =0$
* $\nabla \times \mathbf {B} =\mu _{0}\mathbf {J}$

#### Coulomb's Law

Solving for the electric field from $\nabla \cdot \mathbf {E} ={\frac {\rho }{\epsilon _{0}}}$ and $\nabla \times \mathbf {E} =\mathbf {0}$ gives: $\mathbf {E} (\mathbf {q} )={\frac {1}{4\pi \epsilon _{0}}}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {\rho (\mathbf {q} ')(\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}dV'$. This formula is referred to as "Coulomb's Law".

#### Biot-Savart Law

Solving for the magnetic field from $\nabla \times \mathbf {B} =\mu _{0}\mathbf {J}$ and $\nabla \cdot \mathbf {B} =0$ gives: $\mathbf {B} (\mathbf {q} )={\frac {\mu _{0}}{4\pi }}\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {\mathbf {J} (\mathbf {q} ')\times (\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}dV'$. This formula is referred to as the "Biot-Savart Law".

### Gravitational Fields

According to Newton's inverse square law of gravity, the attractive force between masses $m_{1}$ and $m_{2}$ is $F={\frac {Gm_{1}m_{2}}{r^{2}}}$ where $r$ is the distance between $m_{1}$ and $m_{2}$, and $G=6.67408\times 10^{-11}{\text{Nm}}^{2}/{\text{kg}}^{2}$ is the gravitational constant. The gravitational field (force per unit mass or acceleration) $\mathbf {g}$ generated by a point mass of $M$ located at $\mathbf {q} '$ is $\mathbf {g} (\mathbf {q} )={\frac {GM}{|\mathbf {q} -\mathbf {q} '|^{2}}}(-{\frac {\mathbf {q} -\mathbf {q} '}{|\mathbf {q} -\mathbf {q} '|}})=-{\frac {GM(\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}$
Given a scalar field $\rho$ that denotes the density of mass at each point $\mathbf {q}$, then the total gravitational field is: $\mathbf {g} (\mathbf {q} )=-G\iiint _{\mathbf {q} '\in \mathbb {R} ^{3}}{\frac {\rho (\mathbf {q} ')(\mathbf {q} -\mathbf {q} ')}{|\mathbf {q} -\mathbf {q} '|^{3}}}dV'$. From this formula, it can be seen that the gravitational field satisfies the following equations: $\nabla \cdot \mathbf {g} =-4\pi G\rho$ and $\nabla \times \mathbf {g} =\mathbf {0}$. This form of characterizing the gravitational field can yield interesting results.
Gauss's divergence theorem implies that given a closed surface $\sigma$, that the total inwards flux $\Phi _{g}$ of a gravitational field is related to the total mass $M$ contained by $\sigma$ by the equation $\Phi _{g}=4\pi GM$.
Consider a spherically symmetric mass distribution centered on the origin: $\rho (\mathbf {q} )=\rho _{r}(|\mathbf {q} |)$. It is clear that the gravitational field is also spherically symmetric: $\mathbf {g} (\mathbf {q} )=-g_{r}(|\mathbf {q} |){\frac {\mathbf {q} }{|\mathbf {q} |}}$. Given a sphere of radius $R$ centered on the origin, the total inwards gravitational flux is $\Phi _{g}=4\pi R^{2}g_{r}(R)$, and the total mass contained is $M=\int _{r=0}^{R}4\pi r^{2}\rho _{r}(r)dr=4\pi \int _{r=0}^{R}r^{2}\rho _{r}(r)dr$. Therefore:
$\Phi _{g}=4\pi GM\implies 4\pi R^{2}g_{r}(R)=(4\pi )^{2}G\int _{r=0}^{R}r^{2}\rho _{r}(r)dr\implies g_{r}(R)={\frac {4\pi G}{R^{2}}}\int _{r=0}^{R}r^{2}\rho _{r}(r)dr$
The gravitational field strength at a distance of $R$ from the origin is $g_{r}(R)={\frac {4\pi G}{R^{2}}}\int _{r=0}^{R}r^{2}\rho _{r}(r)dr$.
Given a solid spherical mass with a uniform density of $\rho _{0}$ and a radius of $R$, the radial density function is $\rho _{r}(r)=\left\{{\begin{array}{cc}\rho _{0}&(r\leq R)\\0&(r>R)\end{array}}\right.$ so the gravitational field strength at a radius of $r$ is:
$g_{r}(r)={\frac {4\pi G\rho _{0}}{3r^{2}}}\left\{{\begin{array}{cc}r^{3}&(r\leq R)\\R^{3}&(r>R)\end{array}}\right.$ $={\frac {4\pi G\rho _{0}}{3}}\left\{{\begin{array}{cc}r&(r\leq R)\\R^{3}/r^{2}&(r>R)\end{array}}\right.$
The gravitational field reaches its maximum of ${\frac {4\pi G\rho _{0}R}{3}}={\frac {G({\frac {4}{3}}\pi R^{3}\rho _{0})}{R^{2}}}={\frac {GM}{R^{2}}}$ at the sphere's surface.
Given a hollow spherical mass of $M$ and a radius of $R$, the radial density function is $\rho _{r}(r)={\frac {M}{4\pi R^{2}}}\delta (r-R)$, where $\delta$ is the 1-dimensional Dirac delta function. The gravitational field strength at a radius of $r$ is:
$g_{r}(r)={\frac {GM}{r^{2}}}\left\{{\begin{array}{cc}0&(r\leq R)\\1&(r>R)\end{array}}\right.$
Contrary to intuition (and many examples from science fiction), the gravitational field inside a hollow spherical shell or Dyson sphere does not pull towards the inner surface, but is instead 0, as shown in the image to the right. If there is any mass nested inside the spherical shell, the gravitational field of the nested mass takes over completely inside the shell.
* 1 2 3 4 5 Griffiths, D. J., Introduction to Electrodynamics, 3rd Edition, Prentice Hall, 1999.