# Brownian Motion Programming Assignment

## Creative Programming Assignments

Below are links to a number of creative programming assignments that we've used at Princeton. Some are from COS 126: Introduction to Computer Science; others are from COS 226: Data Structures and Algorithms. The main focus is on scientific, commercial, and recreational applications. The assignments are posed in terms of C or Java, but they could easily be adapted to C++, C#, Python, or Fortran 90.

Assignment | Description | Concepts | Difficulty |

SCIENTIFIC COMPUTING | |||
---|---|---|---|

Guitar Hero [checklist] | Simulate the plucking of a guitar string using the Karplus-Strong algorithm. | objects, ring buffer data type, simulation | 5 |

Digital Signal Processing [checklist] | Generate sound waves, apply an echo filter to an MP3 file, and plot the waves. | data abstraction, arrays | 5 |

Percolation [checklist] | Monte Carlo simulation to estimate percolation threshold. | union-find, simulation | 5 |

Global Sequence Alignment [checklist] | Compute the similarity between two DNA sequences. | dynamic programming, strings | 5 |

N-Body Simulation [checklist] | Simulate the motion of N bodies, mutually affected by gravitational forces, in a two dimensional space. | simulation, standard input, arrays | 3 |

Barnes-Hut [checklist] | Simulate the motion of N bodies, mutually affected by gravitational forces when N is large. | quad-tree, analysis of algorithms, data abstraction | 8 |

Particle Collision Simulation | Simulate the motion of N colliding particles according to the laws of elastic collision. | priority queue, event-driven simulation | 7 |

Atomic Nature of Matter [checklist] | Estimate Avogadro's number using video microscopy of Brownian motion. | depth-first search, image processing, data abstraction, data analysis | 8 |

Root Finding [checklist] | Compute square roots using Newton's method. | loops, numerical computation | 2 |

Cracking the Genetic Codes [checklist] | Find the genetic encoding of amino acids, given a protein and a genetic sequence known to contain that protein. | strings, file input | 5 |

RECREATION | |||

Mozart Waltz Generator | Create a two-part waltz using Mozart's dice game. | arrays | 3 |

Rogue [checklist] | Given a dungeon of rooms and corridors, and two players (monster and rogue) that alternate moves, devise a strategy for the monster to intercept the rogue, and devise a strategy for the rogue to evade the monster. | graph, breath first search, depth first search, bridges | 8 |

8 Slider Puzzle [checklist] | Solve Sam Loyd's 8 slider puzzle using AI. | priority queue, A* algorithm | 5 |

GRAPHICS AND IMAGE PROCESSING | |||

Mandelbrot Set [checklist] | Plot the Mandelbrot set. | functions, arrays, graphics | 3 |

H-tree [checklist] | Draw recursive patterns. | recursion, graphics | 3 |

Sierpinski Triangle [checklist] | Draw recursive patterns. | recursion, graphics | 3 |

Collinear Points [checklist] | Given a set of Euclidean points, determine any groups of 4 or more that are collinear. | polar sorting, analysis of algorithms | 4 |

Smallest Enclosing Circle [checklist] | Given a set of Euclidean points, determine the smallest enclosing circle. | computational geometry, randomized algorithm | 8 |

Planar Point Location [checklist] | Read in a set of lines and determine whether two query points are separated by any line. | computational geometry, binary tree | 6 |

COMBINATORIAL OPTIMIZATION | |||

Small World Phenomenon | Use the Internet Movie Database to compute Kevin Bacon numbers. | graph, breadth-first search, symbol table | 7 |

Map Routing | Read in a map of the US and repeatedly compute shortest paths between pairs of points. | graph, Dijkstra's algorithm, priority queue, A* algorithm. | 7 |

Bin Packing | Allocate sound files of varying sizes to disks to minimize the number of disks. | priority queue, binary search tree, approximation algorithm | 5 |

Traveling Salesperson Problem | Find the shortest route connecting 13,509 US cities. | linked list, heuristics | 5 |

Open Pit Mining | Given an array of positive and negative expected returns, find a contiguous block that maximizes the expected profit. | divide-and-conquer, analysis of algorithms | 5 |

Baseball Elimination | Given the standings of a sports league, determine which teams are mathematically eliminated. | reduction, max flow, min cut | 3 |

Assignment Problem | Solve the assignment problem by reducing it to min cost flow. | reduction, min cost flow | 3 |

Password Cracking | Crack a subset-sum password authentication scheme. | hashing, space-time tradeoff | 7 |

TEXT PROCESSING | |||

Natural Language Modeling | Create a Markov model of an input text and use it to automatically generate stylized pseudo-random text. | suffix sorting or hashing | 6 |

