Explanation of why this derivation of Schmidt decomposition works I'm following Preskill's notes and he derives the Schmidt decomposition in the following way:
Let a bipartite state be $\psi_{AB} = \sum_{i,j}\lambda_{ij}\vert i\rangle\vert j\rangle = \sum_{i} \vert i\rangle\vert \tilde{i}\rangle$, where I simply choose $\sum_j \lambda_{ij}\vert j\rangle = \vert \tilde{i}\rangle$.
I choose a set of basis vectors $\vert i\rangle$ such that the partial state is diagonal, that is $\rho_A = \sum_i p_i\vert i\rangle\langle i\vert$. But I can also obtain $\rho_A = Tr_B(\rho_{AB}) = Tr_B\sum_{i,j} \vert i\rangle\langle j\vert \otimes \vert \tilde{i}\rangle\langle \tilde{j}\vert = \sum_{ij} \langle \tilde{j}\vert\tilde{i}\rangle \vert i\rangle\langle j\vert$. The last part can be computed by explicitly writing out the trace over $B$ and using the properties of an orthonormal basis.
Thus, we have $\rho_{A} = \sum_i p_i\vert i\rangle\langle i\vert = \sum_{ij} \langle \tilde{j}\vert\tilde{i}\rangle \vert i\rangle\langle j\vert$. That is $\langle \tilde{j}\vert \tilde{i}\rangle = p_i\delta_{ij}$. Suddenly, the $\vert\tilde{i}\rangle$ are all orthogonal to each other. 
Why does choosing the basis where $\rho_A$ is diagonal also give you orthogonal vectors in $B$? This seemed to drop out of the sky for me although the math is clear. What is the physical meaning of this?
 A: 
Why does choosing the basis where $\rho_A$ is diagonal also give you orthogonal vectors in $B$? 

The answer is in the proof shown in the question. I'll write it out here in a slightly different way to try to help highlight what's happening:
Suppose that the state
$$
\psi_{AB}=\sum_n |A_n\rangle |B_n\rangle
\tag{1}
$$
is such that the reduced state
$$
\rho_A = \text{Trace}_{B}(\psi_{AB})
\tag{2}
$$
is diagonal in the $A_n$ basis. More explicitly, the reduced state is defined by
$$
\rho_A 
= \sum_k 
 \big(\sum_n |A_n\rangle \langle \hat B_k|B_n\rangle\big)
 \big(\sum_m \langle B_m|\hat B_k\rangle \langle A_m|\,\big)
\tag{3}
$$
where the vectors $|\hat B_k\rangle$ are orthonormal by definition (because we're using them to compute the trace). This implies
$$
\rho_A 
= \sum_{n,m} |A_n\rangle
 \langle B_m|B_n\rangle \langle A_m|.
\tag{4}
$$
We assumed that $\rho_A$ is diagonal in the $A_n$ basis, and the terms in the sum in (4) are all linearly independent, this is only possible if the coefficient of each individual off-diagonal term is zero:
$$
  \langle B_m|B_n\rangle = 0.
$$
Thus equation (1) is in Schmidt form.
A: Let us start from the Schmidt decomposition $|\psi\rangle = \sum s_i |a_i\rangle |b_i\rangle$.
Now consider the reduced state of $A$: $\rho_A=\sum s_i^2 |a_i\rangle\langle a_i|$.  This is, the eigenbasis of A is exactly the basis you need for the Schmidt decomposition!
Thus, if you write your state using that eigenbasis of Alice,
$$
|\psi\rangle = \sum_i |a_i\rangle \Big(\sum_j \lambda_{ij}|j\rangle\Big)\ ,
$$
the part $|\tilde b_i\rangle=\sum_j \lambda_{ij}|j\rangle$ must be equal to $s_i|b_i\rangle$, since the Schmidt decomposition is unique (modulo degeneracies).
