One realization of such a phenomenon can be achieved via the Kortewegde Vries (KdV) equation. It is a nonlinear thirdorder partial differential equation for a function $u\left(t\,x\right)$ of the form
>

$\mathrm{KdV}\u2254\frac{\partial}{\partial t}u\left(tcomma;x\right)+6u\left(tcomma;x\right)\cdot \frac{\partial}{\partial x}\mathit{}u\left(tcomma;x\right)plus;\frac{{\partial}^{3}}{{\partial}^{\phantom{\rule[0.0ex]{0.2em}{0.0ex}}}{x}^{3}}u\left(tcomma;x\right)=0$

${\mathrm{KdV}}{\u2254}\frac{{\partial}}{{\partial}{t}}\phantom{\rule[0.0ex]{0.4em}{0.0ex}}{u}\left({t}{\,}{x}\right){+}{6}{}{u}\left({t}{\,}{x}\right){}\left(\frac{{\partial}}{{\partial}{x}}\phantom{\rule[0.0ex]{0.4em}{0.0ex}}{u}\left({t}{\,}{x}\right)\right){+}\frac{{{\partial}}^{{3}}}{{\partial}{{x}}^{{3}}}\phantom{\rule[0.0ex]{0.4em}{0.0ex}}{u}\left({t}{\,}{x}\right){=}{0}$
 (1.1) 
with $t$ being the time and $x$ the space coordinate.
A specific type of solutions to this equation are solitary waves. A single soliton is described by the solution of (1.1) of the form
>

$\mathrm{Soliton}\u2254v\cdot {\mathrm{sech}}^{2}\left(\sqrt{\frac{v}{2}}\cdot \left(x2\cdot v\cdot t\mathrm{x\_\_0}\right)\right)\:$

as can be easily verified by pdetest
>

$\mathrm{pdetest}\left(u\left(t\,x\right)\=\mathrm{Soliton}\,\mathrm{KdV}\right)$

returning $0$, i.e. the KdVequation is annihilated by the soliton solution. This solution contains two constants: $v$ corresponding to the amplitude of this wave and $\mathrm{x\_\_0}$ that is related to a spatial shift (a third one is already fixed in such a way that $\mathrm{lim\_\_x\→\±\∞}\=0$). Note that the width and the velocity of this traveling wave are related to its amplitude $v$ and are not independent parameters. The following plot animates the propagation of a soliton with adjustable $v$.
Note that the soliton maintains its shape while propagating. In many other models, wave packets are a superposition of plane waves of different wavelengths. If their propagation velocity depends on their frequency, then the wave packet disperses as it travels through space. In the KdV model nonlinear effects exactly cancel these dispersive effects and the waves maintain their shape.
Due to the nonlinearity of the KdV equation, the superposition of two solitons is not necessarily another solution. Indeed, you can easily check using $\mathrm{pdetest}$ that, for example, $2\cdot \mathrm{Soliton}$ is not a solution. Instead, there exists a method to obtain multisoliton solutions called classical inverse scattering method, see, for example, Dunajski. Applying this method, the twosoliton solution turns out to be given by
>

$\mathrm{TwoSolitons}\u22542\left(\mathrm{v\_\_1}\mathrm{v\_\_2}\right)\cdot \frac{\mathrm{v\_\_1}\cdot {\mathrm{sech}}^{2}\left(\sqrt{\frac{\mathrm{v\_\_1}}{2}}\left(x2\cdot \mathrm{v\_\_1}\cdot t\mathrm{x\_\_1}\right)\right)\+\mathrm{v\_\_2}\cdot {\mathrm{csch}}^{2}\left(\sqrt{\frac{\mathrm{v\_\_2}}{2}}\left(x2\cdot \mathrm{v\_\_2}\cdot t\mathrm{x\_\_2}\right)\right)}{{\left(\sqrt{2\cdot \mathrm{v\_\_1}}\mathrm{tanh}\left(\sqrt{\frac{\mathrm{v\_\_1}}{2}}\left(x2\cdot \mathrm{v\_\_1}\cdot t\mathrm{x\_\_1}\right)\right)\sqrt{2\cdot \mathrm{v\_\_2}}\mathrm{coth}\left(\sqrt{\frac{\mathrm{v\_\_2}}{2}}\left(x2\cdot \mathrm{v\_\_2}\cdot t\mathrm{x\_\_2}\right)\right)\right)}^{2}}\:\phantom{\rule[0.0ex]{0.0em}{0.0ex}}$

which can easily be verified using $\mathrm{pdetest}$ (note that the evaluation takes a few seconds)
>

$\mathrm{pdetest}\left(u\left(t\,x\right)\=\mathrm{TwoSolitons}\,\mathrm{KdV}\right)$

In the following you can study this solution and identify important properties of solitary waves.