Natural Language Modeling | Create a Markov model of an input text and use it to automatically generate stylized pseudo-random text. | Markov chains, graph | 4 |

Markovian Candidate [checklist] | Create a Markov model of an input text to perform speech attribution. | artificial intelligence, symbol table | 6 |

Word Searching | Search for words horizontally, vertically and diagonally in a 2D character array | tries | 7 |

Redundancy Detector | Find the longest repeated sequence in a given text. | suffix sorting, strings | 4 |

Text Indexing | Build an inverted index of a text corpus and find the position of query strings in the text. | suffix sorting or binary search tree | 4 |

COMMUNICATION | |||

Linear Feedback Shift Register | Encrypt images using a linear feedback shift register. | objects, encryption | 4 |

Pictures from Space | Detect and fix data errors in transmission using a Hadamard code. | 2D arrays, error-correcting codes | 3 |

Prefix Free Codes | Decode a message compressed using Huffman codes. | binary trees, data compression | 4 |

Burrows-Wheeler | Implement a novel text compression scheme that out-compresses PKZIP. | suffix sorting, arrays, data compression | 7 |

RSA Cryptosystem | Implement the RSA cryptosystem. | big integers, repeated squaring, analysis of algorithms | 8 |

DISCRETE MATH | |||

Linked List Sort | Shellsort a linked list. | linked list, shellsort | 4 |

Batcher Sort | Implement Batcher's even-odd mergesort. | divide-and-conquer, parallel sorting hardware | 6 |

Rational Arithmetic | Implement a Rational number data type. | struct, data abstraction, Euclid's algorithm | 3 |

Factoring | Factor large integers using Pollard's rho method. | big integers, Euclid's algorithm | 5 |

Deques and Randomized Queues | Create deque and randomized queue ADTs. | abstract data types, generics | 5 |

Linear Congruential Random Number Generator | Find the cycle length of a pseudo-random number generator using Floyd's algorithm. | loops, mod | 2 |

Stock Market | Predict the performance of a stock using Dilbert's rule. | loops | 2 |

Subset Sum | Partition the square roots of 1 to 100 into two subsets so that their sum is as close as possible to each other. | various | 6 |

Loops and Conditionals | Binary logarithm, checkerboard pattern, random walk, Gaussian distribution. | loops and conditionals | 1 |

Here are some Nifty Assignments created by instructors at other universities. They are more oriented towards recreational applications, but are fun and creative.

The purpose of this project is to re-affirm the atomic nature of matter by tracking the motion of particles undergoing Brownian motion, fitting this data to Einstein's model, and estimating Avogadro's number.

Perspective The atom played a significant role in 20th century physics and chemistry, but prior to 1908 the reality of atoms and molecules was not universally accepted. In 1827, the botanist Robert Brown observed the random erratic motion of wildflower pollen grains immersed in water using a microscope. This motion would later become known as Brownian motion. Einstein hypothesized that this Brownian motion was the result of millions of tiny water molecules colliding with the larger pollen grain particles.

In one of his "miracle year" (1905) papers, Einstein formulated a quantitative theory of Brownian motion to justify the "existence of atoms of definite finite size". His theory provided experimentalists with a method to count molecules with an ordinary microscope by observing their collective effect on a larger immersed particle. In 1908, Jean Baptiste Perrin used the recently invented ultra-microscope to experimentally validate Einstein's kinetic theory of Brownian motion, thereby providing the first direct evidence supporting the atomic nature of matter. His experiment also provided one of the earliest estimates of Avogadro's number. For this work, Perrin won the 1926 Nobel Prize in physics.

The Problem In this project, you will redo a version of Perrin's experiment. Your job is greatly simplified because with modern video and computer technology (in conjunction with your programming skills), it is possible to accurately measure and track the motion of an immersed particle undergoing Brownian motion. We supply video microscopy data of polystyrene spheres ("beads") suspended in water, undergoing Brownian motion. Your task is to write a program to analyze this data, determine how much each bead moves between observations, fit this data to Einstein's model, and estimate Avogadro's number.

Here is a movie (avi, mov) of several beads undergoing Brownian motion. Below is a typical raw image (left) and a cleaned-up version (right) using thresholding, as described below.

Page 2 of 8

Each image shows a two-dimensional cross section of a microscope slide. The beads move in and

out of the microscope's field of view (the jQuery20008039107078491861_1512323871465??- and ????-directions). Beads also move in the ????-

direction, so they can move in and out of the microscope's depth of focus; this results in halos,

and it can also result in beads completely disappearing from the image.

Problem 1. (Particle Identification) The first challenge is to identify the beads amidst the noisy

data. Each image is 640-by-480 pixels, and each pixel is represented by a Color object which

needs to be converted to a luminance value ranging from 0.0 (black) to 255.0 (white). Whiter

