  A discrete non-linear $\sigma$-model is obtained by triangulate both the
space-time $M^{d+1}$ and the target space $K$. If the path integral is given by
the sum of all the complex homomorphisms $\phi: M^{d+1} \to K$, with an
partition function that is independent of space-time triangulation, then the
corresponding non-linear $\sigma$-model will be called topological non-linear
$\sigma$-model which is exactly soluble. Those exactly soluble models suggest
that phase transitions induced by fluctuations with no topological defects
(i.e. fluctuations described by homomorphisms $\phi$) usually produce a
topologically ordered state and are topological phase transitions, while phase
transitions induced by fluctuations with all the topological defects give rise
to trivial product states and are not topological phase transitions. If $K$ is
a space with only non-trivial first homotopy group $G$ which is finite, those
topological non-linear $\sigma$-models can realize all 3+1D bosonic topological
orders without emergent fermions, which are described by Dijkgraaf-Witten
theory with gauge group $\pi_1(K)=G$. Here, we show that the 3+1D bosonic
topological orders with emergent fermions can be realized by topological
non-linear $\sigma$-models with $\pi_1(K)=$ finite groups, $\pi_2(K)=Z_2$, and
$\pi_{n>2}(K)=0$. A subset of those topological non-linear $\sigma$-models
corresponds to 2-gauge theories, which realize and classify bosonic topological
orders with emergent fermions that have no emergent Majorana zero modes at
triple string intersections. The classification of 3+1D bosonic topological
orders may correspond to a classification of unitary fully dualizable fully
extended topological quantum field theories in 4-dimensions.
