The plot in the question must have been obtained with `stepSize = 15`, not `stepSize = 2`.  Using the latter value gives a smooth plot,

[![enter image description here][1]][1]

The computation takes about 78 sec on my PC.  To address the specific issue in the question, the run time can be reduced by two orders of magnitude using

    Block[{stepSize = 2, end = TMax, tt, rd}, 
        tSolpbc = Table[uSolpbc[tt, 0], {tt, 0, end, stepSize}]; 
        rd = ParallelTable[
            UnitStep[0.01 - Norm[tSolpbc[[nt]] - tSolpbc[[nτ]], 2]], 
            {nt, end/stepSize}, {nτ, end/stepSize}]; MatrixPlot[rd]]

which produces the same plot.  Evidently, most of the run time used in the original computation was consumed by computing `uSolpbc` `(end/stepSize)^2` times.  The revised computation computes it only `end/stepSize` times.

  [1]: https://i.sstatic.net/adb7C.png