The other parts, other than the inverse square, were clear already before Newton, or at least were easy to guess. That the force of gravity is proportionality to mass of a small object responding to the field of another comes from Galileo's observation of the universal acceleration of free fall. If the acceleration is constant, the force is proportional to the mass. By Newton's third law, the force is equal and opposite on the two objects, so you can conclude that it should be proportional to the second mass too.

The model which gives you this is if you assume that everything is made from some kind of universal atom, and this atom feels an inverse square attraction of some magnitude. If you sum over all the pairwise attractions in two bodies, you get an attraction which is proportional to the number of atoms in body one times the number of atoms in body two.

So the only part that was not determined by simple considerations like this was the falloff rate. I should point out that if you look at two sources of a scalar field, and look at the force, it is always proportional to $g_1$ times $g_2$, where $g_1$ and $g_2$ are the propensity of each source to make a field by itself. Further, if you put two noninteracting sources next to each other, this g is additive, if the field is noninteracting, essentially for the reasons described above--- the independent attractions are independent. So that the proportionality to an additive body constant you multiply over the two bodies is clear. That for gravity, the g is the mass, this was established by Galileo.

### More mathemtically

Let's call the force law between the objects $F(m_1,m_2,r)$. We know that if we put the body m_1 in free fall, the acceleration doesn't depend on the mass, so

$$ F(m_1,m_2,r) = m_1 G(m_2,r) $$

So that the mass will cancel in Newton's law to give a universal acceleration. This gives you the relation

$$ F( a m_1 , m_2, r ) = a F(m_1,m_2,r) $$

We know that if we put body 2 in free fall, the same cancellation happens, but we also know Newton's third law: $F(m_1,m_2,r)= F(m_2,m_1,r)$ so that

$$ F( m_1, a m_2, r) = a F(m_1,m_2,r) $$

So you now write

$$ F( m_1 \times 1 , m_2 \times 1, r) = m_1 F( 1, m_2\times 1 , r) = m_1 m_2 F(1,1,r) $$

And this tells you that the force is proportional to the masses times a function of r. The form of the function is undetermined.

An independent argument for the scaling is that if you consider the object m_1 as composed of two nearby independent objects of mass $m_1/2$, then

$$ F(m_1/2 , m_2 , r) + F(m_1/2 , m_2 , r) = F(m_1,m_2,r)$$

Then the same conclusion follows.

These types of scaling arguments are second nature by now, and they are automatically done by matching units. So if you have a force per unit mass, the force between two massive particles must be per unit mass 1 and per unit mass 2.

This general argument fails for direct three-body forces, where the force between 3 bodies is not decomposable as a sum of forces between the pairs bodies individually. There are no macroscopic examples, since the pairwise additivity is true for linear fields, but the force between nucleons has a 3-body component.

Just consider the water, due to gravitational attraction(which I'm not sure how effective it is in this case); water molecules like to be as close to each other as possible. This means that they like to push the bubbles as close to each other as possible. Since the air has negligible mass, its gravitational forces can be neglected compared to the water ones.

The other way to go through this reasoning is by what has been suggested in the question, i.e. assuming the bubbles have negative mass.
This solution has few steps, as following:

Say we have a huge spherical lump of water(with the radius $R$), without any bubbles inside. The gravity potential of this sphere is $\frac{-3GM^2}{5R}$, but we are not going to use that.

Now, say we don't remove a small spherical part(radius $r$ and mass $m$) of the water and replace it with a same-sized sphere with density $\rho'$(or mass $m'$) but rather add it to the current sphere(a ghost like sphere which can only interact through gravity with the world). To calculate the gravitational force acting on this sphere using Shell's theorem, we also need to know the distance from the center; assume it's $x$. Since the sphere had been in equilibrium before, the new net force will be(note the force should be proportional to the mass of each object):

$$ -\frac{GM'(x)}{x^2}m' \tag{1}$$
where $M'(x)$ is the mass of water inside a sphere with radius $x$.

- The next step is to do the same with another small sphere of water. To make relations more simplified, I will put this one at $-x$. Now the net force on the first and second sphere will be:

$$F_1=-\frac{GM'(x)}{x^2} m' - \frac{G m' m'}{(2x)^2} \\
F_2=\frac{GM'(x)}{(x)^2} m' + \frac{G m' m'}{(2x)^2}$$

Or the accelerations:

$$a_1=-\frac{GM'(x)}{x^2} - \frac{G m' }{(2x)^2} \\
a_2=\frac{GM'(x)}{(x)^2} + \frac{G m' }{(2x)^2}$$

- In the case of our problem $m' = -m$, therefore:

$$a_1=-\frac{GM'(x)}{x^2} + \frac{G m }{(2x)^2},$$

($M'(x) \gg m$)which is towards the center(and the other sphere). So it looks like the two spheres are attracting each other.

Now there are some ambiguities here:

- Do negative masses behave mathematically consistent(we can use $F=m a$), which I have assumed to be the case.
- Is this attraction duo to the big sphere of water or the other sphere? Looking at equation $(1)$, it seems to be the first case; unless my previous sign convention is wrong which simply means that the big sphere of water will try blow the bubbles away, although an attraction force between them(if they are closer than a certain distance they will attract).

Also I should point, if the bubbles get in touch; they will immediately collapse into a single bubble. This is due to the surface tension, not the gravitational effects for sure.

## Best Answer

You need to multiply by four not by two. To see why let's draw the situation:

You are assuming the situation is as shown in the top diagram. So the two $M_1$s attract each other and the two $M_2$s attract each other. Those are the forces shown by the red lines.

But you also need to include the force between $M_1$ on one side and $M_2$ on the other. Those are the green lines in the second diagram. So the total force will be:

$$ F = F_{1-1} + F_{2-2} + F_{1-2} + F_{2-1} $$

Suppose we make the mass of each of the balls $M$, so if we combine them the total mass on each side is $2M$. If we combine the masses first then calculate the force we get:

$$ F = \frac{G(2M)(2M)}{d^2} = 4\frac{GM^2}{d^2} $$

If we calculate the forces treating the masses separately we get${}^1$:

$$\begin{align} F &= F_{1-1} + F_{2-2} + F_{1-2} + F_{2-1} \\ &= \frac{GMM}{d^2} + \frac{GMM}{d^2} + \frac{GMM}{d^2} + \frac{GMM}{d^2} \\ &= 4\frac{GM^2}{d^2} \end{align}$$

So the forces are the same.

This also applies to the example of the electrostatic force that you mention.

^{${}^1$ Strictly speaking the distances $M_1 \rightarrow M_2'$ and $M_2 \rightarrow M_1'$ are slightly greater than the distance $M_1 \rightarrow M_1'$ and $M_2 \rightarrow M_2'$ but we'll assume this difference is negligibly small.}