How can gravity affect light?
Newton's law does predict the bending of light. However it predicts a value that is a factor of two smaller than actually observed.
The Newtonian equation for gravity produces a force:
$$ F = \frac{GMm}{r^2} $$
so the acceleration of the smaller mass, $m$, is:
$$ a = \frac{F}{m} = \frac{GM}{r^2}\frac{m}{m} $$
If the particle is massless then $m/m = 0/0$ and this is undefined, however if we take the limit of $m \rightarrow 0$ it's clear that the acceleration for a massless object is just the usual $a = GM/r^2$. That implies a photon will be deflected by Newtonian gravity, and you can use this result to calculate the deflection due to a massive object with the result:
$$ \theta_{Newton} = \frac{2GM}{c^2r} $$
The calculation is described in detail in this paper. The relativistic calculation gives:
$$ \theta_{GR} = \frac{4GM}{c^2r} $$
The point of Eddington's 1919 expedition was not to show that light was bent when no bending was expected, but rather to show that the bending was twice as great as expected.
One can in principle consider a Schwarzschild spacetime:
$ds^2 = -\left(1- \frac{2M}{r}\right)dt^2 + \frac{dr^2}{1- \frac{2M}{r}} + r^2 \left(d\theta^2 + \sin^2 \theta d \phi^2\right)$
The Lagrangian of geodesics is then given by:
$\mathcal{L} = \frac{1}{2} \left[- \left(1 - \frac{2M}{r}\right)\dot{t}^2 + \frac{\dot{r}^2}{1-\frac{2M}{r}} + r^2 \dot{\theta}^2 + r^2 \sin^2 \theta \dot{\phi}^2\right]$
After applying the Euler-Lagrange equations, and exploiting the fact that the S-metric is spherically symmetric and static, one obtains the orbital equation for light as (After defining $u = 1/r$) as:
$\frac{d^2 u}{d\phi^2} + u = 3 M u^2$.
It 's pretty difficult to solve this ODE. In fact, I don't think a closed-form solution exists. One can apply a perturbation approach. Defining an impact parameter $b$, one can obtain an ansatz solution to this ODE as:
$u = \frac{1}{b} \left[\cos \phi + \frac{M}{b} \left(1 + \sin^2 \phi\right)\right]$.
One can derive the following relationship:
$u\left(\frac{\pi}{2} + \frac{\delta \phi}{2}\right) = 0$,
where $\delta \phi$ is the deflection angle.
Now, finally Taylor expanding $u$ above around $\pi/2$, one can show that, in fact:
$\delta \phi = \frac{4M}{b}$,
which is the required result.
In this answer we derive the formula for the angle of deflection
$$ \theta ~=~\frac{2GM}{b}\left(\frac{1}{v_0^2} + \frac{1}{c^2}\right) +{\cal O}(M^2) \tag{1}$$
of a (massive or massless) particle in a Schwarzschild spacetime. Here $b$ is the impact parameter and $v_0$ is the asymptotic speed (which for a massless particle is $c$). The Newtonian limit is $c\to \infty$.
Sketched proof:
By spherical symmetry, we can take the orbit plane = equatorial plane: $\theta=\frac{\pi}{2}$. We start from the geodesic equations $$ \hspace{-5em}\begin{align} E~=~&n^{-1}\frac{dt}{d\lambda} , \cr n^{-1} ~:=~&1 - \frac{r_s}{r}, \cr r_s~:=~&\frac{2GM}{c^2}, \end{align}\tag{5.61/7.43/6.3.12} $$ $$ \hspace{-5em} L~=~r^2\frac{d\phi}{d\lambda}, \tag{5.62/7.44/6.3.13} $$ $$ \hspace{-5em}\begin{align} \epsilon c^2~=~&n^{-1}c^2 \left(\frac{dt}{d\lambda}\right)^2 \cr ~-~&n\left(\frac{dr}{d\lambda}\right)^2 -r^2\left(\frac{d\phi}{d\lambda}\right)^2, \end{align}\tag{5.55/7.39/6.3.10} $$ cf. Refs. 1-3. Here $$ [\lambda]=\text{Time}, \qquad [E] ~=~\text{dimensionless}, \qquad [L] ~=~\frac{\text{Length}^2}{\text{Time}}. \tag{2}$$
Massive particle: $\epsilon=1$ and $\lambda=\tau$ proper time.
Massless particle: $\epsilon=0$ and $\lambda$ not proper time.
This leads to $$\hspace{-5em} \begin{align} \underbrace{ (cE)^2 }_{\text{"energy"}} ~=~& \underbrace{ \left(\frac{dr}{d\lambda}\right)^2 }_{\text{"kinetic energy"}}\cr ~+~& \underbrace{n^{-1}\left(\frac{L^2}{r^2}+\epsilon c^2 \right) }_{\text{"effective potential energy"}}. \end{align}\tag{5.64/7.46/6.3.14} $$ The reader may ponder whether the $\epsilon$ parameter leads to a discontinuity in the angle $\theta$ of deflection (1) between the massive and massless case? We shall see below that it does not.
Q: How do we identify the constants of motion $E$ and $L$ with the observables $b$ and $v_0$ at spatial infinity $r=\infty$? A: Note that $$\begin{align}r v_{\phi}~:=~ r^2\frac{d\phi}{dt}\quad\longrightarrow&\quad bv_0~=~h~=~\frac{L}{E}\cr &\quad\text{for}\quad r~\to~\infty, \end{align}\tag{3}$$ and $$\begin{align}\frac{dr}{dt} \quad\longrightarrow&\quad v_0~=~c\frac{\sqrt{E^2-\epsilon}}{E}\cr &\quad\text{for}\quad r~\to~\infty.\end{align}\tag{4}$$ Eq. (4) means that the energy constant is $E=\gamma_0=\left(1-\frac{v_0^2}{c^2}\right)^{-1/2}$ in the massive case, and it is undetermined in the massless case.
If we define the reciprocal radial coordinate $$ u~:=~\frac{1}{r},\tag{5}$$ we get a 3rd-order polynomial$^1$ $$\begin{align} \left(\frac{du}{d\phi}\right)^2~=~&P(u)\cr ~:=~&\left(\frac{cE}{L}\right)^2 - n^{-1}\left(u^2+\frac{\epsilon c^2}{L^2}\right) \cr ~=~&r_s(u-u_+)(u-u_-)(u-u_0),\end{align}\tag{6}$$ with 3 roots $$\begin{align}u_{\pm}~=~&\pm\frac{c}{L}\sqrt{E^2-\epsilon} ~+~\frac{r_s}{2}\left(\frac{cE}{L}\right)^2~+~{\cal O}(r_s^2)\cr ~=~&\pm \frac{1}{b}~+~\frac{GM}{h^2}~+~{\cal O}(M^2) ,\end{align} \tag{7}$$ and
$$u_0~=~\frac{1}{r_s} +{\cal O}(r_s) .\tag{8}$$During the scattering process the reciprocal radial coordinate $u$ goes from 0 to the root $u_+$ and then back again to 0. The half-angle is then $$\begin{align}\phi&~~~\stackrel{(6)}{=}~\int_0^{u_+} \! \frac{du}{\sqrt{P(u)}} \cr &~~~=~\int_0^{u_+} \! \frac{du}{\sqrt{(u_+-u)(u-u_-)\left(1-r_su\right)}} \cr &\stackrel{u=u_+x}{=}~\int_0^1 \! \frac{dx}{\sqrt{(1-x)(x+\alpha)}}\left(1+\beta x\right) +{\cal O}(r_s^2) \cr &~~~=~\beta\sqrt{\alpha}+(\alpha\beta+\beta+2)\arctan\frac{1}{\sqrt{\alpha}}+{\cal O}(r_s^2) \cr &~~~=~\frac{r_s }{2b}+2\arctan\left(1+\frac{r_s }{2b}\frac{c^2}{v_0^2}\right)+{\cal O}(r_s^2)\cr &~~~=~\frac{r_s }{2b}+2\left(\frac{\pi}{4}+\frac{r_s }{4b}\frac{c^2}{v_0^2}\right)+{\cal O}(r_s^2), \end{align}\tag{9}$$ where we have defined $$\begin{align} \alpha~:=~&-\frac{u_-}{u_+}\cr ~=~&1-\frac{r_s}{L}\frac{E^2}{\sqrt{E^2-\epsilon}} +{\cal O}(r_s^2) \cr ~=~&1-\frac{r_s }{b}\frac{c^2}{v_0^2} +{\cal O}(r_s^2)\end{align} \tag{10}$$ and $$\begin{align} \beta~:=~&\frac{r_s}{2}u_+ ~=~\frac{r_s}{2}\frac{\sqrt{E^2-\epsilon}}{L}+{\cal O}(r_s^2)\cr ~=~&\frac{r_s}{2b} +{\cal O}(r_s^2). \end{align}\tag{11}$$ The constant $\beta$ vanishes in the Newtonian limit $c\to \infty$.
Finally, we can calculate the angle of deflection $$\begin{align}\theta ~=~&2\left(\phi-\frac{\pi}{2}\right)\cr ~\stackrel{(9)}{=}~&\frac{r_s}{b}\left(1 + \frac{c^2}{v_0^2}\right) +{\cal O}(r_s^2),\end{align}\tag{12}$$ which is the sought-for formula (1). $\Box$
References:
Sean Carroll, Spacetime and Geometry: An Introduction to General Relativity, 2003; Section 5.4.
Sean Carroll, Lecture Notes on General Relativity, Chapter 7. The pdf file is available here.
R. Wald, GR, 1984; Section 6.3.
--
$^1$ In the Newtonian limit $c\to \infty\leftrightarrow r_s \to 0$ the 3rd-order polynomial (6) is replaced with a 2nd-order polynomial $$ \left(\frac{du}{d\phi}\right)^2 ~=~(u_+-u)(u-u_-).\tag{13}$$ Differentiation leads to the Binet equation $$ \frac{d^2u}{d\phi^2}+u ~=~\frac{r_s}{2}\left(\frac{cE}{L}\right)^2~=~\frac{GM}{h^2} .\tag{14}$$ The solutions are hyperbolas $$u~=~\frac{GM}{h^2}(1+ e\cos(\phi\!-\!\phi_0)), \tag{15}$$ where $e>1$ is eccentricity and $\phi_0$ is a phase offset.