Derivation of the behavior of solutions to $\tan x = x$ This question is related to Chapter IV, Note IV.36 of Flajolet & Sedgewick's Analytic Combinatorics, and The question: Sum of the squares of the reciprocals of the fixed points of the tangent function as well:
Let $x_k$ be the $k^{th}$ positive root of the equation $\tan z = z$.  Then, the sum $S_r = \sum_k x_k^{-2r}$ are rational numbers for $r \ge 1$.  For instances, $S_1 = 1/10$, $S_2 = 1/350$, $S_3 = 1/7875$ (from Note IV.36, pp. 269.)
I follow the approach in the book for Bernoulli numbers (Chapter IV, IV.6.1, pp. 268.)  Consider the function:
\begin{align}
f(z) = \frac{1}{\tan z - z} \tag{1}
\end{align}
It is obvious that $\{x_k\}$ are poles of $f(z)$, with $k \in \mathbb{Z} \backslash \{0\}$.  Further,
\begin{align}
Res[f(z);z=x_k] &= \frac{1}{\frac{d}{dz}(\tan z - z)|_{z=x_k}} \\
&= \frac{1}{\tan^2 x_k} \\
&= \frac{1}{x_k^2} \tag{2}
\end{align}
Hence,
\begin{align}
\frac{1}{\tan z - z} \sim \frac{1}{x_k^2}\cdot \frac{1}{z-x_k}, \text{ for } z \to x_k \tag{3}
\end{align}
By the Cauchy's coefficient formula, we have
\begin{align}
f_n = [z^n]f(z) &= \int_C \frac{f(z)}{z^{n+1}} dz \tag{4} \\
&= -\sum_{k \in \mathbb{Z} \backslash \{0\}} \frac{1}{x_k^2}\cdot \frac{1}{x_k^{n+1}}
\end{align}
where $C$ is a contour encircling all poles.  Since $f(z)$ is an odd function, the poles are in pairs like $\pm x_k$.  Then, $f_n = 0$, when $n$ is even.  As a result,
\begin{align}
f_{2n-1} &= -2 \sum_{k=1}^\infty x_k^{-2(n+1)}, \text{ for } n \ge 1 \tag{5}
\end{align}
Therefore, the following relation is established:
\begin{align}
S_r = \sum_{k=1}^\infty x_k^{-2r} = -\frac{1}{2}f_{2r-3}, \text{ for } r > 1 \tag{6}
\end{align}
From the expansion of (1), that is,
\begin{align}
\frac{1}{\tan z - z} = \frac{3}{z^3} - \frac{6}{5z} - \frac{1}{175}z - \frac{2}{7875}z^3-\frac{37}{3031875}z^5+O(z^6) \tag{7}
\end{align}
It shows that $S_2$ and $S_3$ calculated from (6) are correct.  But $S_1$ cannot be calculated from (6) since $f_{-1}$ is related to the pole at $z=0$, which is excluded from (6).
Then, for $S_1$, I follow the approach in an answer to The question to consider the integral in the region including $z=0$.  That is,
\begin{align}
\int_C f(z)dz = 2i\pi \sum_{k \in \mathbb{Z}} Res[f(z); z=x_k] \tag{8}
\end{align}
From the expansion in (7), we have $Res[f(z);z=0]=-\frac{6}{5}$.  Together with (2), (8) becomes
\begin{align}
\int_C f(z)dz = 2i\pi \left(-\frac{6}{5} + 2\sum_{k=1}^\infty x_k^{-2}\right) \tag{9}
\end{align}
Since it is found that (details in an answer to The question)
\begin{align}
\int_C f(z)dz = -2i\pi \tag{10}
\end{align}
(9) becomes
\begin{align}
-\frac{6}{5} + 2\sum_{k=1}^\infty x_k^{-2} = -1
\end{align}
Or
\begin{align}
S_1 = \sum_{k=1}^\infty x_k^{-2} = \frac{1}{10} \tag{11}
\end{align}
This result is well known as proved in The question.
It looks the flow of logic leading to the results (6) and (11) are both correct.  However, the point I don't understand is: both (4) and (9) include an integral around the same contour and both give the same result ($-2i\pi$).  Nevertheless, by Cauchy's residue theorem, they equal to the sum of different sets of residues - one with the pole at $0$ (9) while the other without (4).  The difference in residue set is due to the difference in the regions in concern ($\mathbb{C}$ vs. $\mathbb{C} \backslash \{0\}$.)  But why the contour integrals of a function at these different regions give the same result?
Could you please indicate what is wrong in the above logic?  Thank you.

You don't say precisely what contour $C$ you are using.  For definiteness, I'm going to assume it's the square $C_m$ with side length $2m\pi$ ($m$ a positive integer), center 0, and sides parallel to the axes given in robjohn's answer to your linked question.  Then by the residue theorem, for all nonnegative integers $s$,
$$\int_{C_m} \frac{f(z)}{z^s} \, dz =2\pi i\left(f_{s-1}+\sum_{\hbox{$j\ne 0$, $x_j$ inside $C_m$}} \frac{1}{x_j^{s+2}}\right),$$
where $x_{-1}$, $x_{-2}$, $\dots$ are the negative real roots of $\tan x=x$, and $f_k=[z^k]f(z)$.
So, since the contour is the same, regardless of $s$, the pole at 0 is always included; it corresponds to $f_{s-1}$.  The difference between the $s=0$ and $s>0$ cases is in the evaluation of the integral.  As robjohn showed in his answer, for large $m$, $f(z)$ is approximately $-z^{-1}$ on $C_m$, so $f(z)/z^s$ is approximately $-z^{-(s+1)}$ and
$$\lim_{m\to\infty} \int_{C_m} \frac{f(z)}{z^s} \, dz =\left\{\begin{array}{cl}
-2\pi i, & \hbox{if $s=0$,}\\
0, & \hbox{if $s>0$,}\end{array}\right.$$
so, letting $m\to\infty$,
$$\sum_{j\ne 0} \frac{1}{x_j^{s+2}} = -f_{s-1}-\left\{\begin{array}{cl}
1, & \hbox{if $s=0$,}\\
0, & \hbox{if $s>0$.}\end{array}\right.
$$
