The integral form of Gauss's law is:

${\displaystyle {\scriptstyle _{S}}}$ ${\displaystyle \mathbf {E} \cdot \mathrm {d} \mathbf {A} }$${\displaystyle ={\frac {Q}{\varepsilon _{0}}}}$

for any closed surface S containing charge Q. By the divergence theorem, this equation is equivalent to:

$${\displaystyle \iiint _{V}\nabla \cdot \mathbf {E} \,\mathrm {d} V={\frac {Q}{\varepsilon _{0}}}}$$

for any volume V containing charge Q. By the relation between charge and charge density, this equation is equivalent to: $${\displaystyle \iiint _{V}\nabla \cdot \mathbf {E} \,\mathrm {d} V=\iiint _{V}{\frac {\rho }{\varepsilon _{0}}}\,\mathrm {d} V}$$ for any volume V. In order for this equation to be simultaneously true for every possible volume V, it is necessary (and sufficient) for the integrands to be equal everywhere. Therefore, this equation is equivalent to:

$${\displaystyle \nabla \cdot \mathbf {E} ={\frac {\rho }{\varepsilon _{0}}}.}$$ Thus the integral and differential forms are equivalent.

In this proof, we will show that the equation $${\displaystyle \nabla \cdot \mathbf {E} ={\dfrac {\rho }{\varepsilon _{0}}}}$$ is equivalent to the equation $${\displaystyle \nabla \cdot \mathbf {D} =\rho _{\mathrm {free} }}$$ Note that we are only dealing with the differential forms, not the integral forms, but that is sufficient since the differential and integral forms are equivalent in each case, by the divergence theorem.

We introduce the polarization density P, which has the following relation to E and D: $${\displaystyle \mathbf {D} =\varepsilon _{0}\mathbf {E} +\mathbf {P} }$$ and the following relation to the bound charge: $${\displaystyle \rho _{\mathrm {bound} }=-\nabla \cdot \mathbf {P} }$$ Now, consider the three equations: $${\displaystyle {\begin{aligned}\rho _{\mathrm {bound} }&=\nabla \cdot (-\mathbf {P} )\\\rho _{\mathrm {free} }&=\nabla \cdot \mathbf {D} \\\rho &=\nabla \cdot (\varepsilon _{0}\mathbf {E} )\end{aligned}}}$$ The key insight is that the sum of the first two equations is the third equation. This completes the proof: The first equation is true by definition, and therefore the second equation is true if and only if the third equation is true. So the second and third equations are equivalent, which is what we wanted to prove.

Coulomb's law states that the electric field due to a stationary point charge is: $${\displaystyle \mathbf {E} (\mathbf {r} )={\frac {q}{4\pi \varepsilon _{0}}}{\frac {\mathbf {e} _{r}}{r^{2}}}}$$ where

• e_r is the radial unit vector,
• r is the radius, |r|,
• ε₀ is the electric constant,
• q is the charge of the particle, which is assumed to be located at the origin.

Using the expression from Coulomb's law, we get the total field at r by using an integral to sum the field at r due to the infinitesimal charge at each other point s in space, to give $${\displaystyle \mathbf {E} (\mathbf {r} )={\frac {1}{4\pi \varepsilon _{0}}}\int {\frac {\rho (\mathbf {s} )(\mathbf {r} -\mathbf {s} )}{|\mathbf {r} -\mathbf {s} |^{3}}}\,\mathrm {d} ^{3}\mathbf {s} }$$ where ρ is the charge density. If we take the divergence of both sides of this equation with respect to r, and use the known theorem^[9]

$${\displaystyle \nabla \cdot \left({\frac {\mathbf {r} }{|\mathbf {r} |^{3}}}\right)=4\pi \delta (\mathbf {r} )}$$ where δ(r) is the Dirac delta function, the result is $${\displaystyle \nabla \cdot \mathbf {E} (\mathbf {r} )={\frac {1}{\varepsilon _{0}}}\int \rho (\mathbf {s} )\,\delta (\mathbf {r} -\mathbf {s} )\,\mathrm {d} ^{3}\mathbf {s} }$$

Using the "sifting property" of the Dirac delta function, we arrive at $${\displaystyle \nabla \cdot \mathbf {E} (\mathbf {r} )={\frac {\rho (\mathbf {r} )}{\varepsilon _{0}}},}$$ which is the differential form of Gauss's law, as desired.

