Why is the boundary term in the Einstein-Hilbert action, the Gibbons-Hawking-York term, generally "missing" in General Relativity courses, IMPORTANT from the variational viewpoint, geometrical setting and the needs of Black Hole Thermodynamics? Shouldn't it be also included in modern courses of General Relativity despite its global effect on the equations of motion is irrelevant (at least in the classical theory of relativistic gravity)?
[Physics] Boundary term in Einstein-Hilbert action
actionboundary-termsgeneral-relativitylagrangian-formalismvariational-principle
Related Solutions
I've never seen a paper where the calculation is performed in a manifestly covariant manner. However, I've posted a set of reference notes on my website (http://jacobi.luc.edu/notes.html) that contains the variations needed to carry out the calculation. Let me summarize the calculation here.
The action for gravity on a compact region $M$ with boundary $\partial M$ is $$I_{EH} + I_{GHY} = \frac{1}{2 \kappa^2} \int_{M}d^{d+1}x \sqrt{-g} R + \frac{1}{\kappa^2} \int_{\partial M} d^{d}x \sqrt{-h} K ~.$$ The metric on $M$ is $g_{\mu\nu}$, and $R = g^{\mu\nu} R_{\mu\nu}$ is the Ricci Scalar. The induced metric on the boundary $\partial M$ is $h_{\mu\nu} = g_{\mu\nu} - n_{\mu} n_{\nu}$, where $n^{\mu}$ is the (spacelike) unit vector normal to $\partial M \subset M$. Now consider a small variation in the metric: $g_{\mu\nu} \to g_{\mu\nu} + \delta g_{\mu\nu}$. The quantities appearing in the Einstein-Hilbert part of the action change in the following manner: $$ \delta \sqrt{-g} = \frac{1}{2} \sqrt{-g} g^{\mu\nu} \delta g_{\mu\nu}$$ $$ \delta R = -R^{\mu\nu} \delta g_{\mu\nu} + \nabla^{\mu}\left(\nabla^{\nu} \delta g_{\mu\nu} - g^{\nu\lambda} \nabla_{\mu} \delta g_{\nu\lambda} \right)$$ Thus, the change in $I_{EH}$ is $$\begin{aligned}\delta I_{EH} = & \frac{1}{2\kappa^{2}}\int_{M} d^{d+1}x \sqrt{-g} \left(\frac{1}{2} g^{\mu\nu} R - R^{\mu\nu} \right)\delta g_{\mu\nu}\\ & + \frac{1}{\kappa^2} \int_{\partial M} d^{d}x \sqrt{-h} \frac{1}{2} n^{\mu} \left(\nabla^{\nu} \delta g_{\mu\nu} - g^{\nu\lambda} \nabla_{\mu} \delta g_{\nu\lambda}\right)~,\end{aligned}$$ with the boundary term coming from the volume integral of the total derivative in $\delta R$. The variations of the quantities in the GHY term are a bit more complicated to work out, but they all basically follow from standard definitions and this result for the variation of the normal vector: $$\delta n_{\mu} = \frac{1}{2} n_{\mu} n^{\nu} n^{\lambda} \delta g_{\nu\lambda} = \frac{1}{2} \delta g_{\mu\nu} n^{\nu} + c_{\mu}~.$$ In the second equality I've introduced a vector $c_{\mu}$ that is orthogonal to $n^{\mu}$; it is given by $$c_{\mu} = - \frac{1}{2} h_{\mu}{}^{\lambda} \delta g_{\nu\lambda} n^{\nu} ~.$$ The reason I've introduced this vector is that the variation in the trace of the extrinsic curvature can be written as $$\delta K= - \frac{1}{2} K^{\mu\nu} \delta g_{\mu\nu} - \frac{1}{2} n^{\mu}\left(\nabla^{\nu} \delta g_{\mu\nu} - g^{\nu\lambda} \nabla_{\mu} \delta g_{\nu\lambda} \right) + D_{\mu} c^{\mu}$$ where $D_{\mu}$ is the covariant derivative along $\partial M$ that is compatible with the induced metric $h_{\mu\nu}$. So, the change in the GHY part of the action is $$\delta I_{GHY} = \frac{1}{\kappa^2} \int_{\partial M} d^{d}x \sqrt{-h}\left(\frac{1}{2}h^{\mu\nu} \delta g_{\mu\nu} K + \delta K \right)~.$$ Combining this with $\delta I_{EH}$ we see that the several terms cancel, leaving $$\begin{aligned} \delta I = & \frac{1}{2\kappa^2}\int_{M} d^{d+1}x \sqrt{-g}\left(\frac{1}{2} g^{\mu\nu} R - R^{\mu\nu} \right)\delta g_{\mu\nu}\\ & + \frac{1}{\kappa^2}\int_{\partial M} d^{d}x \sqrt{-h}\left(\frac{1}{2}(h^{\mu\nu} K - K^{\mu\nu})\delta g_{\mu\nu} + D_{\mu} c^{\mu} \right) ~.\end{aligned}$$ We discard the term $D_{\mu} c^{\mu}$, which is a total boundary derivative.
A long time ago, when I was in graduate school, a friend and I went to Bryce DeWitt with a similar question. He thought about it for a while, and then said to us:
"Do you know what your problem is? Too much book learning."
He explained that the answer could almost certainly be found in a book or paper, if we looked hard enough, but we would learn much more by working out the answer ourselves. This may be the best advice I ever received as a graduate student.
Why am I telling you this? Because the answer to your question about the first variation of the action (and the notes I linked to) provides you with everything you need to solve this problem. Give it a try! If you get stuck on a specific step in the derivation, come back here and ask about it -- I'll be happy to work through the details with you. But this is a great chance for you to apply what you learned from your last question.
Best Answer
Two points. First, the variation of the GHY term is discussed in detail in this post:
Explicit Variation of Gibbons-Hawking-York Boundary Term
Second, the GHY term on its own is enough to render the boundary value problem well-defined, but it is not sufficient for a physically interesting variational principle. This requires additional surface terms in the action, discussed by Regge and Teitelboim in the Hamiltonian formulation
http://adsabs.harvard.edu/abs/1974AnPhy..88..286R
and by Mann and Marolf in the Lagrangian case
http://arxiv.org/abs/hep-th/0511096