33. Monte Carlo Methods 2 by MIT OpenCourseWare

Summary by www.lecturesummary.com: 1. Introduction, Optimization Problems (MIT 6.0002 Intro to Computational Thinking and Data Science) by MIT OpenCourseWare

• 1. Introduction and Homework Questions

  • Lecture starts with thanks to Creative Commons License and sponsors of MIT Open Courseware. The subject remains Fantastics.
  • Instructor asks whether there are homework questions, hoping students are considering it.
  • Teacher interprets the silence as the homework being "entirely in hand."

  2. Time-Dependent Probability Distributions

  • Lecture transitions to time-dependent probability distributions.
  • Experimental predictions often involve knowledge of how the probability of observables evolves with time, expressed as dP/dt.
  • Example: A box containing a radioactive uranium atom. The chances of the atom being inside the box vary with time.

  3. Macroscopic vs. Microscopic Thinking

  • A plausible equation for the probability of the uranium atom being in the box resembles exponential decay (P(t) = P_initial * e^(-t/tau)).
  • This equation is incorrect for a single atom as it treats the number of atoms as a continuous variable.
  • We are accustomed to macroscopic thinking, dealing with variables like mass as if they were continuous, despite materials being made of discrete atoms.
  • With one atom, it is either present or not (1 or 0). The decay influences the probability that it remains present.
  • In an experiment with a million boxes, each starting with one uranium atom, the observed number of atoms would show discrete drops as atoms decay.
  • This observed macroscopic behavior is essentially sampling from the probability distribution for a single box a million times.

  4. The Need for Discrete (Probability) Equations in Small Systems

  • Equations like dP/dt for complex systems will still include the subtlety between continuous probabilities and discrete experimental results.
  • Each experiment samples from the time-varying probability distribution, leading to different outcomes due to the system's intrinsic randomness.
  • Instructor advises reading Joe Scott's notes on chemical kinetics for further understanding.
  • Chemical kinetic equations under standard circumstances are also approximations since they treat discrete variables as continuous.
  • Continuum equations are used extensively since, in large systems, the contribution of one or two atoms is negligible.
  • Quantum mechanics accounts for everything through wave functions, and actual equations focus on probabilities.
  • This issue is more critical for smaller systems containing fewer molecules.
  • This applies in biology (within cells), single molecule spectroscopy, and emulsion polymerization.

  5. The Kinetic Master Equation

  • To describe such systems, we must consider graininess and set up the equation for dP/dt.
  • The probability P of interest often represents the probability of having a particular number of molecules for each type in the system.
  • The kinetic master equation models the evolution of the probability of being in a particular state as a result of reactions.
  • For a reaction A + B → C, the chance of state (NA, NB, NC) decreases due to the forward reaction.
  • The likelihood of state (NA, NB, NC) increases when another state can interact to yield this state.
  • This collection of equations is the kinetic master equation, representing actual probability evolution.
  • The master equation is a linear differential equation for the probabilities P.
  • The size of the vector P can be enormous, scaling with the maximum number of molecules considered.
  • If the system has over 10^8 possible states, exact solution of the master equation is not feasible.

  6. Stochastic Solutions (Kinetic Monte Carlo)

  • If the master equation is too extensive, a stochastic solution method should be utilized to sample from the probability distribution P(t).
  • Chemical kinetics generally entails one step at a time; the chance of two reactions occurring simultaneously is negligible.
  • For reactions like A + A → B, the rate in continuum kinetics goes as the square of the number of A molecules.
  • The units of rate constants in the master equation differ from usual macroscopic kinetics.
  • The master equation rate constant for a bimolecular reaction is the normal K/volume.
  • Transport may also be incorporated into the master equation by reducing macroscopic transport rates to microscopic time constants.

  • 7. Kinetic Monte Carlo Method

    The method to solve large systems is called kinetic Monte Carlo, invented by Gillespie.

    8. Goal of Kinetic Monte Carlo

    The goal is to sample from the true probability density P(t) without explicitly computing it. This is analogous to sampling from the Boltzmann distribution in Monte Carlo Metropolis without explicitly writing down the density.

    9. Approach

    The approach is to trace individual trajectories of the system. At any moment, for the actual state (NA, NB, NC), we examine what may occur in a small time interval. We select the time interval so that either nothing happens or only one event does (e.g., one reaction).

    10. Uranium Example

    Example: Begin with a single atom. Compute the probability that it decays during a small time interval. Generate a random number to determine whether it decayed or not. If not, try again. If yes, it's lost. This forms a single possible trajectory. Repeating this many times provides different possible paths (e.g., lived 3 years and 14 minutes), which sample from the actual probability density.

    11. The Gillespie Algorithm

    The Gillespie algorithm calculates the:

    • Time until the next event
    • What that event is

    12. Steps of the Gillespie Algorithm

    1. Calculate the rates of all conceivable events that may occur from the current state.
    2. Add up these rates to obtain the total rate 'a', the expected rate per second for any event to occur. In the example of uranium with 100 atoms, 'a' is 100 times as fast as the decay of one atom. The greater the number of molecules you have, the greater the overall rate 'a', and the sooner something will be predicted to occur. This overall rate 'a' must be recalculated every time an event occurs since the state (counts of molecules) has altered.
    3. Determine the time step (delta t, also referred to as tau) to the next event occurrence. This duration is calculated using a random number R1 (0 to 1) and the sum rate 'a' through the equation: delta t = (1/a) * log(1/R1). This equation yields a longer period if 'a' is small (less probable for an event occurrence) and shorter if 'a' is big.
    4. Find which particular event occurred by taking another random number R2 (between 0 and 1). Enumerate all possible events (e.g., forward reaction A+B->C, reverse reaction C->A+B) and their rates.
    5. The probability of a particular process (Process 1) occurring is its rate divided by the overall rate 'a'. If R2 is below this ratio of probabilities for Process 1, then Process 1 occurs. If you have more than one process, you divide up the range proportionally with regard to relative probabilities (rates/a) and determine into which range R2 belongs.
    6. Change the system state (the counts of molecules NA, NB, NC) based on the selected event. For instance, if the forward reaction A+B->C occurred, reduce NA and NB by one and increase NC by one.
    7. Repeat the process: Recalculate the total rate 'a' in terms of the new state, calculate the time to the next event, find which event occurs, change the state, and so on. This produces one trajectory or simulation run.

    13. With Kinetic Monte Carlo

    • One simulation run gives one sample trajectory of what may occur.
    • To interpret the probability distribution P(t), one needs to perform multiple simulations (trajectories) from the initial conditions. Storing the outcome of the simulations (states at various times) permits the construction of histograms and the study of the probability distribution. These trajectories are draws from P(t).
    • Kinetic Monte Carlo (the Gillespie algorithm) is relatively straightforward (a three-step process carried out by the computer) and widely applied, even in problems where the master equation might be solved analytically, as it is simpler to do.