Weird matrix identity I stumbled upon a weird equality that I can't seem to explain:
For a $N\times N$ (symmetric) matrix $\mathbf{X}$ that satisfies the eigenvalue problem
$$\mathbf{X} \mathbf{v}_k = \lambda_k \mathbf{v}_k$$
and a $N\times1$ matrix $\mathbf{u}$ that is split in 2 (elements in the first half are equal to 1, equal to -1 for the second half), i.e.
$$\mathbf{u}^T := \frac{1}{\sqrt{N}}\bigl[\underbrace{1 ...\ 1}_{\times N/2}\ \underbrace{-1\ ....-1}_{\times N/2}\bigr]$$
the following (apparantly) holds
$$ \mathbf{u}^T(z \mathbf{I} - \mathbf{X})^{-1} \mathbf{u} = \sum_{k=1}^N \frac{(\mathbf{u}^T \mathbf{v}_k)^2}{z-\lambda_k}.$$
$z$ is complex but vanishingly close to the real axis.

I'm guessing it has something to do with the resolvent
$$ (z \mathbf{I} - \mathbf{X})^{-1} \equiv \sum_{k=0}^\infty \frac{\mathbf{X}^k}{z^{k+1}}$$
but I can't figure it...
 A: You have an unstated assumption that $v_k$ are ortho-normal (ortho is ensured if the matrix has distinct eigenvalues). One can always choose $v_k$ so they are ortho-normal.
That said, we can write $u$ as
$$ 
u = \sum_1^k \sigma_k v_k
$$
where $$ \sigma_k = u^T v_k$$
Hence
$$ \begin{align} u^T (z I - A)^{-1} u &= u^T\sum_1^k \sigma_k(z I - A)^{-1} v_k \\
&=  \sum_1^k \sigma_k \frac{1}{z-\lambda_k} u^T v_k \\
&=  \sum_1^k  \frac{\left(u^T v_k\right)^2 }{z-\lambda_k} \\
\end{align}
$$
The above result is true for any real vector $u$.
A: I think the identity
$\mathbf{u}^T(z \mathbf{I} - \mathbf{X})^{-1} \mathbf{u} = \sum_{k=1}^N \frac{(\mathbf{u}^T \mathbf{v}_k)^2}{z-\lambda_k} \tag{1}$
holds for any $\mathbf u$, provided we take the $\mathbf v_k$ to be normalized to unity, which always a possible choice in the present context.  The matrix $\mathbf X$, being symmetric, is possessed of an orthonormal eigenbasis $\mathbf v_k$, and we can expand $\mathbf u$ in this basis, to wit:
$\mathbf u = \sum_1^N (\mathbf u^T \mathbf v_k) \mathbf v_k; \tag{2}$
now observe that $z$, being complex, satisfies $z \ne \lambda_k$ for all $k$, $1 \le k \le N$, whence $1 /(z - \lambda_k)$ is well-defined.  This being the case, we have, for each $\mathbf v_k$,
$(z - \mathbf X) \mathbf v_k = (z - \lambda_k) \mathbf v_k, \tag{3}$
which immediately leads to
$\mathbf v_k = (1 / (z - \lambda_k)) (z - \mathbf X)\mathbf v_k.  \tag{4}$
Substituting (4) into (2) yields
$\mathbf u = \sum_1^N (\mathbf u^T \mathbf v_k) (1 / (z - \lambda_k)) (z - \mathbf X)\mathbf v_k = (z - \mathbf X) \sum_1^N (\mathbf u^T \mathbf v_k)(1 / (z - \lambda_k)) \mathbf v_k, \tag{5}$
and since $z - \mathbf X$ is nonsingular by virtue of the fact that $z$ cannot be a (real) eigenvalue of $X$, we may write
$(z - \mathbf X)^{-1} \mathbf u = \sum_1^N (\mathbf u^T \mathbf v_k)(1 / (z - \lambda_k)) \mathbf v_k; \tag{6}$
taking the product of each side with $\mathbf u^T$ now gives
$\mathbf u^T (z - \mathbf X)^{-1} \mathbf u = \sum_1^N (\mathbf u^T \mathbf v_k)(1 / (z - \lambda_k)) (\mathbf u^T \mathbf v_k) =  \sum_1^N (1 / (z - \lambda_k)) (\mathbf u^T \mathbf v_k)^2, \tag{7}$
the desired result.
Note that if the $\mathbf v_k$ are not normalized, (2) would fail, and the ensuing equations would (presumably) have to be modified.
Hope this helps.  Holiday Greetings,
and as always,
Fiat Lux!!!