pixels correspond to beads (foreground) and blacker pixels to water (background). We break the

problem into three pieces:

1. Read the image. Use the Picture data type to read in the image.

2. Classify the pixels as foreground or background. We use a simple, but effective, technique

known as thresholding to separate the pixels into foreground and background components:

all pixels with monochrome luminance values strictly below some threshold ????(tau) are

considered background, and all others are considered foreground. The two pictures in figure

above illustrate the original frame (left) and the same frame after thresholding (right), using

???? = 180.0. This value of ???? results in an effective cutoff for the supplied data.

3. Find the blobs. A polystyrene bead is typically represented by a disc-like shape of at least

some minimum number ???? (typically 25) of connected foreground pixels. A blob or connected

component is a maximal set of connected foreground pixels, regardless of its shape or size.

Page 3 of 8

We will refer to any blob containing at least ???? pixels as a bead. The center-of-mass of a blob

(or bead) is the average of the x- and y-coordinates of its constituent pixels.

Define a helper data type Blob in blob.py that has the following API:

Method Description

Blob() An empty blob ????

b.add(i,j) Add pixel (????, ????) to the ????

b.mass() The number of pixels in ????, i.e., its mass

b.distanceto(c) The distance between the centers of ???? and ????

str(b) String representation of ????'s mass and center of mass

Next, define a data type BlobFinder in blob_finder.py that has the following API. Use

depth-first search to efficiently identify the blobs.

Method Description

BlobFinder(pic, tau) A blob finder bf to find blobs in the picture pic using a

luminance threshold tau

bf.getBeads(P) List of all beads with ≥ ???? pixels

Finally, implement a test client _main() in blob_finder that takes an integer ????, a float tau,

and the name of a JPEG file as command-line arguments; then creates a BlobFinder object

using a luminance threshold of tau; prints out all of the beads with at least P pixels; and finally,

prints out all of the blobs (beads with at least 1 pixel). Note: For full credit you must write

_main() as if it were a client. That means that you must not access any private methods or

private instance variables of the BlobFinder data type.

Page 4 of 8

The program identifies 13 blobs in the sample frame, 7 of which are beads. Our string

representation of a blob specifies its mass (number of pixels) and its center of mass (in the 640-

by-480 picture). By convention, pixels are measured from left-to-right, and from top-to-bottom

(instead of bottom-to-top).

Problem 2. (Particle Tracking) The next step is to determine how far a bead moved from onetime

step ???? to the next ???? + Δ????. For our data, Δ???? = 0.5 seconds per frame. We assume the data is

such that each bead moves a relatively small amount, and that two beads do not collide. However,

we must account for the possibility that the bead disappears from the frame, either by

departing the microscope's field of view in the ???? - or ???? -directions, or moving out of the

microscope's depth of focus in the ????-direction. Thus, for each bead at time ???? + Δ????, we calculate

the closest bead at time ???? (in Euclidean distance) and identify these two as the same beads.

$ python blob_finder.py 25 180.0 data/run_1/frame00000.jpg

7 Beads:

37 (220.0270, 122.8919)

36 (297.8333, 394.5000)

39 (312.3077, 215.8205)

31 (433.7742, 375.4839)

32 (475.5000, 44.5000)

31 (525.2903, 443.2903)

35 (632.7714, 154.5714)

13 Blobs:

37 (220.0270, 122.8919)

1 (254.0000, 223.0000)

17 (255.4118, 233.8824)

23 (265.8261, 316.4348)

36 (297.8333, 394.5000)

39 (312.3077, 215.8205)

23 (373.0000, 357.1739)

19 (390.8421, 144.8421)

31 (433.7742, 375.4839)

32 (475.5000, 44.5000)

31 (525.2903, 443.2903)

24 (591.0000, 399.5000)

35 (632.7714, 154.5714)

Page 5 of 8

However, if the distance is too large, i.e., greater than Δ(delta) pixels, we assume that one of

the beads has either just begun or ended its journey. We record the displacement that each

bead travels in the Δ???? units of time.

Implement a client program bead_tracker.py that takes an integer ????, a float tau, a float

delta, and a sequence of JPEG filenames as command-line arguments, identifies the beads in

each JPEG image using BlobFinder, and prints out (one per line, formatted with 4 decimal

places to the right of decimal point) the radial distance that each bead moves from one frame

to the next (assuming it is no more than delta). Note that it is not necessary to explicitly track a

bead through a sequence of frames - you only need to worry about identifying the same bead in

two consecutive frames.

Problem 3. (Data Analysis) Einstein's theory of Brownian motion connects microscopic properties

(e.g., radius, diffusivity) of the beads to macroscopic properties (e.g., temperature, viscosity) of

