I am trying to solve the kdv equation numerically.

Whatever initial condition I give, I finally find the solution to break into some high amplitude solitonic parts moving to the right, and some very small amplitude parts moving to the left. Is there any way to analytically show that the initial function, under evolution, will always break into two parts, moving to right and left?
"It can be shown that any sufficiently fast decaying smooth solution will eventually split into a finite superposition of solitons travelling to the right plus a decaying dispersive part travelling to the left. This was first observed by Zabusky & Kruskal (1965) and can be rigorously proven using the nonlinear steepest descent analysis for oscillatory Riemann–Hilbert problems."
Wikipedia says this about the asymptotes. Is there any way to prove this? I do not have any background in the method(steepest descent analysis of pde) wikipedia says, or on how they are done actually. Any preliminary help would be helpful.
If you find an invariant family of infinite groups and tease out the stabilizers of those groups, the PDE will fall apart into said stabilizers (because they partition the manifold). You can then use the Lie algebra between the stabilizers to solve the equation. What you are finding are the separatrices that lie between the singularities and saddle points in the manifold, which do not change under group transformation. Said separatrices are your asymptotic solutions. From scratch it's a complicated thing, but happily there is an algorithm:
$$\frac{\partial \nu}{\partial \tau}+ \frac{\nu}{1}\frac{\partial \nu}{\partial \xi}+\frac{\delta^2}{1}\frac{\partial^3 \nu}{\partial \xi^3}=0$$
Count the number of variables. If it's in the numerator, it's positive. If in the denominator, negative. Calculate this for each monomial and set them equal to each other, as shown.
$$\nu-\tau=2\nu-\xi=2\delta+\nu-3\xi$$
$\nu_1-\tau=2\nu_1-\xi \rightarrow \nu_1=\xi-\tau$ Thus, $\nu_1(\xi,\tau)=\frac{C\xi}{\tau}$Taking appropriate derivatives and substituting them back into KdV gives the value of C=1. Therefore the first asymptotic solution is $$\nu_1(\tau,\xi)=\frac{\xi}{\tau}$$
$2\nu_2-\xi=2\delta+\nu_2-3\xi \rightarrow \nu_2=2\delta-2\xi$ Thus, $\nu_2(\delta,\xi)=\frac{D\delta^2}{\xi^2}$ Substituting derivatives into KdV gives the value of D=-12. Thus the second asymptotic solution is $$\nu_2(\delta, \xi)=\frac{-12\delta^2}{\xi^2}$$ These solutions are easily checked. Interestingly, if you put them together as $\nu=\frac{\xi}{\tau}-\frac{12\delta^2}{\xi^2}$ they ALMOST make a complete solution, off by only one value of $\frac{12\delta^2}{\xi^2\tau}$.
I doubt if an instructor will accept a group theory application for the solution of your problem, since it is unlikely that they know anything about it. It came out of efforts to solve superconductor stability problems in the Tokamak, and isn't widely known or used outside of the Magnetics Division at Y-12. For more information I suggest Lawrence Dresner's books on the subject (who was my instructor) or the original texts by Sophus Lie published in 1884. The algorithm was my contribution, but it's obvious once you understand the principles involved.