For a domain $Omega subset mathbb{R}^n$ and a small number $frak{T} > 0$, let [ mathcal{E}_0(Omega) = lambda_1(Omega) + {frak{T}} {text{tor}}(Omega) = inf_{u, w in H^1_0(Omega)setminus {0}} frac{int | abla u|^2}{int u^2} + {frak{T}} int frac{1}{2} | abla w|^2 - w ] be a modification of the first Dirichlet eigenvalue of $Omega$. It is well-known that over all $Omega$ with a given volume, the only sets attaining the infimum of $mathcal{E}_0$ are balls $B_R$; this is the Faber-Krahn inequality. The main result of this paper is that, if for all $Omega$ with the same volume and barycenter as $B_R$ and whose boundaries are parametrized as small $C^2$ normal graphs over $partial B_R$ with bounded $C^2$ norm, [ int |u_{Omega} - u_{B_R}|^2 + |Omega triangle B_R|^2 leq C [mathcal{E}_0(Omega) - mathcal{E}_0(B_R)] ] (i.e. the Faber-Krahn inequality is linearly stable), then the same is true for any $Omega$ with the same volume and barycenter as $B_R$ without any smoothness assumptions (i.e. it is nonlinearly stable). Here $u_{Omega}$ stands for an $L^2$-normalized first Dirichlet eigenfunction of $Omega$. Related results are shown for Riemannian manifolds. The proof is based on a detailed analysis of some critical perturbations of Bernoulli-type free boundary problems. The topic of when linear stability is valid, as well as some applications, are considered in a companion paper.