Gram-Schmidt orthogonalization for a symmetric bilinear form (source code)

= Gram-Schmidt orthogonalization for a symmetric bilinear form

For a symmetric bilinear form $B$, the same projection formula
$$
u_k=v_k-\sum_{j<k}\frac{B(v_k,u_j)}{B(u_j,u_j)}u_j
$$
produces a $B$-orthogonal basis whenever the chosen pivots $B(u_j,u_j)$ are nonzero. Reordering the input vectors can avoid a zero pivot when a suitable anisotropic vector remains.