Do roots of these polynomials approach the negative of the Euler-Mascheroni constant? Let $p_n$ be the $n$th degree polynomial that sends $\frac{k(k-1)}{2}$ to $\frac{k(k+1)}{2}$ for $k=1,2,...,n+1$. E.g., $p_2(x) = (6+13x -x^2)/6$ is the unique quadratic polynomial $p(x)$ satisfying $p(0) = 1$, $p(1) = 3$, and $p(3) = 6$. Then it appears that $p_n(x)-x$ always has precisely one negative real root, and moreover this root (as a function of $n$) appears to approach $-0.577$ as $n$ gets large. Do these roots indeed approach a limit, and is this limit indeed the negative of the Euler-Mascheroni constant?
For $n=60$, the root is about $-0.580$; for $n=120$, the root is about $-0.577$. That's as far as I've gone.
(If the question seems unmotivated, I'll give some background, though it might not help anyone find an answer. One way to think about the "values" of certain divergent series is to view them as fixed points of an associated function; e.g., for the series $1+2+4+8+\dots$, the linear map $x \mapsto 2x+1$ sends the sum of the first $n$ terms of the series to the sum of the first $n+1$ terms, and the fixed point of this map is $-1$, which is indeed the natural value to assign to the divergent series. The problem I'm posting arose from trying to analyze $1+2+3+4+\dots$ in a similar way.)

According to Mathematica your polynomials satisfy the recurrence relation
$$
(2 n+1) p(n) \left(n^2+3 n-2 x+2\right)+p(n+1) \left(-4 n^3-18 n^2+4 n x-27 n+2 x-14\right)+(2 n+3) (n+2)^2 p(n+2)=0
$$
with initial conditions $p_1=1+2x,p_2=1-x(x-13)/6$. The solution to this equation is
$$
p_n(x)=\frac{1}{4} \left(\frac{\Gamma \left(n+\frac{1}{2}\right) \cos \left(\frac{1}{2} \pi  \sqrt{8 x+1}\right) \Gamma \left(n+\frac{1}{2} \sqrt{8 x+1}+\frac{3}{2}\right) \Gamma \left(n-\frac{1}{2} \sqrt{8 x+1}+\frac{3}{2}\right) \, _4\tilde{F}_3\left(1,n+\frac{1}{2},n-\frac{1}{2} \sqrt{8 x+1}+\frac{3}{2},n+\frac{1}{2} \sqrt{8 x+1}+\frac{3}{2};n+\frac{3}{2},n+2,n+2;1\right)}{\pi }+2 \, _3F_2\left(-\frac{1}{2},\frac{1}{2}-\frac{1}{2} \sqrt{8 x+1},\frac{1}{2} \sqrt{8 x+1}+\frac{1}{2};\frac{1}{2},1;1\right)+4 x+2\right)
$$
For $n\to\infty$, Mathematica claims that this becomes
$$
p_\infty(x)=\frac{1}{2} \left(\, _3F_2\left(-\frac{1}{2},\frac{1}{2}-\frac{1}{2} \sqrt{8 x+1},\frac{1}{2} \sqrt{8 x+1}+\frac{1}{2};\frac{1}{2},1;1\right)+2 x+1\right)
$$
which solves $p(x)=x$ at around
$$
x=-0.573825523080029241015952733\dots
$$
