fluid dynamics - How Reynolds number was derived?

I'm studying fluid dynamics and recently the formula $Re=\frac{\rho vd}{\eta}$ was presented to me. I'm curious to know how Reynolds came up with this relations between this different variables.

Did $Re=\frac{\rho vd}{\eta}$ result from the formula $Re = \frac{\text{Inertial Forces}}{\text{Viscous Forces}}$ or did this last equation came up as an intuition/ physical interpretation after the Reynolds number was first discovered?

I tried to find the history behind Reynolds "scientific procedure", how he found the number, but I wasn't successful.

Answer

There's no magic behind it. It was done by non-dimensionalizing the momentum equation in the Navier-Stokes equations.

Starting with:

$$\frac{\partial u_i}{\partial t} + u_j\frac{\partial u_i}{\partial x_i} = -\frac{1}{\rho}\frac{\partial P}{\partial x_i} + \nu \frac{\partial^2 u_i}{\partial x_i x_j}$$

which is the momentum equation for an incompressible flow. Now you non-dimensionalize things by choosing some appropriate scaling values. Let's look at just the X-direction equation and assume it's 1D for simplicity. Introduce $\overline{x} = x/L$, $\overline{u} = u/U_\infty$, $\tau = tU_\infty/L$, $\overline{P} = P/(\rho U_\infty^2)$ and then substitute those into the equation. You get:

$$\frac{\partial U_\infty \overline{u}}{\partial \tau L/U_\infty} + U_\infty\overline{u}\frac{\partial U_\infty \overline{u}}{\partial L \overline{x}} = - \frac{1}{\rho}\frac{\partial \overline{P}\rho U_\infty^2}{\partial L\overline{x}} + \nu \frac{\partial^2 U_\infty \overline{u}}{\partial L^2 \overline{x}^2}$$

So now, you collect terms and divide both sides by $U_\infty^2/L$ and you get:

$$\frac{\partial \overline{u}}{\partial \tau} + \overline{u}\frac{\partial \overline{u}}{\partial \overline{x}} = -\frac{\partial \overline{P}}{\partial \overline{x}} + \frac{\nu}{U_\infty L}\frac{\partial^2 \overline{u}}{\partial \overline{x}^2}$$

Where now you should see that the parameter on the viscous term is $\frac{1}{Re}$. Therefore, it falls out naturally from the definitions of the non-dimensional parameters.

The intuition

There's some other ways to come up with it. The Buckingham Pi theorem is a popular way (demonstrated in Floris' answer) where you collect all of the units in your problem in this case $L, T, M$ and find a way to combine them into a number without dimension. There is one way to do that, which ends up being the Reynolds number.

The interpretation of inertial to viscous forces comes from looking at the non-dimensional equation. If you inspect the magnitude of the terms, namely the convective (or inertial term) and the viscous term, the role of the number should be obvious. As $Re \rightarrow 0$, the magnitude of the viscous term $\rightarrow \infty$, meaning the viscous term dominates. As $Re \rightarrow \infty$, the viscous term $\rightarrow 0$ and so the inertial terms dominates. Therefore, one can say that the Reynolds number is a measure of the ratio of inertial forces to viscous forces in a flow.