Fix bug in the case when ph maps one wf onto another one.

The old phs_symmetrization() fails to produce ph eigenstates with PH^2 = +1 when one of the wave functions is mapped exactly onto another one by PH, as in this case the PH symmetrized wf's from the two are identical. This is the case with wfs = np.eye(2) and particle_hole = sigma[1], the second vector returned is an eigenstate with eigenvalue -1 not +1. In general this can be even worse, if there are multiple such vector pairs some vectors may not be ph eigenstates at all.

The fix catches this case, projects the wf's on the subspace where ph eigenstates are not found and repeats the procedure. In the generic case it does the same calculation as before. It is guaranteed to finish, as in every iteration at least one new ph eigenstate is found. Also added a check that wfs form a linearly independent system and the spaces spanned by wfs and Pdot(wfs) are identical.

Tests are modified to clear up the random unitary generation process (using kwant.rmt). Tests are now guaranteed to feed orthonormal set of vectors to phs_symmetrization(). Add extra tests for the fully off-diagonal case where the old method fails.