The idea here is to experiment with some simple iterations based on messing with the Toda (and similar) flows.

Toda Flow

To implement the Toda flow I need to split a possibly complex matrix.

SplitMatrix[B_]:=Module[
{BL = LowerTriangularize[B,-1]},
B-BL
]
SplitComplexMatrix[B_]:= Module[
{BL=LowerTriangularize[B,-1],π1,π2},
π1=(1/2)(BL+ Conjugate[BL])+ I Im[ DiagonalMatrix[Diagonal[B]]];
π2=B-π1;
{π1,π2}
]

The evolution 
	B'[t]==π1[B].B-B.π1[[B]
should be iso-spectral and should force force the lower triangular portion of
B to zero. It does seem to work.  For a real matrix the diagonal should become ordered.

Clear[F,B,t]
n=12;
TMax=300;
f[B_]:=B;
B0=RandomReal[{-1,1},{n,n}];
B0=LowerTriangularize[UpperTriangularize[B0,-1],1];
B0=B0 + Transpose[B0];
F[B_List]:= Module[
{π1},
π1=SplitComplexMatrix[B][[2]];
B.π1 - π1.B
]
BSol=B/.NDSolve[{
B'[t]==(10^2) F[B[t]],
B[0]==B0
},
B,
{t, 0, TMax},
AccuracyGoal->2,
PrecisionGoal->2
] [[1]]

The following shows that the evolution is isospectral and that the diagonal becomes ordered.

Needs["DifferentialEquations`InterpolatingFunctionAnatomy`"];
λ0=Eigenvalues[B0];
ListLogPlot[
Transpose[{
Flatten[InterpolatingFunctionGrid[BSol]],
Map[Norm[1-((Eigenvalues[#])/(λ0))]&,InterpolatingFunctionValuesOnGrid[BSol]]
}]
]

ts=Flatten[InterpolatingFunctionGrid[BSol]];
BDias=Map[Diagonal,InterpolatingFunctionValuesOnGrid[BSol]];
Fred=Table[
tt=ts[[i]];
ds=BDias[[i]];
Table[ {tt,ds[[j]]},{j, n}],
{i, Length[ts]}
];
ListLogLinearPlot[Transpose[Fred],
PlotRange->All]

OK so what we have learned.
The Toda flow really does try to set the order. In practice it does wuite well. The last pieces to resolve typically involve switching two entries which is kinda interesting.  I would guess that the knowledge of the final version could be used to advantage in a wave form relaxation algorithm.