Let ${\displaystyle \Omega \subseteq R^{3}}$ be a bounded open set, and $${\displaystyle \mathbf {E} _{0}(\mathbf {r} )={\frac {1}{4\pi \varepsilon _{0}}}\int _{\Omega }\rho (\mathbf {r} '){\frac {\mathbf {r} -\mathbf {r} '}{\left\|\mathbf {r} -\mathbf {r} '\right\|^{3}}}\mathrm {d} \mathbf {r} '\equiv {\frac {1}{4\pi \varepsilon _{0}}}\int _{\Omega }e(\mathbf {r,\mathbf {r} '} ){\mathrm {d} \mathbf {r} '}}$$ be the electric field, with ${\displaystyle \rho (\mathbf {r} ')}$ a continuous function (density of charge).

It is true for all ${\displaystyle \mathbf {r} \neq \mathbf {r'} }$ that ${\displaystyle \nabla _{\mathbf {r} }\cdot \mathbf {e} (\mathbf {r,r'} )=0}$.

Consider now a compact set ${\displaystyle V\subseteq R^{3}}$ having a piecewise smooth boundary ${\displaystyle \partial V}$ such that ${\displaystyle \Omega \cap V=\emptyset }$. It follows that ${\displaystyle e(\mathbf {r,\mathbf {r} '} )\in C^{1}(V\times \Omega )}$ and so, for the divergence theorem:

$${\displaystyle \oint _{\partial V}\mathbf {E} _{0}\cdot d\mathbf {S} =\int _{V}\mathbf {\nabla } \cdot \mathbf {E} _{0}\,dV}$$

But because ${\displaystyle e(\mathbf {r,\mathbf {r} '} )\in C^{1}(V\times \Omega )}$,

$${\displaystyle \mathbf {\nabla } \cdot \mathbf {E} _{0}(\mathbf {r} )={\frac {1}{4\pi \varepsilon _{0}}}\int _{\Omega }\nabla _{\mathbf {r} }\cdot e(\mathbf {r,\mathbf {r} '} ){\mathrm {d} \mathbf {r} '}=0}$$ for the argument above (${\displaystyle \Omega \cap V=\emptyset \implies \forall \mathbf {r} \in V\ \ \forall \mathbf {r'} \in \Omega \ \ \ \mathbf {r} \neq \mathbf {r'} }$ and then ${\displaystyle \nabla _{\mathbf {r} }\cdot \mathbf {e} (\mathbf {r,r'} )=0}$)

Therefore the flux through a closed surface generated by some charge density outside (the surface) is null.

Now consider ${\displaystyle \mathbf {r} _{0}\in \Omega }$, and ${\displaystyle B_{R}(\mathbf {r} _{0})\subseteq \Omega }$ as the sphere centered in ${\displaystyle \mathbf {r} _{0}}$ having ${\displaystyle R}$ as radius (it exists because ${\displaystyle \Omega }$ is an open set).

Let ${\displaystyle \mathbf {E} _{B_{R}}}$ and ${\displaystyle \mathbf {E} _{C}}$ be the electric field created inside and outside the sphere respectively. Then,

${\displaystyle \mathbf {E} _{B_{R}}={\frac {1}{4\pi \varepsilon _{0}}}\int _{B_{R}(\mathbf {r} _{0})}e(\mathbf {r,\mathbf {r} '} ){\mathrm {d} \mathbf {r} '}}$, ${\displaystyle \mathbf {E} _{C}={\frac {1}{4\pi \varepsilon _{0}}}\int _{\Omega \setminus B_{R}(\mathbf {r} _{0})}e(\mathbf {r,\mathbf {r} '} ){\mathrm {d} \mathbf {r} '}}$ and ${\displaystyle \mathbf {E} _{B_{R}}+\mathbf {E} _{C}=\mathbf {E} _{0}}$

$${\displaystyle \Phi (R)=\oint _{\partial B_{R}(\mathbf {r} _{0})}\mathbf {E} _{0}\cdot d\mathbf {S} =\oint _{\partial B_{R}(\mathbf {r} _{0})}\mathbf {E} _{B_{R}}\cdot d\mathbf {S} +\oint _{\partial B_{R}(\mathbf {r} _{0})}\mathbf {E} _{C}\cdot d\mathbf {S} =\oint _{\partial B_{R}(\mathbf {r} _{0})}\mathbf {E} _{B_{R}}\cdot d\mathbf {S} }$$

The last equality follows by observing that ${\displaystyle (\Omega \setminus B_{R}(\mathbf {r} _{0}))\cap B_{R}(\mathbf {r} _{0})=\emptyset }$, and the argument above.

The RHS is the electric flux generated by a charged sphere, and so:

$${\displaystyle \Phi (R)={\frac {Q(R)}{\varepsilon _{0}}}={\frac {1}{\varepsilon _{0}}}\int _{B_{R}(\mathbf {r} _{0})}\rho (\mathbf {r} '){\mathrm {d} \mathbf {r} '}={\frac {1}{\varepsilon _{0}}}\rho (\mathbf {r} '_{c})|B_{R}(\mathbf {r} _{0})|}$$ with ${\displaystyle r'_{c}\in \ B_{R}(\mathbf {r} _{0})}$

Where the last equality follows by the mean value theorem for integrals. Using the squeeze theorem and the continuity of ${\displaystyle \rho }$, one arrives at:

$${\displaystyle \mathbf {\nabla } \cdot \mathbf {E} _{0}(\mathbf {r} _{0})=\lim _{R\to 0}{\frac {1}{|B_{R}(\mathbf {r} _{0})|}}\Phi (R)={\frac {1}{\varepsilon _{0}}}\rho (\mathbf {r} _{0})}$$

Taking S in the integral form of Gauss's law to be a spherical surface of radius r, centered at the point charge Q, we have

$${\displaystyle \oint _{S}\mathbf {E} \cdot d\mathbf {A} ={\frac {Q}{\varepsilon _{0}}}}$$

By the assumption of spherical symmetry, the integrand is a constant which can be taken out of the integral. The result is $${\displaystyle 4\pi r^{2}{\hat {\mathbf {r} }}\cdot \mathbf {E} (\mathbf {r} )={\frac {Q}{\varepsilon _{0}}}}$$ where r̂ is a unit vector pointing radially away from the charge. Again by spherical symmetry, E points in the radial direction, and so we get $${\displaystyle \mathbf {E} (\mathbf {r} )={\frac {Q}{4\pi \varepsilon _{0}}}{\frac {\hat {\mathbf {r} }}{r^{2}}}}$$ which is essentially equivalent to Coulomb's law. Thus the inverse-square law dependence of the electric field in Coulomb's law follows from Gauss's law.