$$\begin{align} &\text{Power collected on earth} \\ &= \int I\ dA \\ &= \int_0^{R_e} \frac{P_0}{4\pi[(R_0 + R_e - \sqrt {R_e^2 - R^2})^2 + R^2]} 2\pi R\ dR \\ &= \frac{P_0}{2} \int_0^{R_e} \frac{R}{ 2R^2 - 2(R_0+R_e)\sqrt{R_e^2 - R^2} + R_e^2 + (R_0 + R_e)^2}\ dR \end{align}$$
The only variable here is R. But would like the answer to contain all the other constants, so I may proceed with further manipulations later. Thank you, and hope to get some help from the mathematicians present here.
For those who are interested in the physics derivation, here it is, step by step (upon request by @hardmath):
Consider the Earth located at a distance of Ro from the sun. Consider the sun as a point source of power, emitting power of $P_0$.
Intensity at every point on the surface of the earth is not the same, as each point on the earth is located at a different distance from the sun. Only concentric rings are located the same distance from the sun, and thus will capture the same intensity.
Intensity expression:
Let the radius of a concentric ring be $R$. Let the depth of the concentric ring from the surface of earth be $x$.
Intensity, $$I = \frac{P_0}{4\pi[(R_0 + x)^2 + R^2)]} \tag{1}$$
Finding an expression for $x$ in terms of $R$:
$$R_e^2 = R^2 + (R_e - x)^2,$$ where $R_e$ is the radius of the Earth
$$x = R_e - \sqrt {R_e^2 - R^2} \tag{2}$$
Subst (2) into (1): $$I = \frac{P_0}{4\pi[(R_0 + R_e - \sqrt {R_e^2 - R^2})^2 + R^2]} \tag{3}$$
Area expression:
Next, I find an expression for area of an infinitesimally thin concentric ring, $dA$. I can slice the concentric ring once, to make it into a rectangle.
Thus, area of concentric ring, $$dA = 2\pi R\ dR \tag{4}$$
Hence, the power collected can be calculated using the expression above.
For the sake of typing and also to make the solution look neat, let $$a=2(R_0+R_e)$$ $$b=R_e$$ $$c=R_e^2+(R_0+R_e)^2$$
So, the required integral becomes $$\int_{0}^{R_e}\frac{R\space dR}{2R^2-a\sqrt{b^2-R^2}+c}$$
Substitute $b^2-R^2=t^2$.
$\Rightarrow R\space dR=-t\space dt$
At $R=0$, we have $t=b=R_e$
At $R=R_e$, we have $t=b^2-R_e^2=0$. \begin{align} \Rightarrow\int_{0}^{R_e}\frac{R\space dR}{2R^2-a\sqrt{b^2-r^2}+c} &= \int_{R_e}^{0}\frac{-t\space dt}{2(b^2-t^2)-at+c}\\ &=\int_{R_e}^{0}\frac{t\space dt}{2t^2+at-(2b^2+c)}\\ \end{align}
This is a pretty standard integral of the form $\int\frac{Linear}{Quadratic}$ which we solve by expressing $\text{Numerator}=l(\text{Derivative of Denominator})+m$ where $l$ and $m$ are constants. This process reduces the original integral into two simpler integrals which can be solved easily.