Robert Millikan's famous oil-drop experiments were carried out at the University of Chicago from 1909 to determine the magnitude of the charge of the electron (since May 2019, this quantity has been fixed by definition at $1.602176634 \times 10^{-19}\;\mathrm{C}$.) In a single experiment, an electrically charged oil droplet was observed to fall a known distance, $d$, between two uncharged plates at its terminal velocity, $v_g$: from the time taken, $t_g$, the droplet's radius, $a$, can be deduced. Next, a voltage was applied to the plates, inducing an electric field between them. As the droplet rises under the resulting net force, the time taken, $t_e$, for it to move back up through the same distance, $d$, can be used to deduce its total charge, $q$, which is observed to be an integer multiple of the same base value, $e$, that is: $q=Ne$.
In this example we adopt a coordinate system in which the droplet's vertical position, $z$, increases in the "up" direction. For the free-fall part of the experiment, when the droplet falls at constant terminal velocity $v_g=-d/t_g$ there is no net force on it: the sum of the gravitational and drag forces is zero: $$ F_g + F_d = 0 \quad \Rightarrow -m'g - 6\pi\eta a v_g = 0, $$ where $m' = \frac{4}{3}\pi a^3\rho' = \frac{4}{3}\pi a^3(\rho_\mathrm{oil} - \rho_\mathrm{air})$ is the effective mass of the droplet (after the mass of air it displaces is taken into account), $g = 9.803\;\mathrm{m\,s^{-1}}$ is the acceleration due to gravity in Chicago, and $\eta = 1.859 \times 10^{-5}\;\mathrm{kg\,m^{-1}\,s^{-1}}$ is the air viscosity under the experimental conditions (temperature, humidity, etc.). Rearranging, we get the following expression for the droplet radius: $$ a = \sqrt{\frac{-9\eta v_g}{2\rho'g}}. $$ When a suitable voltage is applied to the plates and the droplet moves upwards at a constant velocity $v_e=d/t_e$, the force due to the electric field is balanced by gravity and the drag force (at this new velocity): \begin{align*} & F_e + F_g + F_\mathrm{d}' = 0 \quad \Rightarrow qE + 6\pi\eta a v_g - 6\pi\eta a v_e = 0\\ \Rightarrow \quad & q = \frac{6\pi\eta a(v_e - v_g)}{E} = \frac{6\pi\eta ad}{E}\left( \frac{1}{t_g} + \frac{1}{t_e} \right) \end{align*} Each droplet (labeled A – H) was observed three times for each different charge, $q$, acquired by exposure to X-rays (up to seven experiments per droplet).
The data below, which are also provided in the file eg10-millikan-data.txt, give the time data for a number of such experiments conducted with an oil of density $\rho_\mathrm{oil} = 917.3\;\mathrm{kg\,m^{-3}}$ on a day for which $\rho_\mathrm{air} = 1.17\;\mathrm{kg\,m^{-3}}$. The magnitude of the electric field was $E = 322.1\;\mathrm{kN\,C^{-1}}$ and the distance the drops move, $d=11.09\;\mathrm{mm}$. We can use these data to estimate $e$ (assuming it is not fixed by definition) as follows.
drop expt tg te tg te tg te
A 1 13.102 46.822 12.941 46.896 13.086 46.681
A 2 12.938 86.767 13.032 86.952 13.086 86.746
A 3 13.023 61.082 12.958 60.826 12.998 60.860
A 4 12.943 86.747 12.922 86.840 13.054 86.899
B 1 11.434 56.305 11.350 56.097 11.246 56.282
B 2 11.402 75.823 11.584 75.819 11.487 76.063
B 3 11.591 44.717 11.397 44.851 11.364 44.776
B 4 11.443 75.905 11.368 75.975 11.457 76.041
B 5 11.434 75.939 11.414 75.880 11.444 75.929
B 6 11.559 75.892 11.414 75.924 11.292 75.985
B 7 11.394 44.716 11.589 44.753 11.401 44.794
C 1 16.197 100.458 16.010 100.486 16.329 100.461
C 2 16.241 47.727 16.106 47.714 16.177 47.625
C 3 16.133 37.879 16.267 37.746 16.203 37.709
C 4 16.170 64.765 16.136 64.649 16.229 64.508
D 1 16.176 38.017 16.127 37.910 16.282 38.020
D 2 16.275 38.280 16.092 38.208 16.133 38.092
D 3 16.422 48.327 16.073 48.284 16.212 48.184
D 4 16.134 38.202 16.258 38.270 16.105 38.229
D 5 16.164 102.562 16.217 102.673 16.194 102.696
E 1 12.275 55.020 12.116 54.962 12.307 54.978
E 2 12.157 54.772 12.183 54.967 12.046 55.219
E 3 12.146 55.004 12.118 54.938 12.346 54.869
E 4 12.319 43.635 12.243 43.552 12.073 43.582
F 1 14.172 61.946 14.174 61.970 14.069 61.959
F 2 14.145 90.718 13.955 90.707 14.075 90.866
F 3 14.070 62.147 14.074 61.961 14.247 61.892
F 4 14.017 61.968 14.101 61.921 14.106 62.174
G 1 9.723 50.375 9.527 50.482 9.502 50.508
G 2 9.463 63.755 9.670 63.853 9.509 63.827
G 3 9.448 63.804 9.407 63.899 9.563 63.768
G 4 9.327 63.855 9.518 63.967 9.533 63.824
H 1 13.192 73.375 13.167 73.338 13.316 73.449
H 2 13.042 42.642 13.387 42.428 13.334 42.459
H 3 13.389 42.379 13.244 42.373 13.055 42.610
H 4 13.114 73.161 13.226 73.384 13.207 73.257
H 5 13.030 73.295 13.022 73.419 13.438 73.512
First, define the necessary parameters:
eta = 1.859e-5 # air viscosity, kg.m-1.s-1
rho_air = 1.17 # air density, kg.m-3
rho_oil = 917.3 # oil density, kg.m-3
rhop = rho_oil - rho_air
g = 9.803 # acceleration due to gravity, m.s-2
d = 11.09e-3 # rise/fall distance, m
E = -322.1e3 # electric field vector (points down!)
Next, read in the data, assigning the first two columns to a MultiIndex:
In [x]: import pandas as pd
In [x]: df = pd.read_csv('eg10-millikan-data.txt', delim_whitespace=True,
index_col=[0, 1])
In [x]: df.head()
Out[x]:
tg te tg.1 te.1 tg.2 te.2
drop expt
A 1 13.102 46.822 12.941 46.896 13.086 46.681
2 12.938 86.767 13.032 86.952 13.086 86.746
3 13.023 61.082 12.958 60.826 12.998 60.860
4 12.943 86.747 12.922 86.840 13.054 86.899
B 1 11.434 56.305 11.350 56.097 11.246 56.282
Note that pandas has added a counting integer to the column names to make them distinct.
We will start with just a single droplet, taking the transpose of its data:
In [x]: dropA = df.loc['A'].T
In [x]: dropA
Out[x]:
expt 1 2 3 4
tg 13.102 12.938 13.023 12.943
te 46.822 86.767 61.082 86.747
tg.1 12.941 13.032 12.958 12.922
te.1 46.896 86.952 60.826 86.840
tg.2 13.086 13.086 12.998 13.054
te.2 46.681 86.746 60.860 86.899
We would prefer to label each row as simply 'tg' or 'te':
In [x]: dropA.index = dropA.index.str.slice(0, 2)
In [x]: dropA
Out[x]:
expt 1 2 3 4
tg 13.102 12.938 13.023 12.943
te 46.822 86.767 61.082 86.747
tg 12.941 13.032 12.958 12.922
te 46.896 86.952 60.826 86.840
tg 13.086 13.086 12.998 13.054
te 46.681 86.746 60.860 86.899
We require the average of all of the values of $t_g$ (in the absence of the electric field the droplet takes the same time to fall the distance $d$) and the average value of $t_e$ for each column (each experiment may have a different droplet charge, but the fall – rise times are measured three times for each experiment):
In [x]: tg = dropA.loc['tg'].values.mean()
In [x]: te = dropA.loc['te'].mean()
In [x]: tg
Out[x]: 13.006916666666667
In [x]: te
Out[x]:
expt
1 46.799667
2 86.821667
3 60.922667
4 86.828667
dtype: float64
Now use the value of $t_g$ to calculate the droplet's radius:
In [x]: a = np.sqrt(9*eta*d/tg/2/rhop/g)
In [x]: a
Out[x]: 2.8181654881967875e-06
or about 2.82 μm. The charge we deduce for each experiment is:
In [x]: q = 6 * np.pi * eta * a * d / E * (1/tg + 1/te)
In [x]: q
Out[x]:
expt
1 -3.340563e-18
2 -3.005663e-18
3 -3.172143e-18
4 -3.005631e-18
dtype: float64
Repeating this for all the droplets, we can add a column, q to the DataFrame df:
for drop in df.index.levels[0]:
drop_df = df.loc[drop].T
drop_df.index = drop_df.index.str.slice(0, 2)
tg = drop_df.loc['tg'].values.mean()
te = drop_df.loc['te'].mean()
a = np.sqrt(9*eta*d/tg/2/rhop/g)
q = 6 * np.pi * eta * a * d / E * (1/tg + 1/te)
df.loc[drop, 'q'] = q.values
It is now helpful to sort the droplet charges by magnitude and to plot the sorted array and its differences:
In [x]: sorted_q = sorted(-df.loc[:, 'q'])
In [x]: plt.plot(sorted_q)
In [x]: plt.ylabel('$|q|$')
In [x]: plt.twinx()
In [x]: dq = np.diff(sorted_q)
In [x]: plt.plot(dq)
In [x]: plt.ylabel(r'$\Delta q$')
In [x]: plt.show()
It certainly seems possible that the droplet charge is always a multiple of some value between $1\times 10^{-19}\;\mathrm{C}$ and $2\times 10^{-19}\;\mathrm{C}$. We can therefore estimate the value of $|e|$:
In [x]: e_estimate = dq[(dq>1.e-19) & (dq<2.e-19)].mean()
In [x]: e_estimate
Out[x]: 1.5697150510604604e-19
We can now add a column to df for the number of elementary charges we hypothesize for each experiment:
In [x]: df['N'] = (df['q'] / e_estimate).astype(int)
Considering all the data then gives us our estimate for the magnitude of the electron charge:
In [x]: (df['q']/df['N']).mean()
Out[x]: 1.5923552150386455e-19
within 1% of the defined value.