the fluid in which the beads are immersed. This amazing theory enables us to estimate

Avogadro's number with an ordinary microscope by observing the collective effect of millions

of water molecules on the beads.

1. Estimating the self-diffusion constant. The self-diffusion constant ???? characterizes the

stochastic movement of a molecule (bead) through a homogeneous medium (the water

molecules) because of random thermal energy. The Einstein-Smoluchowski equation states

that the random displacement of a bead in one dimension has a Gaussian distribution with

mean zero and variance ????2 = 2????Δ???? , where Δ???? is the time interval between position

measurements. That is, a molecule's mean displacement is zero and its mean square

$ python bead_tracker.py 25 180.0 25.0 data/run_1/frame00000.jpg

data/run_1/frame00001.jpg

7.1833

4.7932

2.1693

5.5287

5.4292

4.3962

Page 6 of 8

displacement is proportional to the elapsed time between measurements, with the constant

of proportionality 2???? . We estimate ???????? by computing the variance of all observed bead

displacements in the ???? and ???? directions. Let (Δ????1, Δ????1), ... , (Δ????????, Δ????????) be the ???? bead

displacements, and let r1, ... , ???????? denote the radial displacements. Then

????2 =

(Δ????1

2 + ⋯ + Δ????????

2) + (Δ????1

2 + ⋯ + Δ????????

2)

2????

=

r1 + ⋯ + ????????

2????

.

For our data, we have Δ???? = 0.5 . So, our estimate for ????2 is an estimate for ???? as well

(????2 = 2 × ???? × Δ???? = 2 × ???? × 0.5 = ????). Note that the radial displacements in the formula

above are measured in meters. The radial displacements output by your

bead_tracker.py program is measured in pixels. To convert from pixels to meters,

multiply by ????. ???????????? × ????????−???? (meters per pixel). The value of ???? is the count of the total

number of displacements read.

2. Estimating the Boltzmann constant. The Stokes-Einstein relation asserts that the selfdiffusion

constant ???? of a spherical particle immersed in a fluid is given by

???? =

????????

6????????????

,

where, for our data ???? (absolute temperature) is 297 degrees Kelvin (room temperature), ????

(eta, viscosity of water) is 9.135 × 10−4 Nsm−2 (at room temperature), ???? (rho, radius of

bead) is 0.5 × 10−6, and ???? is the Boltzmann constant. All parameters are given in SI units. The

Boltzmann constant is a fundamental physical constant that relates the average kinetic

energy of a molecule to its temperature. We estimate ???? by measuring all the parameters in

the Stokes-Einstein equation, and solving for ????.

3. Estimating Avogadro's number. Avogadro's number ???????? is defined to be the number of

particles in a mole. By definition, ???? = ????/???????? , where the universal gas constant ???? is

approximately 8.31457 JK−1mol−1. Use ????/???? as an estimate of Avogadro's number.

Page 7 of 8

For the final part, implement a client program avogadro.py that reads in the displacements

from standard input and computes an estimate of Boltzmann's constant and Avogadro's

number using the formulae described above.

Data We provide ten datasets (they are under the data directory), obtained by William Ryu

(Princeton University) using fluorescent imaging. Each run contains a sequence of two hundred

640-by-480 color JPEG images, frame00000.jpg through frame00199.jpg and is stored

in a subdirectory run_1 through run_10. The directory also contains some reference solutions.

Files to Submit

1. blob.py

2. blob_finder.py

3. bead_tracker.py

4. avogadro.py

5. Report.pdf

Before You Submit:

• Your programs must be fully commented.

• Make sure your programs meet the input and output specifications by running the

following command on the terminal:

where the optional argument <problems> lists the problems (Problem1, Problem2, etc.)

you want to test; all the problems are tested if no argument is given.

• Make sure your programs meet the style requirements by running the following

command on the terminal:

$ python bead_tracker.py 25 180.0 25.0 data/run_1/frame00000.jpg

data/run_1/frame00001.jpg | python avogadro.py

Boltzman = 1.173701e-23

Avogadro = 7.084062e+23

$ python run_tests .py -v [< problems>]

$ pep8 <program>

Page 8 of 8

where <program> is the .py file which style you want to check.

• Create a folder and named as your <first name>_<lat name>. Then, compress that folder

to create <first name>_<lat name>.zip file.

• Open a blank document in MS Word. Save your report as Problem<number>.pdf once you

finish your report. Make sure your report is not too verbose, does not contain

spelling/grammatical mistakes.

Acknowledgements This project is an adaptation of the Atomic Nature of Matter assignment

developed at Princeton University by Thomas Clarke, Robert Sedgewick, Scott Vafai, and Kevin

Wayne.

## 0 thoughts on “Brownian Motion Programming Assignment”

-->