I have set of N m-dimensional vectors $\{\phi_i\}$ which gradually loses mutual orthogonality in an algorithm. => I have to re-orthogonalize them every few iterations. But if I use the Gram–Schmidt algorithm for instance, the first vectors (pivot) are perfectly preserved while the the last vectors are strongly perturbed (losing information about their orientation). i.e. the result is sensitive to the order in which I pick the vectors.
I need something which tries to preserve all vectors as much as possible. Like if each has its own inertia.
I was rather thinking to calculate something like orthogonalization forces:
$\vec f_i = -k \sum_{\forall j \neq i} <\vec \phi_i|\vec \phi_j> \vec \phi_j$
($k$ being some stiffness constant ) and run some minimization e.g. gradient descent ($ \phi_i \leftarrow \phi_i+ \vec f_i dt$) or dynamical quenching, possibly even BFGS.
However this is perhaps quite inefficient as each iteration needs $O(mN^2)$ operations, and convergence is quite slow.
Is there more efficient way?
(either using less ops per iteration or with faster convergence)
... but with this property that it preserve orientation of all vectors equally well.
EDIT:
to illustrate the problem for N=3:
- start with set of vectors $\{ (1,0,0),(0,1,0),(0,0,1) \}$
- the vectors get perturbed by some process resulting in $\{(1,e_1,e_2),(e_3,1,e_4),(e_5,e_6,1) \}$, where small random numbers $e$ is error caused by the perturbation.
- Grand-Schmidt will keep direction of $\phi_1 = (\sqrt{1-e_1^2-e_1^2},e_1,e_2)$ and rotate $\phi_2,\phi_3$ acordingly. This means that the whole coordinate system get rotated by $\sin^{-1}(e_1),\sin^{-1}(e_2)$.
- This can introduces bias or drift, it does not preserve angular momentum. This is a huge problem in many physical applications.
- For example if $\phi_i$ represent exapnsion coeficients of bonding molecular orbital in atomic orbital basiset (LCAO) e.g. ($s,p_x,p_y,p_z$) this drift introduce fake forces which makes the whole molecule rotate without real physical reason.