Variance of exit time for simple symmetric random walk
(Rather stunningly, this was downvoted, thrice. Too much mathematical content, perhaps? :-) Anyway, happy reading!)
Recall that, to compute the mean exit time $E(\tau)$ of the simple random walk $(X_n)_{n\geqslant0}$ from $(a,b)$, one can use two clever martingales, namely, $(X_n)_{n\geqslant0}$ itself and $(Y_n)_{n\geqslant0}$ defined as $Y_n=X_n^2-n$. Then the identity $E(X_\tau)=0$ determines the distribution of $X_\tau$ and the identity $E(Y_\tau)=0$ shows that $E(\tau)=E(X_\tau^2)$, yielding $$E(\tau)=-ab. $$ One can also use some other clever martingales, actually a whole family of them, namely the processes $(Z^{(c)}_n)_{n\geqslant0}$, indexed by some parameter $c$, and defined as $$Z^{(c)}_n=e^{cX_n}(\cosh c)^{-n}.$$ These are indeed martingales and, expanding them as random series $$Z^{(c)}_n=\sum\limits_{k\geqslant1}\varphi_k(X_n,n)\frac1{k!}c^k,$$ one gets some random variables $M^{(k)}_n=\varphi_k(X_n,n)$ such that every $M^{(k)}=(M^{(k)}_n)_{n\geqslant0}$ is itself a martingale, hence such that $$E(\varphi_k(X_\tau,\tau))=\varphi_k(0,0).$$ Expanding $e^{cx}(\cosh c)^{-t}$ into powers of $c$, one gets $\varphi_0(x,t)=1$ (boring), $\varphi_1(x,t)=x$ and $\varphi_2(x,t)=x^2-t$ (one recognizes the martingales $M^{(1)}=X$ and $M^{(2)}=Y$ introduced above), but also the new functions $$\varphi_3(x,t)=x^3-3xt,\qquad \varphi_4(x,t)=x^4-6tx^2+3t^2+2t. $$ Then the identity $E(M^{(1)}_\tau)=0$ yields the distribution of $X_\tau$ as $P(X_\tau=a)=p$ and $P(X_\tau=b)=q$ with $$p=b/(b-a),\qquad q=-a/(b-a).$$ The identity $E(M^{(2)}_\tau)=0$ yields $E(\tau)=E(X_\tau^2)=a^2p+b^2q$ hence all this confirms that $$E(\tau)=-ab. $$ Passing (at last) to unknown territory, note that the identity $E(M^{(3)}_\tau)=0$ yields $$3E(X_\tau \tau)=E(X_\tau^3)=a^3p+b^3q=-ab(a+b),$$ hence $u=E(\tau;X_\tau=a)$ and $v=E(\tau;X_\tau=b)$ solve the system $u+v=E(\tau)$, $au+bv=E(X_\tau \tau)$, which fully determine them. More explicitely, $$u+v=-ab,\qquad 3(au+bv)=-ab(a+b).$$ And now, the identity $E(M^{(4)}_\tau)=0$ yields $3E(\tau^2)=6E(\tau X_\tau^2)-E(X_\tau^4)-2E(\tau)$ with $E(\tau X_\tau^2)=a^2u+b^2v$ and $E(X_\tau^4)=a^4p+b^4q$ hence $$3E(\tau^2)=6(a^2u+b^2v)-(a^4p+b^4q)+2ab.$$ One should mention to conclude that all this is inspired by the well-known analogous exit problem for Brownian motion, which leads to similar constructions, only neater...
One can use a difference equation to solve this problem. By conditioning on the first step you get
$$P_x(\tau=t)=\frac{1}{2} P_{x+1}(\tau=t-1)+\frac{1}{2}P_{x-1}(\tau=t-1).$$
Now multiply by $t^2$ and sum over $t$:
$$\sum_{t=1}^\infty t^2 P_x(\tau=t)=\sum_{t=1}^\infty \frac{1}{2} t^2 P_{x+1}(\tau=t-1) + \frac{1}{2} t^2 P_{x-1}(\tau=t-1).$$
This is not quite what we want yet, because we're summing $t^2$ against the probability of being $t-1$. We'd like the terms on the right side to be in terms of $t-1$. So we add and subtract appropriately, obtaining:
$$E_x[\tau^2]=\sum_{t=1}^\infty \frac{1}{2} ((t-1)^2+2(t-1)+1) P_{x+1}(\tau=t-1) + \frac{1}{2} ((t-1)^2 +2(t-1)+1) P_{x-1}(\tau=t-1).$$
Denoting $E_x[\tau^2]$ by $v(x)$ and $E_x[\tau]$ by $m(x)$ and noting that the sum of a probability distribution over all its elements is $1$, this simplifies to
$$v(x)=\frac{1}{2} v(x+1)+m(x+1)+\frac{1}{2} v(x-1) + m(x-1)+1.$$
Put another way:
$$\frac{1}{2} v(x+1)-v(x)+\frac{1}{2}v(x-1)=-m(x+1)-m(x-1)-1.$$
The left side can be thought of as a discrete approximation of the second derivative, while the right side is exactly a quadratic function (known from the OP). This suggests a quartic ansatz; substituting in this ansatz and the formula for $m(x)$ gives five equations in the five unknown coefficients. The resulting equation has a unique solution which can be found (though it is a bit tedious by hand).