Disclaimer 1: Though the author of the OP was prompted to use the proper format for the question this has been neglected for a significant amount of time. While I do not think that the author should receive any help due to neglecting to play by the rules, this is a very famous paper and perhaps other users can benefit in the future, and I think that should be the main concern of the site.
Disclaimer 2: Since the edit I did was too extreme for the taste of some people, and apparently prefer the unedited version, let's give some background so that anyone who's really interested in this can read in the standards (??) of the site.
Consider a gravitational theory with a dilaton field and the action that is given by Charged black holes in string theory
\begin{align}
S = \int d^4 x \sqrt{-g} (- R + 2 (\nabla \phi)^2 + e^{-2 \phi} F^2)
\end{align}
note that you can find an erratum here; erratum to the paper.
The Einstein field equations are shown to be
\begin{align}
\begin{split}
&\nabla_{\mu} (e^{-2 \phi} F^{\mu \nu}) = 0 \\
&\nabla^2\phi + \frac{1}{2} e^{-2 \phi} F^2 =0 \\
&R_{\mu \nu} = 2 \nabla_{\mu} \phi \nabla_{\nu} \phi + 2 e^{-2\phi} F_{\mu \rho} F_{\nu}^{\rho} - \frac{1}{2} g_{\mu \nu} e^{-2 \phi} F^2
\end{split}
\end{align}
Regarding the computational part of the OP. Before the edit there was a hand-written computation, which was not neat and I could not read it.
I am giving below some results regarding the variations relevant to the problem.
\begin{align}
\delta(\sqrt{-g}) = - \frac{1}{2} \sqrt{-g} g^{\mu \nu} \delta g_{\mu \nu}
\end{align}
having used $\delta g = g^{\mu \nu} \delta g_{\mu \nu}$.
We also have,
\begin{align}
\begin{split}
\delta R &= (\delta g^{\mu \nu}) R_{\mu \nu} + g^{\mu \nu} \delta R_{\mu \nu}\\
&= (\delta g^{\mu \nu}) R_{\mu \nu} + g^{\mu \nu} (\nabla_{\rho} (\delta \Gamma^{\rho}_{\mu \nu})-\nabla_{\nu} (\delta \Gamma^{\sigma}_{\sigma \mu})) \\
&= (\delta g^{\mu \nu}) R_{\mu \nu} + \nabla_{\rho} (g^{\mu \nu} \delta \Gamma^{\rho}_{\mu \nu} - g^{\mu \rho} \delta \Gamma^{\sigma}_{\sigma \mu})
\end{split}
\end{align}
where in the second step we have used the Palatini identity and in the third the metric compatibility and a relabeling of the indices.
Moreover,
\begin{align}
\begin{split}
\delta (g^{\mu \nu} \nabla_{\mu} \phi \nabla_{\nu} \phi) &= \delta (g^{\mu \nu}) \nabla_{\mu} \phi \nabla_{\nu} \phi + 2 \nabla^{\mu} \delta \phi \nabla_{\mu} \phi \\
&= \delta (g^{\mu \nu}) \nabla_{\mu} \phi \nabla_{\nu} \phi - 2 (\nabla^2 \phi) \delta \phi +2 \nabla^{\mu} (\delta \phi \nabla_{\mu} \phi)
\end{split}
\end{align}
where we have used integration by parts to obtain the final line and of course $\nabla^2 = \nabla^{\mu} \nabla_{\mu}$.
Finally, we have
\begin{align}
\begin{split}
\delta(e^{-2 \phi} F^2) &= - 2 e^{-2 \phi} F^2 \delta \phi + e^{-2\phi} \delta F^2 \\
&= -2 e^{-2 \phi} F^2 \delta \phi + 2 e^{-2 \phi} F_{\nu \rho} F_{\mu}^{\rho} \delta g^{\mu \nu} - 4 (\nabla_{\sigma} e^{-2 \phi} F^{\sigma \beta}) \delta A_{\beta} + 4 \nabla_{\sigma} (e^{-2 \phi} F^{\sigma \beta} \delta A_{\beta})
\end{split}
\end{align}
In all of the above relations, the last term is a total divergence. Now, one is at a position to derive the field equations, then take the trace of the $R_{\mu \nu}$ equation and re-write the Ricci scalar in order to obtain the form given in the paper, if I am not mistaken.
Finally, regarding the $\phi \rightarrow 0$ limit. After taking the limit, one can use the first two relations and show the Einstein field equations in a vacuum, again if I have not made any stupid mistakes.